8.1
general documentation
cs_cdofb_monolithic_sles.h File Reference
#include "cs_defs.h"
#include "cs_cdofb_monolithic_priv.h"
+ Include dependency graph for cs_cdofb_monolithic_sles.h:

Go to the source code of this file.

Functions

cs_cdofb_monolithic_sles_tcs_cdofb_monolithic_sles_create (cs_lnum_t n_faces, cs_lnum_t n_cells)
 Create an empty cs_cdofb_monolithic_sles_t structure. More...
 
void cs_cdofb_monolithic_sles_clean (cs_cdofb_monolithic_sles_t *msles)
 Free a part of the structure. More...
 
void cs_cdofb_monolithic_sles_free (cs_cdofb_monolithic_sles_t **p_msles)
 Free memory related to cs_cdofb_monolithic_sles_t structure. More...
 
void cs_cdofb_monolithic_sles_init_sharing (const cs_mesh_t *mesh, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant)
 Set pointers to shared structures. More...
 
void cs_cdofb_monolithic_sles_finalize (void)
 Free if needed structure(s) associated CDO face-based schemes with a monolithic velocity-pressure coupling. More...
 
void cs_cdofb_monolithic_set_sles (cs_navsto_param_t *nsp, void *context)
 Start setting-up the Navier-Stokes equations when a monolithic algorithm is used to couple the system. No mesh information is available at this stage. nsp is not declared as const to avoid compilation warnings but it should be modified at this stage. More...
 
int cs_cdofb_monolithic_solve (const cs_navsto_param_t *nsp, const cs_equation_param_t *eqp, const cs_cdo_system_helper_t *sh, cs_param_sles_t *slesp, cs_cdofb_monolithic_sles_t *msles)
 Solve a linear system arising from the discretization of the Navier-Stokes equation with a CDO face-based approach. The full system is treated as one block and solved as it is. In this situation, PETSc or MUMPS are usually considered. More...
 
int cs_cdofb_monolithic_block_krylov (const cs_navsto_param_t *nsp, const cs_equation_param_t *eqp, const cs_cdo_system_helper_t *sh, cs_param_sles_t *slesp, cs_cdofb_monolithic_sles_t *msles)
 Solve a linear system arising from the discretization of the Navier-Stokes equation with a CDO face-based approach. The system is split into blocks to enable more efficient preconditioning techniques. The main iterative solver is a Krylov solver such as GCR, or MINRES. More...
 
int cs_cdofb_monolithic_gkb_solve (const cs_navsto_param_t *nsp, const cs_equation_param_t *eqp, const cs_cdo_system_helper_t *sh, cs_param_sles_t *slesp, cs_cdofb_monolithic_sles_t *msles)
 Use the GKB algorithm to solve the saddle-point problem arising from CDO-Fb schemes for Stokes and Navier-Stokes with a monolithic coupling. More...
 
int cs_cdofb_monolithic_uzawa_cg_solve (const cs_navsto_param_t *nsp, const cs_equation_param_t *eqp, const cs_cdo_system_helper_t *sh, cs_param_sles_t *slesp, cs_cdofb_monolithic_sles_t *msles)
 Use the preconditioned Uzawa-CG algorithm to solve the saddle-point problem arising from CDO-Fb schemes for Stokes, Oseen and Navier-Stokes with a monolithic coupling. This algorithm is based on Koko's paper "Uzawa conjugate gradient method for the Stokes problem: Matlab implementation with P1-iso-P2/ P1 finite element". More...
 
int cs_cdofb_monolithic_uzawa_al_incr_solve (const cs_navsto_param_t *nsp, const cs_equation_param_t *eqp, const cs_cdo_system_helper_t *sh, cs_param_sles_t *slesp, cs_cdofb_monolithic_sles_t *msles)
 Use the Uzawa algorithm with an Augmented Lagrangian (ALU) technique in an incremental way to solve the saddle-point problem arising from CDO-Fb schemes for Stokes, Oseen and Navier-Stokes with a monolithic coupling. More...
 

