The solver modules consists of several modules:
The easiest way to use the solver module is to setup all configuration and use the generic solve function to solve a stencil-based linear equation system A=b with the help of the CG method or Chebyshev iteration and a preconditioner.
In MPIOM, the solver configuration in generated in mo_solver_hlp.f90
. There, the read_solver_config
function reads the configuration from file SLVCTL and fills the type solver_config with it. SLVCTL is generated by the preparte_job_mpiom_hamocc
script. It comes with detailed comments on its usage and the meaning of the various parameters. Additionally, this configuration is checked (check_solver_config
), broadcast (broadcast_solver_config
) and logged (logwrite_solver_config
). With the help of the functions get_interior
and set_extents
the extents of the solution vector z1o is set.
During the startup of MPIOM in mpiom.f90
, a solver_config
is generated and passed to init_solver to initialize the solver library with the user-set configuration. The actual utilization of the solver takes place in mo_tro.f90
in the function solve_barotro_sys
. There, first the calc_barotro_rhs
function is called to calculate the right-hand-side of the barotropic system then stencil (i.e. the linear system) of the barotropic system is set via set_stencil. Next, the checkprep_params
from mo_solver_hlp
function checkes and prepares all parameters by looking at the parameters set by the solver initialization and complementing them. For instance the eigenvalues lambda_min and max can be set in SLVCTL or are retrieved automatically during this preparation step which is executed only once. After this a call to the solve function solves the system. As parameters this function needs the barotropic stencil operation as given by apply_stencil_dp , the right-hand-side b1o, the starting value z1o with its extents ext_z1o as well as a function for the boundary exchange (bounds_exchange_dp
) which is defined in mo_solver_hlp
in MPIOM. Since almost all functions are defined for single (suffix _sp) and double (suffix _dp) precision we pass the double precision variants of all funtions to solve. The corresponding code is:
The variable z1o now holds the result and kiter the needed iterations. As stated before, the solve function is just a convenient interface to the precond_cg_method and the precond_chebyshev_method as well as all preconditioners in preconditioners module. Direct usage of these two solvers is pretty self-explanatory and can be looked up easily in this reference.
The additive Schwarz method schwarz_method resides in the solvers module. This method uses a local method to solve the global problem first on each partition seperately with fixed Dirichlet boundaries. Then, the extended overlapping boundary halos are exchanged and another local problem is solved. This step is repeated until the solution of the global problem has converged.
To use this method, we first need to set up arrays with extended boundary halos. This is done for instance via the setup_halos
function (which then uses sethaloN
) in mo_tro.f90
of MPIOM. This function enlarges the defining arrays of the barotropic system, i.e., uf, vf, ff, as well as z1o and b1o. The next step is to use the enlarged fields zuf, zvf, zff to initialize the matrix stencil with set_stencil. Then all parameters are checked again via checkprep_params
in MPIOM. After this we can directly call schwarz_method with the function reference apply_stencil. The following code snipped illustrates this.
It is also possible to use the CG method in single precision while obtaining double precision accuracy in the end. The hope is that solving in single precision (sp) is faster than in double precision (dp) because twice as many sp numbers can be transferred in the same time as db numbers. Also, if the VMC (on PowerPC) or SSE units (on Intel) is used, vector operations in sp can be twice as fast as in dp. On the POWER6 architecture the VMX units as well as the compiler is quite limited according to IBM but the future POWER7+ with its VSX units will drastically increase the sp performace compared to dp.
To make future use of this approach, nearly all functions and routines in the solver library have a single and double precision variant. This was accomplished by splitting each module into a declaration file, e.g. solvers.f90
and an implementation file, e.g. solvers_multi.f90
. Preprocessor macros take care that the implementation file is included into the declaration file with the corresponding data kinds set accordingly. Also, each name of a multi-precision function gets the suffix _sp resp. _dp for single precision resp. double precision.
The iter_refinement implements the multi-precision iterative refinement. It takes as parameters the stencil operations in single and double precision. The code below illustrates the usage taken from solve_barotro_sys
.
Das diesem Bericht zugrundeliegende Vorhaben wurde mit Mitteln des Bundesministeriums für Bildung, und Forschung unter dem Förderkennzeichen 01IH08004E gefördert. Die Verantwortung für den Inhalt dieser Veröffentlichung liegt beim Autor.