version 0.6.0
Linear system solvers

Topics

 HYPRE linear system interface
 Massively parallel HYPRE library interface routines.
 
 LIS linear system interface
 Parallel LIS library interface routines.
 
 MUMPS linear solver library
 Interface for the MUMPS linear solver library.
 

Classes

type  type_solver::t_base_solver
 Parent type for all solvers. More...
 
type  type_solver::t_direct_solver
 Parent type for all direct solvers. More...
 
type  type_solver::t_iterative_solver
 Parent type for all iterative solvers. More...
 

Functions

subroutine mod_linear_system_solver_cell_interface::linear_system_solver_cell_interface (matrix, sol, rhs, solver, boundary_condition, stencil, ls_map, equation_name, coupled_equations, sol2, boundary_condition2, is_print_residual)
 Interface routine for cell solvers. More...
 
subroutine mod_linear_system_solver_face_interface::linear_system_solver_face_interface (matrix, sol, rhs, solver, boundary_condition, face_stencil, face_ls_map, cell_ls_map)
 Fortran interface routine for face solvers (HYPRE SSTRUCT interface, Lis, Notus and MUMPS direct solvers). More...
 
subroutine mod_scaled_residual::compute_scaled_residual_cell (matrix, rhs, solution, ls_map, stencil, scaled_residual)
 Compute the scaled residual of a cell-based linear system. More...
 
subroutine mod_scaled_residual::compute_scaled_residual_face (matrix, rhs, solution, ls_map, stencil, scaled_residual)
 Compute the scaled residual of a face-based linear system. More...
 

Detailed Description

This directory provides a consistent interface to linear system solvers, which can be either built-in or third-party backend. A derived type system selects the solver to use, carry solver-specific parameters and carry result information. The routines linear_system_solver_cell_interface and linear_system_solver_face_interface solves a linear system by dispatching to the appropriate solver.

Currently, the built-in backends are:

and the third-party backends are:

The table below presents available solvers for each interface

Backend linear_system_solver_cell_interface linear_system_solver_face_interface
Hypre_Struct Yes No
Hypre_SStruct Yes Yes
Mumps Yes Yes
Notus BiCGStab Yes Yes
Lis Yes Yes

Common interface for solvers

Notus uses inheritance of derived types to determine which solver to use. The two fundamental types are the type_solver::direct_solver type and the type_solver::iterative_solver type, which both are derived from the type_solver::base_solver type.

base_solver
 │
 ├── direct_solver
 │
 └── iterative_solver

Every type_solver::direct_solver has a residual attribute, which receives the solver's result residual.

Every type_solver::iterative_solver has max_iter and tolerance parameters, which holds the maximum number of iterations and the residual tolerance desired. The derived type also has a final_iter and a residual argument, which receives the actual number of iteration and the residual after the last iteration.

Notus' Jacobi preconditioner

Notus has its own preconditioner that can be applied prior to any solver backend. This is a left Jacobi preconditioner that simply divides each row by the diagonal coefficient.

Every solver has a is_jacobi_initial_preconditioner attribute to activate or deactivate the Notus' Jacobi preconditioner.

Function Documentation

◆ compute_scaled_residual_cell()

subroutine mod_scaled_residual::compute_scaled_residual_cell ( double precision, dimension(:), intent(in)  matrix,
double precision, dimension(:), intent(in)  rhs,
double precision, dimension(:,:,:), intent(in)  solution,
type(t_ls_map), intent(in)  ls_map,
type(t_cell_stencil), intent(in)  stencil,
double precision, intent(out)  scaled_residual 
)

The scaled residual is equal to

\[ \frac{|Ax - b|_\infty}{|A|_\infty|x|_\infty} \]

where the infinity norm of the matrix \( A \) is equal to

\[ \max\left{\sum_{j=1}^n |A_{ij}| \mid| i \in [1,n] \right} \]

Parameters
[in]matrixmatrix of the linear system.
[in]rhsright-hand side of the linear system (serialized).
[in]solutionsolution of the linear system in (i,j,k) form.
[in]ls_maplinear system mapping structure.
[in]stencilstencil of the discretization.
[out]scaled_redidualscaled residual.

◆ compute_scaled_residual_face()

subroutine mod_scaled_residual::compute_scaled_residual_face ( double precision, dimension(:), intent(in)  matrix,
double precision, dimension(:), intent(in)  rhs,
type(t_face_field), intent(in)  solution,
type(t_face_ls_map), intent(in)  ls_map,
type(t_face_stencil), intent(in)  stencil,
double precision, intent(out)  scaled_residual 
)

