v0.15.0
Loading...
Searching...
No Matches
ConstrainMatrixCtx.hpp
Go to the documentation of this file.
1/** \file ConstrainMatrixCtx.hpp
2 *
3 * Can be used if constrains are linear, i.e. are not function of time.
4 */
5
6
7
8#ifndef __PROJECTION_MATRIX_CTX_HPP__
9#define __PROJECTION_MATRIX_CTX_HPP__
10
11namespace MoFEM {
12/**
13 * \brief structure for projection matrices
14 * \ingroup projection_matrix
15 *
16 */
17
19
21
22 KSP kSP;
23 Mat C, CT, CCT, CTC, K;
25 Vec X, Qx, KQx;
31
32 // Scatter is created form problem_x to problem_y, or scatter is given
33 // in the constructor
34
35 VecScatter sCatter;
37
38 PetscLogEvent MOFEM_EVENT_projInit;
39 PetscLogEvent MOFEM_EVENT_projQ;
40 PetscLogEvent MOFEM_EVENT_projP;
41 PetscLogEvent MOFEM_EVENT_projR;
42 PetscLogEvent MOFEM_EVENT_projRT;
44
45 /**
46 * Construct data structure to build operators for projection matrices
47 *
48 * User need to set matrix C to make it work
49 *
50 * \param x_problem problem on which vector is projected
51 * \param y_problem problem used to construct projection matrices
52 * \param create_ksp create ksp solver otherwise user need to set it up
53 */
54 ConstrainMatrixCtx(MoFEM::Interface &m_field, string x_problem,
55 string y_problem, bool create_ksp = true,
56 bool own_contrain_matrix = false);
57
58 ConstrainMatrixCtx(MoFEM::Interface &m_field, VecScatter scatter,
59 bool create_ksp = true, bool own_contrain_matrix = false);
60
62 ierr = destroyQorP();
63 CHKERRABORT(mField.get_comm(), ierr);
64 ierr = destroyQTKQ();
65 CHKERRABORT(mField.get_comm(), ierr);
67 ierr = MatDestroy(&C);
68 CHKERRABORT(mField.get_comm(), ierr);
69 }
70 };
71
72 PetscReal rTol, absTol, dTol;
73 PetscInt maxIts;
74
75 /**
76 * \brief initialize vectors and matrices for Q and P shell matrices,
77 * scattering is set based on x_problem and y_problem
78 *
79 * \param x is a vector from problem x
80 */
82
83 /**
84 * \brief initialize vectors and matrices for CTC+QTKQ shell matrices,
85 * scattering is set based on x_problem and y_problem
86 */
88
89 /**
90 * \brief re-calculate CT and CCT if C matrix has been changed since
91 * initialization
92 */
94
95 /**
96 * \brief re-calculate CTC matrix has been changed since initialization
97 */
99
100 /**
101 * \brief destroy sub-matrices used for shell matrices P, Q, R, RT
102 */
104
105 /**
106 * \brief destroy sub-matrices used for shell matrix QTKQ
107 */
109
110 friend MoFEMErrorCode ProjectionMatrixMultOpQ(Mat Q, Vec x, Vec f);
111 friend MoFEMErrorCode ConstrainMatrixMultOpP(Mat P, Vec x, Vec f);
112 friend MoFEMErrorCode ConstrainMatrixMultOpR(Mat R, Vec x, Vec f);
113 friend MoFEMErrorCode ConstrainMatrixMultOpRT(Mat RT, Vec x, Vec f);
114 friend MoFEMErrorCode ConstrainMatrixMultOpCTC_QTKQ(Mat CTC_QTKQ, Vec x,
115 Vec f);
116
119};
120
121/**
122 * \brief Multiplication operator for Q = I-CTC(CCT)^-1C
123 *
124 * \code
125 * Mat Q; //for problem
126 * ConstrainMatrixCtx
127 projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
128 * CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&Q);
129 * CHKERR
130 MatShellSetOperation(Q,MATOP_MULT,(void(*)(void))ProjectionMatrixMultOpQ);
131 * CHKERR
132 MatShellSetOperation(Q,MATOP_DESTROY,(void(*)(void))ConstrainMatrixDestroyOpPorQ);
133 *
134 * \endcode
135
136 * \ingroup projection_matrix
137
138 */
139MoFEMErrorCode ProjectionMatrixMultOpQ(Mat Q, Vec x, Vec f);
140
141/**
142 * \brief Multiplication operator for P = CT(CCT)^-1C
143 *
144 * \code
145 * Mat P; //for problem
146 * ConstrainMatrixCtx
147 projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
148 * CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&P);
149 * CHKERR
150 MatShellSetOperation(P,MATOP_MULT,(void(*)(void))ConstrainMatrixMultOpP);
151 *
152 * \endcode
153
154 * \ingroup projection_matrix
155
156 */
157MoFEMErrorCode ConstrainMatrixMultOpP(Mat P, Vec x, Vec f);
158
159/**
160 * \brief Multiplication operator for R = CT(CCT)^-1
161 *
162 * \code
163 * Mat R; //for problem
164 * ConstrainMatrixCtx
165 projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
166 * CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&R);
167 * CHKERR
168 MatShellSetOperation(R,MATOP_MULT,(void(*)(void))ConstrainMatrixMultOpR);
169 *
170 * \endcode
171
172 * \ingroup projection_matrix
173
174 */
175MoFEMErrorCode ConstrainMatrixMultOpR(Mat R, Vec x, Vec f);
176
177/**
178 * \brief Multiplication operator for RT = (CCT)^-TC
179 *
180 * \code
181 * Mat RT; //for problem
182 * ConstrainMatrixCtx
183 projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
184 * CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&RT);
185 * CHKERR
186 MatShellSetOperation(RT,MATOP_MULT,(void(*)(void))ConstrainMatrixMultOpRT);
187 *
188 * \endcode
189
190 * \ingroup projection_matrix
191
192 */
193MoFEMErrorCode ConstrainMatrixMultOpRT(Mat RT, Vec x, Vec f);
194
195/**
196 * \brief Multiplication operator for RT = (CCT)^-TC
197 *
198 * \code
199 * Mat CTC_QTKQ; //for problem
200 * ConstrainMatrixCtx
201 projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
202 * CHKERR
203 MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&CTC_QTKQ);
204 * CHKERR
205 MatShellSetOperation(CTC_QTKQ,MATOP_MULT,(void(*)(void))ConstrainMatrixMultOpCTC_QTKQ);
206 * CHKERR
207 MatShellSetOperation(CTC_QTKQ,MATOP_DESTROY,(void(*)(void))ConstrainMatrixDestroyOpQTKQ);
208 *
209 * \endcode
210 *
211
212 * \ingroup projection_matrix
213
214 */
215MoFEMErrorCode ConstrainMatrixMultOpCTC_QTKQ(Mat CTC_QTKQ, Vec x, Vec f);
216
217/**
218 * \brief Destroy shell matrix Q
219 *
220 * \code
221 * Mat Q; //for problem
222 * ConstrainMatrixCtx
223 projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
224 * CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&Q);
225 * CHKERR
226 MatShellSetOperation(Q,MATOP_MULT,(void(*)(void))ProjectionMatrixMultOpQ);
227 * CHKERR
228 MatShellSetOperation(Q,MATOP_DESTROY,(void(*)(void))ConstrainMatrixDestroyOpPorQ);
229 *
230 * \endcode
231
232 * \ingroup projection_matrix
233
234 */
236
237/**
238 * \brief Destroy shell matrix
239 *
240 * \code
241 * Mat CTC_QTKQ; //for problem
242 * ConstrainMatrixCtx
243 projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
244 * CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&Q);
245 * CHKERR
246 MatShellSetOperation(Q,MATOP_MULT,(void(*)(void))ConstrainMatrixMultOpCTC_QTKQ);
247 * CHKERR
248 MatShellSetOperation(Q,MATOP_DESTROY,(void(*)(void))mat_destroy_QTKQ);
249 *
250 * \endcode
251
252 * \ingroup projection_matrix
253
254 */
256
257} // namespace MoFEM
258
259#endif // __PROJECTION_MATRIX_CTX_HPP__
260
261/**
262 \defgroup projection_matrix Constrain Projection Matrix
263 \ingroup mofem_basic_finite_elements
264**/
@ R
MoFEMErrorCode ConstrainMatrixMultOpCTC_QTKQ(Mat CTC_QTKQ, Vec x, Vec f)
Multiplication operator for RT = (CCT)^-TC.
MoFEMErrorCode ConstrainMatrixMultOpP(Mat P, Vec x, Vec f)
Multiplication operator for P = CT(CCT)^-1C.
MoFEMErrorCode ConstrainMatrixDestroyOpPorQ(Mat Q)
Destroy shell matrix Q.
MoFEMErrorCode ProjectionMatrixMultOpQ(Mat Q, Vec x, Vec f)
Multiplication operator for Q = I-CTC(CCT)^-1C.
MoFEMErrorCode ConstrainMatrixMultOpR(Mat R, Vec x, Vec f)
Multiplication operator for R = CT(CCT)^-1.
MoFEMErrorCode ConstrainMatrixDestroyOpQTKQ(Mat QTKQ)
Destroy shell matrix.
MoFEMErrorCode ConstrainMatrixMultOpRT(Mat RT, Vec x, Vec f)
Multiplication operator for RT = (CCT)^-TC.
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
structure for projection matrices
MoFEMErrorCode destroyQTKQ()
destroy sub-matrices used for shell matrix QTKQ
friend MoFEMErrorCode ConstrainMatrixMultOpR(Mat R, Vec x, Vec f)
Multiplication operator for R = CT(CCT)^-1.
friend MoFEMErrorCode ConstrainMatrixMultOpP(Mat P, Vec x, Vec f)
Multiplication operator for P = CT(CCT)^-1C.
friend MoFEMErrorCode ConstrainMatrixMultOpRT(Mat RT, Vec x, Vec f)
Multiplication operator for RT = (CCT)^-TC.
friend MoFEMErrorCode ConstrainMatrixDestroyOpQTKQ()
ConstrainMatrixCtx(MoFEM::Interface &m_field, string x_problem, string y_problem, bool create_ksp=true, bool own_contrain_matrix=false)
MoFEMErrorCode recalculateCTC()
re-calculate CTC matrix has been changed since initialization
friend MoFEMErrorCode ConstrainMatrixMultOpCTC_QTKQ(Mat CTC_QTKQ, Vec x, Vec f)
Multiplication operator for RT = (CCT)^-TC.
MoFEMErrorCode destroyQorP()
destroy sub-matrices used for shell matrices P, Q, R, RT
MoFEMErrorCode initializeQTKQ()
initialize vectors and matrices for CTC+QTKQ shell matrices, scattering is set based on x_problem and...
friend MoFEMErrorCode ConstrainMatrixDestroyOpPorQ()
MoFEMErrorCode initializeQorP(Vec x)
initialize vectors and matrices for Q and P shell matrices, scattering is set based on x_problem and ...
MoFEMErrorCode recalculateCTandCCT()
re-calculate CT and CCT if C matrix has been changed since initialization
friend MoFEMErrorCode ProjectionMatrixMultOpQ(Mat Q, Vec x, Vec f)
Multiplication operator for Q = I-CTC(CCT)^-1C.
virtual MPI_Comm & get_comm() const =0
Deprecated interface functions.