Function Documentation

◆ cs_cdofb_monolithic_block_krylov()

int cs_cdofb_monolithic_block_krylov ( const cs_navsto_param_t nsp,
const cs_equation_param_t eqp,
const cs_cdo_system_helper_t sh,
cs_param_sles_t slesp,
cs_cdofb_monolithic_sles_t msles 
)

Solve a linear system arising from the discretization of the Navier-Stokes equation with a CDO face-based approach. The system is split into blocks to enable more efficient preconditioning techniques. The main iterative solver is a Krylov solver such as GCR, or MINRES.

Parameters
[in]nsppointer to a cs_navsto_param_t structure
[in]eqppointer to a cs_equation_param_t structure
[in]shpointer to a cs_cdo_system_helper_t structure
[in,out]slesppointer to a set of parameters to drive the SLES
[in,out]mslespointer to a cs_cdofb_monolithic_sles_t structure
Returns
the (cumulated) number of iterations of the solver

◆ cs_cdofb_monolithic_gkb_solve()

int cs_cdofb_monolithic_gkb_solve ( const cs_navsto_param_t nsp,
const cs_equation_param_t eqp,
const cs_cdo_system_helper_t sh,
cs_param_sles_t slesp,
cs_cdofb_monolithic_sles_t msles 
)

Use the GKB algorithm to solve the saddle-point problem arising from CDO-Fb schemes for Stokes and Navier-Stokes with a monolithic coupling.

Parameters
[in]nsppointer to a cs_navsto_param_t structure
[in]eqppointer to a cs_equation_param_t structure
[in]shpointer to a cs_cdo_system_helper_t structure
[in,out]slesppointer to a set of parameters to drive the SLES
[in,out]mslespointer to a cs_cdofb_monolithic_sles_t structure
Returns
the cumulated number of iterations of the solver

◆ cs_cdofb_monolithic_set_sles()

void cs_cdofb_monolithic_set_sles ( cs_navsto_param_t nsp,
void *  context 
)

Start setting-up the Navier-Stokes equations when a monolithic algorithm is used to couple the system. No mesh information is available at this stage. nsp is not declared as const to avoid compilation warnings but it should be modified at this stage.

Parameters
[in]nsppointer to a cs_navsto_param_t structure
[in,out]contextpointer to a context structure cast on-the-fly

◆ cs_cdofb_monolithic_sles_clean()

void cs_cdofb_monolithic_sles_clean ( cs_cdofb_monolithic_sles_t msles)

Free a part of the structure.

Parameters
[in,out]mslespointer to the structure to clean

◆ cs_cdofb_monolithic_sles_create()

cs_cdofb_monolithic_sles_t* cs_cdofb_monolithic_sles_create ( cs_lnum_t  n_faces,
cs_lnum_t  n_cells 
)

Create an empty cs_cdofb_monolithic_sles_t structure.

Parameters
[in]n_facesnumber of faces (interior + border)
[in]n_cellsnumber of cells
Returns
a pointer to a newly allocated structure

◆ cs_cdofb_monolithic_sles_finalize()

void cs_cdofb_monolithic_sles_finalize ( void  )

Free if needed structure(s) associated CDO face-based schemes with a monolithic velocity-pressure coupling.

◆ cs_cdofb_monolithic_sles_free()

void cs_cdofb_monolithic_sles_free ( cs_cdofb_monolithic_sles_t **  p_msles)

Free memory related to cs_cdofb_monolithic_sles_t structure.

Parameters
[in,out]p_mslesdouble pointer to the structure to free

◆ cs_cdofb_monolithic_sles_init_sharing()

void cs_cdofb_monolithic_sles_init_sharing ( const cs_mesh_t mesh,
const cs_cdo_connect_t connect,
const cs_cdo_quantities_t quant 
)