The scaled residual is equal to

\[ \frac{|Ax - b|_\infty}{|A|_\infty|x|_\infty} \]

where the infinity norm of the matrix \( A \) is equal to

\[ \max\left{\sum_{j=1}^n |A_{ij}| \mid| i \in [1,n] \right} \]

Parameters
[in]matrixmatrix of the linear system.
[in]rhsright-hand side of the linear system (serialized).
[in]solutionsolution of the linear system in (i,j,k) form.
[in]ls_maplinear system mapping structure.
[in]stencilstencil of the discretization.
[out]scaled_redidualscaled residual.

◆ linear_system_solver_cell_interface()

subroutine mod_linear_system_solver_cell_interface::linear_system_solver_cell_interface ( double precision, dimension(:), intent(inout)  matrix,
double precision, dimension(:,:,:), intent(inout)  sol,
double precision, dimension(:), intent(inout)  rhs,
class(t_base_solver), intent(inout)  solver,
type(t_boundary_condition), intent(in)  boundary_condition,
type(t_cell_stencil), intent(in)  stencil,
type(t_ls_map), intent(in)  ls_map,
character(len=*), intent(in)  equation_name,
logical, intent(in), optional  coupled_equations,
double precision, dimension(:,:,:), intent(inout), optional  sol2,
type(t_boundary_condition), intent(in), optional  boundary_condition2,
logical, intent(in), optional  is_print_residual 
)

This routine accepts the following classes of solvers:

  • t_hypre_solver_Struct_interface
  • t_hypre_solver_SStruct_interface
  • lis_solver
  • t_mumps_solver
  • notus_solver_mg

The following steps are done:

  • sol (array of dimension 3) is put in one dimension array solution
  • Point Jacobi preconditioner is applied (left preconditioning)
  • solvers are called
  • solution is put in sol
  • boundary condition are applied on boundaries ghost cells
  • exchange beetween processors is done
Parameters
[in,out]matrixmatrix of the linear system (vectorized)
[in,out]solsolution of the linear system (not vectorized)
[in,out]rhsright-hand-side of the linear system (vectorized)
[in,out]solversolver structure associated to the equation
[in]boundary_conditionboundary condition
[in]stencilstencil structure of the equation
[in]ls_mapmapping associated to the cell grid
[in]equation_namename of the equation (log message)
[in]coupled_equationsoptional logical to solve a 2-fields coupled linear system
[in,out]sol2optional solution of the second unknown
[in]boundary_condition2optional boundary condition associated to the second unknown
Note
The sol argument does not need to be vectorized at upper level, so it should be provided as a rank-3 array, and the vectorization occurs in the routine. On the opposite, rhs array is vectorized at a upper level.

◆ linear_system_solver_face_interface()

subroutine mod_linear_system_solver_face_interface::linear_system_solver_face_interface ( double precision, dimension(:), intent(inout)  matrix,
type(t_face_field), intent(inout)  sol,
double precision, dimension(:), intent(inout)  rhs,
class(t_base_solver), intent(inout)  solver,
type(t_boundary_condition_face), intent(in)  boundary_condition,
type(t_face_stencil), intent(in)  face_stencil,
type(t_face_ls_map), intent(in)  face_ls_map,
type(t_ls_map), intent(in)  cell_ls_map 
)

The following steps are performed:

  • The rank 3 solution 'sol' (array of dimension 3) is copied into a rank 1 array.
  • Jacobi preconditioner is applied (left preconditioning).
  • Linear solvers are called.
  • The rank 1 solution is copied into the 'sol' rank 3 array.
  • Boundary conditions are applied on ghost boundary faces.
  • The solution is exchanged between processes.
Parameters
[in,out]matrixmatrix of the linear system (vectorized)
[in,out]solsolution of the linear system (not vectorized)
[in,out]rhsright-hand side of the linear system (vectorized)
[in,out]solverlinear solver
[in]boundary_conditionboundary conditions on faces
[in]face_stencilstencil structure
[in]face_ls_mapmapping associated to the face grid for HYPRE solver
[in]cell_ls_mapmapping associated to cell grid

Note: solution do not need to be vectorized at a upper level. This is the reason why the corresponding argument is a rank 3 array, that is vectorized inside this routine. On the opposite, rhs array is vectorized at a upper level.