8.0
general documentation
Loading...
Searching...
No Matches
cs_saddle_itsol.h
Go to the documentation of this file.
1#ifndef __CS_SADDLE_ITSOL_H__
2#define __CS_SADDLE_ITSOL_H__
3
4/*============================================================================
5 * In-house iterative solvers defined by blocks and associated to CDO
6 * discretizations for saddle-point system
7 *============================================================================*/
8
9/*
10 This file is part of code_saturne, a general-purpose CFD tool.
11
12 Copyright (C) 1998-2023 EDF S.A.
13
14 This program is free software; you can redistribute it and/or modify it under
15 the terms of the GNU General Public License as published by the Free Software
16 Foundation; either version 2 of the License, or (at your option) any later
17 version.
18
19 This program is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
21 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
22 details.
23
24 You should have received a copy of the GNU General Public License along with
25 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
26 Street, Fifth Floor, Boston, MA 02110-1301, USA.
27*/
28
29/*----------------------------------------------------------------------------*/
30
31#include "cs_defs.h"
32
33/*----------------------------------------------------------------------------
34 * Standard C library headers
35 *----------------------------------------------------------------------------*/
36
37/*----------------------------------------------------------------------------
38 * Local headers
39 *----------------------------------------------------------------------------*/
40
41#include "cs_iter_algo.h"
42#include "cs_matrix.h"
43#include "cs_mesh_adjacencies.h"
44#include "cs_param_sles.h"
45#include "cs_range_set.h"
46#include "cs_sles.h"
47
48/*----------------------------------------------------------------------------*/
49
51
52/*============================================================================
53 * Macro definitions
54 *============================================================================*/
55
56/*============================================================================
57 * Type definitions
58 *============================================================================*/
59
60typedef struct {
61
62 /*
63 * A hybrid storage is used to represent the system described below
64 *
65 * | M11 | M12 | | x1 | |rhs1 |
66 * M = |-----------| |----| = |-----|
67 * | M21 | 0 | | x2 | |rhs2 |
68 *
69 * One assumes that M12 = M21^T
70 */
71
72 int n_m11_matrices; /* 1, 3, 6 or 9 is expected */
74 cs_lnum_t x1_size; /* scatter view */
75 cs_lnum_t max_x1_size; /* max(x1_size, n_m11_cols) */
76
78
79 cs_lnum_t x2_size; /* scatter view */
80 cs_lnum_t max_x2_size; /* max(x2_size, n_m22_cols) */
82
85
86 /* Indexed list used to scan the unassembled m21 operator */
87
89
90 /* Structure used for synchronisation (parallel or periodic). Enable to
91 switch from a scatter view (the mesh view) to a gather view (the algebraic
92 view). This structure is shared. */
93
95
97
98/* Structure handling the block preconditioning */
99
100typedef struct {
101
102 /* Parameters */
103
106
107 /* Block 11 settings (mandatory) */
108
111
112 /* Schur complement settings (optional) */
113
118
119 /* Native arrays for the Schur matrix (optional) */
120
123
124 /* Diagonal approximations of block matrices (optional) */
125
128
130
131/*============================================================================
132 * Public function prototypes
133 *============================================================================*/
134
135/*----------------------------------------------------------------------------*/
146/*----------------------------------------------------------------------------*/
147
150 cs_param_schur_approx_t schur_type,
151 cs_param_sles_t *m11_slesp,
152 cs_sles_t *m11_sles);
153
154/*----------------------------------------------------------------------------*/
160/*----------------------------------------------------------------------------*/
161
162void
164
165/*----------------------------------------------------------------------------*/
179/*----------------------------------------------------------------------------*/
180
181void
184 cs_real_t *x1,
185 cs_real_t *x2,
186 cs_iter_algo_t *algo);
187
188/*----------------------------------------------------------------------------*/
205/*----------------------------------------------------------------------------*/
206
207void
208cs_saddle_gcr(int restart,
209 cs_saddle_system_t *ssys,
211 cs_real_t *x1,
212 cs_real_t *x2,
213 cs_iter_algo_t *algo);
214
215/*----------------------------------------------------------------------------*/
231/*----------------------------------------------------------------------------*/
232
233void
235 const cs_matrix_t *mat,
236 cs_real_t *vec,
237 cs_real_t *matvec);
238
239/*----------------------------------------------------------------------------*/
256/*----------------------------------------------------------------------------*/
257
258void
260 const cs_matrix_t *mat,
261 cs_lnum_t vec_len,
262 cs_real_t *vec,
263 cs_real_t **p_matvec);
264
265/*----------------------------------------------------------------------------*/
266
268
269#endif /* __CS_SADDLE_ITSOL_H__ */
#define BEGIN_C_DECLS
Definition cs_defs.h:509
double cs_real_t
Floating-point value.
Definition cs_defs.h:319
#define END_C_DECLS
Definition cs_defs.h:510
int cs_lnum_t
local mesh entity id
Definition cs_defs.h:313
@ x2
Definition cs_field_pointer.h:223
struct _cs_matrix_t cs_matrix_t
Definition cs_matrix.h:110
Structure and routines handling the SLES settings stored inside a cs_param_sles_t structure.
cs_param_precond_block_t
Definition cs_param_types.h:818
cs_param_schur_approx_t
Strategy to build the Schur complement approximation. This appears in block preconditioning or Uzawa ...
Definition cs_param_types.h:752
void cs_saddle_block_precond_free(cs_saddle_block_precond_t **p_sbp)
Free a cs_saddle_block_precond_t structure.
Definition cs_saddle_itsol.c:1423
void cs_saddle_gcr(int restart, cs_saddle_system_t *ssys, cs_saddle_block_precond_t *sbp, cs_real_t *x1, cs_real_t *x2, cs_iter_algo_t *algo)
Apply the GCR algorithm to a saddle point problem (the system is stored in a hybrid way)....
Definition cs_saddle_itsol.c:1672
void cs_matrix_vector_multiply_gs(const cs_range_set_t *rset, const cs_matrix_t *mat, cs_lnum_t vec_len, cs_real_t *vec, cs_real_t **p_matvec)
Perform a matrix-vector multiplication in case of scatter-view array as input parameter....
Definition cs_saddle_itsol.c:1922
void cs_saddle_minres(cs_saddle_system_t *ssys, cs_saddle_block_precond_t *sbp, cs_real_t *x1, cs_real_t *x2, cs_iter_algo_t *algo)
Apply the MINRES algorithm to a saddle point problem (the system is stored in a hybrid way)....
Definition cs_saddle_itsol.c:1459
cs_saddle_block_precond_t * cs_saddle_block_precond_create(cs_param_precond_block_t block_type, cs_param_schur_approx_t schur_type, cs_param_sles_t *m11_slesp, cs_sles_t *m11_sles)
Create and initialize a cs_saddle_block_precond_t structure.
Definition cs_saddle_itsol.c:1376
void cs_matrix_vector_multiply_gs_allocated(const cs_range_set_t *rset, const cs_matrix_t *mat, cs_real_t *vec, cs_real_t *matvec)
Perform a matrix-vector multiplication in case of scatter-view array as input parameter....
Definition cs_saddle_itsol.c:1867
struct _cs_sles_t cs_sles_t
Definition cs_sles.h:68
Definition cs_mesh_adjacencies.h:68
Structure to handle the convergence of an iterative algorithm.
Definition cs_iter_algo.h:60
Structure storing all metadata related to the resolution of a linear system with an iterative solver.
Definition cs_param_sles.h:91
Definition cs_range_set.h:57
Definition cs_saddle_itsol.h:100
cs_real_t * schur_diag
Definition cs_saddle_itsol.h:121
cs_real_t * m11_inv_diag
Definition cs_saddle_itsol.h:126
cs_param_schur_approx_t schur_type
Definition cs_saddle_itsol.h:105
cs_sles_t * m11_sles
Definition cs_saddle_itsol.h:110
cs_param_sles_t * m11_slesp
Definition cs_saddle_itsol.h:109
cs_param_precond_block_t block_type
Definition cs_saddle_itsol.h:104
cs_matrix_t * schur_matrix
Definition cs_saddle_itsol.h:114
cs_param_sles_t * schur_slesp
Definition cs_saddle_itsol.h:115
double schur_scaling
Definition cs_saddle_itsol.h:117
cs_sles_t * schur_sles
Definition cs_saddle_itsol.h:116
cs_real_t * mass22_diag
Definition cs_saddle_itsol.h:127
cs_real_t * schur_xtra
Definition cs_saddle_itsol.h:122
Definition cs_saddle_itsol.h:60
cs_lnum_t max_x1_size
Definition cs_saddle_itsol.h:75
cs_real_t * m21_unassembled
Definition cs_saddle_itsol.h:83
cs_real_t * rhs1
Definition cs_saddle_itsol.h:77
cs_real_t * rhs2
Definition cs_saddle_itsol.h:81
const cs_range_set_t * rset
Definition cs_saddle_itsol.h:94
int m21_stride
Definition cs_saddle_itsol.h:84
cs_matrix_t ** m11_matrices
Definition cs_saddle_itsol.h:73
cs_lnum_t x1_size
Definition cs_saddle_itsol.h:74
const cs_adjacency_t * m21_adjacency
Definition cs_saddle_itsol.h:88
cs_lnum_t x2_size
Definition cs_saddle_itsol.h:79
int n_m11_matrices
Definition cs_saddle_itsol.h:72
cs_lnum_t max_x2_size
Definition cs_saddle_itsol.h:80