Set pointers to shared structures.

Parameters
[in]meshpointer to the mesh structure
[in]connectpointer to additional CDO connectivities
[in]quantpointer to additional CDO mesh quantities

◆ cs_cdofb_monolithic_solve()

int cs_cdofb_monolithic_solve ( const cs_navsto_param_t nsp,
const cs_equation_param_t eqp,
const cs_cdo_system_helper_t sh,
cs_param_sles_t slesp,
cs_cdofb_monolithic_sles_t msles 
)

Solve a linear system arising from the discretization of the Navier-Stokes equation with a CDO face-based approach. The full system is treated as one block and solved as it is. In this situation, PETSc or MUMPS are usually considered.

Parameters
[in]nsppointer to a cs_navsto_param_t structure
[in]eqppointer to a cs_equation_param_t structure
[in]shpointer to a cs_cdo_system_helper_t structure
[in,out]slesppointer to a set of parameters to drive the SLES
[in,out]mslespointer to a cs_cdofb_monolithic_sles_t structure
Returns
the (cumulated) number of iterations of the solver

◆ cs_cdofb_monolithic_uzawa_al_incr_solve()

int cs_cdofb_monolithic_uzawa_al_incr_solve ( const cs_navsto_param_t nsp,
const cs_equation_param_t eqp,
const cs_cdo_system_helper_t sh,
cs_param_sles_t slesp,
cs_cdofb_monolithic_sles_t msles 
)

Use the Uzawa algorithm with an Augmented Lagrangian (ALU) technique in an incremental way to solve the saddle-point problem arising from CDO-Fb schemes for Stokes, Oseen and Navier-Stokes with a monolithic coupling.

Parameters
[in]nsppointer to a cs_navsto_param_t structure
[in]eqppointer to a cs_equation_param_t structure
[in]shpointer to a cs_cdo_system_helper_t structure
[in,out]slesppointer to a set of parameters to drive the SLES
[in,out]mslespointer to a cs_cdofb_monolithic_sles_t structure
Returns
the cumulated number of iterations of the solver

◆ cs_cdofb_monolithic_uzawa_cg_solve()

int cs_cdofb_monolithic_uzawa_cg_solve ( const cs_navsto_param_t nsp,
const cs_equation_param_t eqp,
const cs_cdo_system_helper_t sh,
cs_param_sles_t slesp,
cs_cdofb_monolithic_sles_t msles 
)

Use the preconditioned Uzawa-CG algorithm to solve the saddle-point problem arising from CDO-Fb schemes for Stokes, Oseen and Navier-Stokes with a monolithic coupling. This algorithm is based on Koko's paper "Uzawa conjugate gradient method for the Stokes problem: Matlab implementation with P1-iso-P2/ P1 finite element".

Parameters
[in]nsppointer to a cs_navsto_param_t structure
[in]eqppointer to a cs_equation_param_t structure
[in]shpointer to a cs_cdo_system_helper_t structure
[in,out]slesppointer to a set of parameters to drive the SLES
[in,out]mslespointer to a cs_cdofb_monolithic_sles_t structure
Returns
the cumulated number of iterations of the solver

Use the preconditioned Uzawa-CG algorithm to solve the saddle-point problem arising from CDO-Fb schemes for Stokes, Oseen and Navier-Stokes with a monolithic coupling. This algorithm is based on Koko's paper "Uzawa conjugate gradient method for the Stokes problem: Matlab implementation with P1-iso-P2/ P1 finite element".

Parameters
[in]nsppointer to a cs_navsto_param_t structure
[in]eqppointer to a cs_equation_param_t structure
[in]shpointer to a cs_cdo_system_helper_t structure
[in,out]slesppointer to a set of parameters to drive the SLES
[in,out]mslespointer to a cs_cdofb_monolithic_sles_t structure
Returns
the cumulated number of iterations of the solver