escript Revision_
Preconditioner.h
Go to the documentation of this file.
1
2/*****************************************************************************
3*
4* Copyright (c) 2003-2020 by The University of Queensland
5* http://www.uq.edu.au
6*
7* Primary Business: Queensland, Australia
8* Licensed under the Apache License, version 2.0
9* http://www.apache.org/licenses/LICENSE-2.0
10*
11* Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12* Development 2012-2013 by School of Earth Sciences
13* Development from 2014-2017 by Centre for Geoscience Computing (GeoComp)
14* Development from 2019 by School of Earth and Environmental Sciences
15**
16*****************************************************************************/
17
18#ifndef __PASO_PRECONDITIONER_H__
19#define __PASO_PRECONDITIONER_H__
20
21#include "Paso.h"
22#include "SystemMatrix.h"
23
24namespace paso {
25
26struct MergedSolver;
27struct Preconditioner;
28typedef boost::shared_ptr<Preconditioner> Preconditioner_ptr;
29typedef boost::shared_ptr<const Preconditioner> const_Preconditioner_ptr;
30
32struct Solver_ILU;
33struct Solver_RILU;
34
35// general preconditioner interface
49
52void Preconditioner_solve(Preconditioner* prec, SystemMatrix_ptr<double> A, double*, double*);
53
54
55// GAUSS SEIDEL & Jacobi
57{
58 bool Jacobi;
59 double* diag;
60 double* buffer;
62};
63
69
72
74 SystemMatrix_ptr<double> A, bool jacobi, bool is_local, bool verbose);
75
77 SparseMatrix_ptr<double> A, bool jacobi, bool verbose);
78
80 Preconditioner_Smoother* gs, double* x, const double* b,
81 dim_t sweeps, bool x_is_initial);
82
84 Preconditioner_LocalSmoother* gs, double* x, const double* b,
85 dim_t sweeps, bool x_is_initial);
86
88 Preconditioner_Smoother* gs, double* x, const double* b,
89 double atol, dim_t* sweeps, bool x_is_initial);
90
92 Preconditioner_LocalSmoother* gs, double* x);
93
96 double* x);
97
99 Preconditioner_LocalSmoother* gs, double* x);
100
102 Preconditioner_LocalSmoother* gs, double* x);
103
106{
107 double* factors;
108};
109
131
132void Solver_ILU_free(Solver_ILU * in);
134void Solver_solveILU(SparseMatrix_ptr<double> A, Solver_ILU* ilu, double* x, const double* b);
135
138void Solver_solveRILU(Solver_RILU* rilu, double* x, double* b);
139
141 SparseMatrix_ptr<double> A_CF, double* invA_FF, index_t* A_FF_pivot,
143
144} // namespace paso
145
146#endif // __PASO_PRECONDITIONER_H__
147
index_t dim_t
Definition DataTypes.h:66
int index_t
type for array/matrix indices used both globally and on each rank
Definition DataTypes.h:61
Definition BiCGStab.cpp:25
void Preconditioner_Smoother_solve(SystemMatrix_ptr< double > A, Preconditioner_Smoother *gs, double *x, const double *b, dim_t sweeps, bool x_is_initial)
Definition Smoother.cpp:100
void Preconditioner_solve(Preconditioner *prec, SystemMatrix_ptr< double > A, double *x, double *b)
Definition Preconditioner.cpp:111
void Preconditioner_LocalSmoother_Sweep(SparseMatrix_ptr< double > A, Preconditioner_LocalSmoother *gs, double *x)
Definition Smoother.cpp:201
Preconditioner * Preconditioner_alloc(SystemMatrix_ptr< double > A, Options *options)
Definition Preconditioner.cpp:46
boost::shared_ptr< SystemMatrix< T > > SystemMatrix_ptr
Definition SystemMatrix.h:42
void Preconditioner_free(Preconditioner *in)
Definition Preconditioner.cpp:35
void Solver_ILU_free(Solver_ILU *in)
Definition ILU.cpp:35
void Preconditioner_LocalSmoother_Sweep_sequential(SparseMatrix_ptr< double > A, Preconditioner_LocalSmoother *gs, double *x)
inplace Gauss-Seidel sweep in sequential mode
Definition Smoother.cpp:221
void Preconditioner_LocalSmoother_Sweep_tiled(SparseMatrix_ptr< double > A, Preconditioner_LocalSmoother *gs, double *x)
void Preconditioner_Smoother_free(Preconditioner_Smoother *in)
Definition Smoother.cpp:34
Preconditioner_Smoother * Preconditioner_Smoother_alloc(SystemMatrix_ptr< double > A, bool jacobi, bool is_local, bool verbose)
constructs the symmetric Gauss-Seidel preconditioner
Definition Smoother.cpp:54
boost::shared_ptr< Preconditioner > Preconditioner_ptr
Definition Preconditioner.h:28
boost::shared_ptr< const Preconditioner > const_Preconditioner_ptr
Definition Preconditioner.h:29
SolverResult
Definition Paso.h:44
void Solver_solveILU(SparseMatrix_ptr< double > A, Solver_ILU *ilu, double *x, const double *b)
Definition ILU.cpp:316
void Preconditioner_LocalSmoother_Sweep_colored(SparseMatrix_ptr< double > A, Preconditioner_LocalSmoother *gs, double *x)
Definition Smoother.cpp:338
void Solver_solveRILU(Solver_RILU *rilu, double *x, double *b)
Definition RILU.cpp:289
void Preconditioner_LocalSmoother_free(Preconditioner_LocalSmoother *in)
Definition Smoother.cpp:42
Solver_RILU * Solver_getRILU(SparseMatrix_ptr< double > A, bool verbose)
Definition RILU.cpp:71
boost::shared_ptr< SparseMatrix< T > > SparseMatrix_ptr
Definition SparseMatrix.h:37
Solver_ILU * Solver_getILU(SparseMatrix_ptr< double > A, bool verbose)
constructs the incomplete block factorization
Definition ILU.cpp:44
void Solver_RILU_free(Solver_RILU *in)
Definition RILU.cpp:35
Preconditioner_LocalSmoother * Preconditioner_LocalSmoother_alloc(SparseMatrix_ptr< double > A, bool jacobi, bool verbose)
Definition Smoother.cpp:64
void Solver_updateIncompleteSchurComplement(SparseMatrix_ptr< double > A_CC, SparseMatrix_ptr< double > A_CF, double *invA_FF, index_t *A_FF_pivot, SparseMatrix_ptr< double > A_FC)
Definition SchurComplement.cpp:28
void Preconditioner_LocalSmoother_solve(SparseMatrix_ptr< double > A, Preconditioner_LocalSmoother *gs, double *x, const double *b, dim_t sweeps, bool x_is_initial)
Definition Smoother.cpp:161
SolverResult Preconditioner_Smoother_solve_byTolerance(SystemMatrix_ptr< double > A, Preconditioner_Smoother *gs, double *x, const double *b, double atol, dim_t *sweeps, bool x_is_initial)
Definition Smoother.cpp:127
#define PASO_DLL_API
Definition paso/src/system_dep.h:29
Definition Options.h:80
Definition Preconditioner.h:57
double * buffer
Definition Preconditioner.h:60
bool Jacobi
Definition Preconditioner.h:58
index_t * pivot
Definition Preconditioner.h:61
double * diag
Definition Preconditioner.h:59
Definition Preconditioner.h:65
bool is_local
Definition Preconditioner.h:67
Preconditioner_LocalSmoother * localSmoother
Definition Preconditioner.h:66
Definition Preconditioner.h:37
dim_t type
Definition Preconditioner.h:38
Solver_RILU * rilu
RILU preconditioner.
Definition Preconditioner.h:47
Preconditioner_Smoother * jacobi
Jacobi preconditioner.
Definition Preconditioner.h:41
Solver_ILU * ilu
ILU preconditioner.
Definition Preconditioner.h:45
dim_t sweeps
Definition Preconditioner.h:39
Preconditioner_Smoother * gs
Gauss-Seidel preconditioner.
Definition Preconditioner.h:43
ILU preconditioner.
Definition Preconditioner.h:106
double * factors
Definition Preconditioner.h:107
RILU preconditioner.
Definition Preconditioner.h:112
double * b_F
Definition Preconditioner.h:126
index_t * mask_F
Definition Preconditioner.h:123
SparseMatrix_ptr< double > A_FC
Definition Preconditioner.h:119
index_t * A_FF_pivot
Definition Preconditioner.h:118
double * inv_A_FF
Definition Preconditioner.h:117
double * x_C
Definition Preconditioner.h:127
index_t * rows_in_C
Definition Preconditioner.h:122
dim_t n_block
Definition Preconditioner.h:114
SparseMatrix_ptr< double > A_CF
Definition Preconditioner.h:120
double * b_C
Definition Preconditioner.h:128
dim_t n
Definition Preconditioner.h:113
dim_t n_C
Definition Preconditioner.h:116
index_t * rows_in_F
Definition Preconditioner.h:121
index_t * mask_C
Definition Preconditioner.h:124
Solver_RILU * RILU_of_Schur
Definition Preconditioner.h:129
dim_t n_F
Definition Preconditioner.h:115
double * x_F
Definition Preconditioner.h:125