KASKADE 7 development version
|
A class representing small KKT systems on a single subdomain. More...
#include <bddc.hh>
A class representing small KKT systems on a single subdomain.
This is intended to work together with BDDCSolver and KKTSolver, and represents subdomains, providing subdomain-local storage and processing capabilities. This is the basis for distributed domain decomposition solvers.
A Subdomain object represents the system
B = \begin{bmatrix} A & C^T \\ C \end{bmatrix},
and provides in particular access to the "Schur complement"
S = \begin{bmatrix} 0 & I \end{bmatrix} \begin{bmatrix} A & C^T \\ C \end{bmatrix}^{-1} \begin{bmatrix} 0 \\ I \end{bmatrix}.
Since A is only positive semidefinite in general, the true Schur complement - C A^{-1} C^T cannot be formed, but is equal in case of invertible A .
The subdomain holds A_i , f_i , and the constraints block C . It also holds a local solution vector u_i . For the overall picture, we refer to KKTSolver.
m | the number of vectorial components in x (usually 1 for scalar systems, and 2 or 3 for vectorial ones) |
Scalar | the field type K , usually double |
CoarseScalar | the field type K_c to be used for representing the coarse solver factorization, usually double , but could also be float |
Transfer | a space transfer policy class |
float
factorization. Extract L and U factors and convert to float?This is a reference implementation.
Public Types | |
using | Scalar = Scalar_ |
The scalar field type used in coefficients. More... | |
using | XEntry = Dune::FieldVector< Scalar, components > |
using | XVector = Dune::BlockVector< XEntry > |
using | Vector = Dune::BlockVector< Dune::FieldVector< Scalar, 1 > > |
using | TransmissionType = typename Transfer::TransmissionType |
A container type used for exchanging coefficients between subdomains. More... | |
using | AMatrix = NumaBCRSMatrix< Dune::FieldMatrix< Scalar, components, components > > |
using | SMatrix = DynamicMatrix< Dune::FieldMatrix< Scalar, 1, 1 > > |
A matrix type with scalar entries, used for the Schur complement matrix. More... | |
Public Member Functions | |
std::vector< int > | neighbors () const |
Reports a list of its neighbor subdomains. More... | |
void | updateSolution (Scalar alpha) |
Updates the solution approximation. More... | |
void | getSolution (XVector &x) const |
Provides the current solution vector x . More... | |
XVector | getSolution () const |
Provides the current solution vector x . More... | |
void | set_solver_type (bool cg_solver) |
Set solver type. More... | |
void | updateIteration (int cg_itr) |
void | getCorrection (XVector &dx) const |
Provides the current correction vector dx . More... | |
void | getRawCorrection (XVector &dx) const |
void | getResidual (XVector &rr) const |
Vector const & | getRestrictedResidual () const |
std::vector< bool > | getActiveSubdomains () |
void | setActiveSubdomains_sub (std::vector< bool > &activeSubdomains_) |
Construction. | |
template<int cComponents> | |
Subdomain (int id, AMatrix const &A, Interfaces const &interfaces, NumaBCRSMatrix< Dune::FieldMatrix< Scalar, cComponents, components > > const &C) | |
Constructor. More... | |
template<class Interface > | |
Subdomain (int id, AMatrix const &A, Interface const &ifs) | |
void | imposeDC_inActiveInterface (Vector &Aq, Vector &q) |
void | imposeDCq (Vector &q) |
void | update_rhs (XVector const &f) |
update rhs for the current subdomain More... | |
void | getInterfaceAveragingWeights (int s, XVector &rho) const |
Reports own interface averaging weights. More... | |
void | setInterfaceAveragingWeights (int s, XVector const &rho) |
Accumulates interface weights of adjacent subdomains. More... | |
Restriction support. | |
These methods implement the "restriction" of the residual to the "coarse" space. The restriction is defined as the transpose of the prolongation, i.e. \Pi^T = \pi^T E^T with averaging \pi and extension E . The extension operator E^T acts locally, realizing the dual harmonic extension \begin{bmatrix} r_i \\ r_D \\ \bar r_D \end{bmatrix}_{\rm new} = \begin{bmatrix} I & \\ H^T & \\ -H^T & I \end{bmatrix} \begin{bmatrix} r_i \\ r_D \end{bmatrix} working on interior residual r_i , interface residual r_D , and averaged interface residual \bar r_D . Here, H = A_{ii}^{-1} A_{iD} denotes the harmonic extension of interface values into the interior. This local operation is performed by the preRestrict() method. The averaging (for a single degree of freedom shared by k subdomains) is \pi^T = \frac{1}{k} \begin{bmatrix} 1 & \dots & 1 \\ \vdots & \ddots & \vdots \\ 1 & \dots & 1 \end{bmatrix}, therefore symmetric, and acts on the interface residual values r_D only. It involves communication between subdomains, exchanging their interface residual r_D , which is realized by the getInterfaceResidual() and setInterfaceResidual() methods. Finally, restrict() must be called as a sentinel to denote the end of communication, and to perform the actual averaging. | |
void | preRestrict () |
First step of the restriction: Dirichlet solve. More... | |
void | getInterfaceResidualUpdate (int s, TransmissionType &dres) const |
Gets interface residual that shall be shared with neighboring subdomain. More... | |
void | setInterfaceResidualUpdate (int s, TransmissionType const &dres) |
Sets interface residual of the correction from neighboring subdomains. More... | |
void | restrict () |
Performs the restriction of the residual vector. More... | |
Coarse solve support. | |
These methods implement the subdomain-local parts of the coarse space solver. | |
void | getS (SMatrix &S) const |
Compute S . More... | |
void | getSchurResidual (Vector &lambda) const |
Provides the Schur residual. More... | |
void | setCorrection (Vector const &c) |
Computes the correction vector [\delta x, \delta\lambda] . More... | |
Prolongation support. | |
These methods implement the "prolongation" of "coarse" space corrections to the "fine" space. The prolongation is defined as \Pi = E\pi with averaging \pi and extension E . The averaging (for a single degree of freedom shared by k subdomains) is \pi = \frac{1}{k} \begin{bmatrix} 1 & \dots & 1 \\ \vdots & \ddots & \vdots \\ 1 & \dots & 1 \end{bmatrix} and acts on the interface degrees of freedom x_D only. It involves communication between subdomains, exchanging their interface values x_D , which is realized by the getInterfaceValues() and setInterfaceValues() methods. The extension operator E acts locally, realizing the harmonic extension \begin{bmatrix} x_i \\ x_D \end{bmatrix}_{\rm new} = \begin{bmatrix} I & H & - H \\ & & I \end{bmatrix} \begin{bmatrix} x_i \\ x_D \\ \bar x_D \end{bmatrix} working on interior dofs x_i , interface dofs x_D , and averaged interface values \bar x_D . Here, H = A_{ii}^{-1} A_{iD} denotes the harmonic extension of interface values into the interior. This local operation is performed by the prolongate() method. | |
void | preProlongate () |
Prepares the exchange of interface values. More... | |
void | getInterfaceValues (int s, TransmissionType &values) const |
Gets boundary values that shall be shared with neighboring subdomain. More... | |
void | setInterfaceValues (int s, TransmissionType const &values) |
Sets boundary values of the correction from neighboring subdomains. More... | |
void | setActiveId (bool active) |
Dune::FieldVector< double, 2 > | prolongate () |
Perform averaging on the shared interfaces, projecting the solution to the subspace of functions continuous across the interfaces. More... | |
Dune::FieldVector< double, 1 > | updateSearchDirectionNorm () |
void | updateSearchDirection (double beta) |
Update search direction q for cg solver. More... | |
Static Public Attributes | |
static const int | components = m |
Number of (vectorial) components of u . More... | |
using Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::AMatrix = NumaBCRSMatrix<Dune::FieldMatrix<Scalar,components,components> > |
using Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::Scalar = Scalar_ |
using Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::SMatrix = DynamicMatrix<Dune::FieldMatrix<Scalar,1,1> > |
using Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::TransmissionType = typename Transfer::TransmissionType |
using Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::Vector = Dune::BlockVector<Dune::FieldVector<Scalar,1> > |
using Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::XEntry = Dune::FieldVector<Scalar,components> |
using Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::XVector = Dune::BlockVector<XEntry> |
Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::Subdomain | ( | int | id, |
AMatrix const & | A, | ||
Interfaces const & | interfaces, | ||
NumaBCRSMatrix< Dune::FieldMatrix< Scalar, cComponents, components > > const & | C | ||
) |
Constructor.
This initializes the solution approximation to zero, which implies an admissible state over all subdomains.
id | the subdomain number |
A | the subdomain-local stiffness matrix |
f | the subdomain-local right hand side |
interfaces | a list of pairs (s,ids) of neighboring subdomain number and the degrees of freedom shared with this subdomain (in the given order) |
C | the constraint matrix |
Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::Subdomain | ( | int | id, |
AMatrix const & | A, | ||
Interface const & | ifs | ||
) |
|
inline |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::getCorrection | ( | XVector & | dx | ) | const |
Provides the current correction vector dx .
[out] | dx | on exit, contains the current correction |
|
inline |
|
inline |
|
inline |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::getRawCorrection | ( | XVector & | dx | ) | const |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::getResidual | ( | XVector & | rr | ) | const |
|
inline |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::getS | ( | SMatrix & | S | ) | const |
Compute S .
This computes
S = \begin{bmatrix} 0 & I \end{bmatrix} \begin{bmatrix} A & C^T \\ C \end{bmatrix}^{-1} \begin{bmatrix} 0 \\ I \end{bmatrix}
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::getSchurResidual | ( | Vector & | lambda | ) | const |
Provides the Schur residual.
This solves the system
\begin{bmatrix} A & C^T \\ C \end{bmatrix} \begin{bmatrix} x+\delta x \\ \lambda +\delta\lambda \end{bmatrix} = \begin{bmatrix} f\\ 0 \end{bmatrix}
and returns the solution component \delta\lambda in the provided vector lambda
. The vectors x and \lambda are maintained internally as solution approximation.
[out] | lambda |
XVector Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::getSolution | ( | ) | const |
Provides the current solution vector x .
Note that the result is returned by value, which may be less efficient than providing a result vector to be filled.
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::getSolution | ( | XVector & | x | ) | const |
Provides the current solution vector x .
[out] | x | on exit, contains the current solution approximation |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::imposeDC_inActiveInterface | ( | Vector & | Aq, |
Vector & | q | ||
) |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::imposeDCq | ( | Vector & | q | ) |
std::vector< int > Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::neighbors | ( | ) | const |
Reports a list of its neighbor subdomains.
|
inline |
Prepares the exchange of interface values.
Has to be called before getInterfaceValues() or setInterfaceValues().
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::preRestrict | ( | ) |
First step of the restriction: Dirichlet solve.
This computes the transpose of the prolongation, i.e., with a suitable sorting of interior and bounday degrees of freedom, the
Dune::FieldVector< double, 2 > Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::prolongate | ( | ) |
Perform averaging on the shared interfaces, projecting the solution to the subspace of functions continuous across the interfaces.
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::restrict | ( | ) |
Performs the restriction of the residual vector.
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::set_solver_type | ( | bool | cg_solver | ) |
Set solver type.
|
inline |
|
inline |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::setCorrection | ( | Vector const & | c | ) |
Computes the correction vector [\delta x, \delta\lambda] .
|
inline |
|
inline |
|
inline |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::update_rhs | ( | XVector const & | f | ) |
update rhs for the current subdomain
f | rhs |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::updateIteration | ( | int | cg_itr | ) |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::updateSearchDirection | ( | double | beta | ) |
Update search direction q for cg solver.
Dune::FieldVector< double, 1 > Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::updateSearchDirectionNorm | ( | ) |
void Kaskade::BDDC::Subdomain< m, Scalar_, CoarseScalar, Transfer >::updateSolution | ( | Scalar | alpha | ) |
Updates the solution approximation.
|
static |