官方建议的ARKODE编程步骤

原创于 2022-03-10 11:13:38 发布 · 344 阅读

0 ·

CC 4.0 BY-SA版权

文章标签：

#sundials

SUNDIAL 专栏收录该内容

27 篇文章

订阅专栏

本文详细介绍了如何初始化并配置多尺度模拟的环境，包括设置问题尺寸、创建初始值向量、创建并配置快速和慢速积分器，以及设置步长、指定根查找问题、时间推进和资源释放等关键步骤。内容涵盖了MPI并行环境的初始化、NVECTOR库的使用，以及ARKStep和MRIStep这两个数值积分器的接口调用和参数设置。

                    
                    Initialize parallel or multi-threaded environment, if appropriate.

For example, call MPI_Init to initialize MPI if used, or set num_threads, the number of threads to use

within the threaded vector functions, if used.

Set problem dimensions, etc.

This generally includes the problem size, N, and may include the local vector length Nlocal.

Note: The variables N and Nlocal should be of type sunindextype.

Set vector of initial values

To set the vector y0 of initial values, use the appropriate functions defined by the particular NVECTOR implementation.

For native SUNDIALS vector implementations (except the CUDA and RAJA based ones), use a call of the form

y0 = N_VMake_***(…, ydata);

if the realtype array ydata containing the initial values of 𝑦 already exists. Otherwise, create a new vector

by making a call of the form

y0 = N_VNew_***(…);

and then set its elements by accessing the underlying data where it is located with a call of the form

ydata = N_VGetArrayPointer_***(y0);

See the sections The NVECTOR_SERIAL Module through The NVECTOR_PTHREADS Module for details.

For the HYPRE and PETSc vector wrappers, first create and initialize the underlying vector, and then create the

NVECTOR wrapper with a call of the form

y0 = N_VMake_***(yvec);

where yvec is a HYPRE or PETSc vector. Note that calls like N_VNew_***(…) and

N_VGetArrayPointer_***(…) are not available for these vector wrappers. See the sections The

NVECTOR_PARHYP Module and The NVECTOR_PETSC Module for details.

If using either the CUDA- or RAJA-based vector implementations use a call of the form

y0 = N_VMake_***(…, c);

where c is a pointer to a suncudavec or sunrajavec vector class if this class already exists. Otherwise,

create a new vector by making a call of the form

N_VGetDeviceArrayPointer_***

or

N_VGetHostArrayPointer_***

Note that the vector class will allocate memory on both the host and device when instantiated. See the sections

The NVECTOR_CUDA Module and The NVECTOR_RAJA Module for details.

Create an ARKStep object for the fast (inner) integration

Call inner_arkode_mem = ARKStepCreate(…) to create the ARKStep memory block.

ARKStepCreate() returns a void* pointer to this memory structure. See the section ARKStep initialization

and deallocation functions for details.

Configure the fast integrator

Specify tolerances, create and attach matrix and/or solver objects, or call ARKStepSet* functions to configure

the fast integrator as desired. See sections A skeleton of the user’s main program and Optional input functions

for details on configuring ARKStep.

Notes on using ARKStep as a fast integrator:

If the inner method is not explicitly specified then the default method in the ARKStep module will be used. If

a particular fast method is desired it should be set in this step. The slow method can be set when configuring the

slow integrator in the following steps.

By default the fast integrator will use adaptive step sizes. To use a fixed fast step a call to

ARKStepSetFixedStep() should be made in this step otherwise fast integration tolerances should be set

in this step as described in A skeleton of the user’s main program.

If attaching a user_data pointer, it should be attached to the slow integrator in the following steps with

MRIStepSetUserData(). This pointer will subsequently be passed to user-provided functions during the

fast integration.

Specifying a rootfinding problem for the fast integration is not supported. Rootfinding problems should be

created and initialized with the slow integrator. See the steps below and MRIStepRootInit() for more

details.

The ARKStep module used for the fast time scale must be configured with an identity mass matrix.

Create an MRIStep object for the slow (outer) integration

Call arkode_mem = MRIStepCreate(…) to create the MRIStep memory block.

MRIStepCreate() returns a void* pointer to this memory structure. See the section MRIStep initialization and deallocation functions for details.

Set the slow step size

Call MRIStepSetFixedStep() to specify the slow time step size.

Set optional inputs

Call MRIStepSet* functions to change any optional inputs that control the behavior of MRIStep from their

default values. See the section Optional input functions for details.

Specify rootfinding problem

Optionally, call MRIStepRootInit() to initialize a rootfinding problem to be solved during the integration

of the ODE system. See the section Rootfinding initialization function for general details, and the section

Optional input functions for relevant optional input calls.

Advance solution in time

For each point at which output is desired, call

ier = MRIStepEvolve(arkode_mem, tout, yout, &tret, itask);

Here, itask specifies the return mode. The vector yout (which can be the same as the vector y0 above) will

contain 𝑦(𝑡out). See the section MRIStep solver function for details.

Get optional outputs

Call MRIStepGet* and/or ARKStepGet* functions to obtain optional output from the slow or fast integrators

respectively. See the section Optional output functions and Optional output functions for details.

Deallocate memory for solution vector

Upon completion of the integration, deallocate memory for the vector y (or yout) by calling the NVECTOR

destructor function:

N_VDestroy(y);

Free solver memory

Call ARKStepFree(&inner_arkode_mem) and MRIStepFree(&arkode_mem) to free

the memory allocated for fast and slow integration modules respectively.

Free linear solver and matrix memory

Call SUNLinSolFree() and (possibly) SUNMatDestroy() to free any memory allocated for

any linear solver and/or matrix objects created above for the fast integrator.

Finalize MPI, if used

Call MPI_Finalize to terminate MPI.