1# ifndef _RHEO_MIXED_SOLVER_H
2# define _RHEO_MIXED_SOLVER_H
23#include "rheolef/uzawa.h"
24#include "rheolef/cg.h"
25#include "rheolef/minres.h"
28template <
class Matrix,
class Vector,
class Solver>
38template <
class Matrix,
class Vector,
class Solver>
48template <
class Matrix,
class Vector,
class Solver,
class Preconditioner>
50 const Vector& Mf,
const Vector& Mg,
const Preconditioner& S1,
51 const Solver& inner_solver_A,
const Float& rho,
54 sopt.label = (sopt.label !=
"" ? sopt.label :
"uzawa_abtbc");
56 Vector v = inner_solver_A.solve (Mf);
58 int status =
uzawa (S,
p, Mh, S1, sopt);
59 u = inner_solver_A.solve (Mf - B.trans_mult(
p));
62template <
class Matrix,
class Vector,
class Solver,
class Preconditioner,
class Real>
64 const Vector& Mf,
const Vector& Mg,
const Preconditioner& S1,
65 const Solver& inner_solver_A,
const Real& rho,
68 sopt.label = (sopt.label !=
"" ? sopt.label :
"uzawa_abtb");
70 Vector v = inner_solver_A.solve (Mf);
72 int status =
uzawa (S,
p, Mh, S1, rho, sopt);
73 u = inner_solver_A.solve (Mf - B.trans_mult(
p));
150template <
class Matrix,
class Vector,
class VectorExpr1,
class VectorExpr2,
151 class Solver,
class Preconditioner>
153 const VectorExpr1& Mf,
const VectorExpr2& Mg,
const Preconditioner& S1,
156 sopt.label = (sopt.label !=
"" ? sopt.label :
"cg_abtbc");
158 Vector v = inner_solver_A.solve (Mf);
160 int status =
cg (S,
p, Mh, S1, sopt);
161 u = inner_solver_A.solve (Mf - B.trans_mult(
p));
164template <
class Matrix,
class Vector,
class VectorExpr1,
class VectorExpr2,
165 class Solver,
class Preconditioner>
167 const VectorExpr1& Mf,
const VectorExpr2& Mg,
const Preconditioner& S1,
170 sopt.label = (sopt.label !=
"" ? sopt.label :
"cg_abtb");
172 Vector v = inner_solver_A.solve (Mf);
174 int status =
cg (S,
p, Mh, S1, sopt);
175 u = inner_solver_A.solve (Mf - B.trans_mult(
p));
178template <
class Matrix,
class Vector,
class Solver,
class Preconditioner>
180 const Vector& Mf,
const Vector& Mg,
const Preconditioner& S1,
183 sopt.label = (sopt.label !=
"" ? sopt.label :
"minres_abtbc");
185 Vector v = inner_solver_A.solve (Mf);
187 int status =
minres(S,
p, Mh, S1, sopt);
188 u = inner_solver_A.solve (Mf - B.trans_mult(
p));
191template <
class Matrix,
class Vector,
class Solver,
class Preconditioner>
193 const Vector& Mf,
const Vector& Mg,
const Preconditioner& S1,
196 sopt.label = (sopt.label !=
"" ? sopt.label :
"minres_abtb");
198 Vector v = inner_solver_A.solve (Mf);
200 int status =
minres(S,
p, Mh, S1, sopt);
201 u = inner_solver_A.solve (Mf - B.trans_mult(
p));
see the Float page for the full documentation
see the solver_option page for the full documentation
This file is part of Rheolef.
int cg_abtbc(const Matrix &A, const Matrix &B, const Matrix &C, Vector &u, Vector &p, const VectorExpr1 &Mf, const VectorExpr2 &Mg, const Preconditioner &S1, const Solver &inner_solver_A, const solver_option &sopt=solver_option())
int cg_abtb(const Matrix &A, const Matrix &B, Vector &u, Vector &p, const VectorExpr1 &Mf, const VectorExpr2 &Mg, const Preconditioner &S1, const Solver &inner_solver_A, const solver_option &sopt=solver_option())
int minres(const Matrix &A, Vector &x, const Vector &Mb, const Preconditioner &M, const solver_option &sopt=solver_option())
int uzawa_abtb(const Matrix &A, const Matrix &B, Vector &u, Vector &p, const Vector &Mf, const Vector &Mg, const Preconditioner &S1, const Solver &inner_solver_A, const Real &rho, const solver_option &sopt=solver_option())
int minres_abtbc(const Matrix &A, const Matrix &B, const Matrix &C, Vector &u, Vector &p, const Vector &Mf, const Vector &Mg, const Preconditioner &S1, const Solver &inner_solver_A, const solver_option &sopt=solver_option())
int minres_abtb(const Matrix &A, const Matrix &B, Vector &u, Vector &p, const Vector &Mf, const Vector &Mg, const Preconditioner &S1, const Solver &inner_solver_A, const solver_option &sopt=solver_option())
int cg(const Matrix &A, Vector &x, const Vector2 &Mb, const Preconditioner &M, const solver_option &sopt=solver_option())
int uzawa(const Matrix &A, Vector &x, const Vector &Mb, const Preconditioner &M, const Real2 &rho, const solver_option &sopt=solver_option())
int uzawa_abtbc(const Matrix &A, const Matrix &B, const Matrix &C, Vector &u, Vector &p, const Vector &Mf, const Vector &Mg, const Preconditioner &S1, const Solver &inner_solver_A, const Float &rho, const solver_option &sopt=solver_option())
Vector operator*(const Vector &p) const
abtb_schur_complement(const Solver &iA1, const Matrix &B1)
Vector operator*(const Vector &p) const
abtbc_schur_complement(const Solver &iA1, const Matrix &B1, const Matrix &C1)