v0.12.1
Public Member Functions | Public Attributes | Friends | List of all members
ConstrainMatrixCtx Struct Reference

structure for projection matrices More...

#include <users_modules/basic_finite_elements/src/ConstrainMatrixCtx.hpp>

Collaboration diagram for ConstrainMatrixCtx:
[legend]

Public Member Functions

 ConstrainMatrixCtx (MoFEM::Interface &m_field, string x_problem, string y_problem, bool create_ksp=true, bool own_contrain_matrix=false)
 
 ConstrainMatrixCtx (MoFEM::Interface &m_field, VecScatter scatter, bool create_ksp=true, bool own_contrain_matrix=false)
 
virtual ~ConstrainMatrixCtx ()
 
MoFEMErrorCode initializeQorP (Vec x)
 initialize vectors and matrices for Q and P shell matrices, scattering is set based on x_problem and y_problem More...
 
MoFEMErrorCode initializeQTKQ ()
 initialize vectors and matrices for CTC+QTKQ shell matrices, scattering is set based on x_problem and y_problem More...
 
MoFEMErrorCode recalculateCTandCCT ()
 re-calculate CT and CCT if C matrix has been changed since initialization More...
 
MoFEMErrorCode recalculateCTC ()
 re-calculate CTC matrix has been changed since initialization More...
 
MoFEMErrorCode destroyQorP ()
 destroy sub-matrices used for shell matrices P, Q, R, RT More...
 
MoFEMErrorCode destroyQTKQ ()
 destroy sub-matrices used for shell matrix QTKQ More...
 

Public Attributes

MoFEM::InterfacemField
 
KSP kSP
 
Mat C
 
Mat CT
 
Mat CCT
 
Mat CTC
 
Mat K
 
Vec Cx
 
Vec CCTm1_Cx
 
Vec CT_CCTm1_Cx
 
Vec CTCx
 
Vec X
 
Vec Qx
 
Vec KQx
 
bool initQorP
 
bool initQTKQ
 
bool createKSP
 
bool createScatter
 
bool cancelKSPMonitor
 
bool ownConstrainMatrix
 
VecScatter sCatter
 
string xProblem
 
string yProblem
 
PetscLogEvent MOFEM_EVENT_projInit
 
PetscLogEvent MOFEM_EVENT_projQ
 
PetscLogEvent MOFEM_EVENT_projP
 
PetscLogEvent MOFEM_EVENT_projR
 
PetscLogEvent MOFEM_EVENT_projRT
 
PetscLogEvent MOFEM_EVENT_projCTC_QTKQ
 
PetscReal rTol
 
PetscReal absTol
 
PetscReal dTol
 
PetscInt maxIts
 

Friends

MoFEMErrorCode ProjectionMatrixMultOpQ (Mat Q, Vec x, Vec f)
 Multiplication operator for Q = I-CTC(CCT)^-1C. More...
 
MoFEMErrorCode ConstrainMatrixMultOpP (Mat P, Vec x, Vec f)
 Multiplication operator for P = CT(CCT)^-1C. More...
 
MoFEMErrorCode ConstrainMatrixMultOpR (Mat R, Vec x, Vec f)
 Multiplication operator for R = CT(CCT)^-1. More...
 
MoFEMErrorCode ConstrainMatrixMultOpRT (Mat RT, Vec x, Vec f)
 Multiplication operator for RT = (CCT)^-TC. More...
 
MoFEMErrorCode ConstrainMatrixMultOpCTC_QTKQ (Mat CTC_QTKQ, Vec x, Vec f)
 Multiplication operator for RT = (CCT)^-TC. More...
 
MoFEMErrorCode ConstrainMatrixDestroyOpPorQ ()
 
MoFEMErrorCode ConstrainMatrixDestroyOpQTKQ ()
 

Detailed Description

structure for projection matrices

Definition at line 28 of file ConstrainMatrixCtx.hpp.

Constructor & Destructor Documentation

◆ ConstrainMatrixCtx() [1/2]

ConstrainMatrixCtx::ConstrainMatrixCtx ( MoFEM::Interface m_field,
string  x_problem,
string  y_problem,
bool  create_ksp = true,
bool  own_contrain_matrix = false 
)

Construct data structure to build operators for projection matrices

User need to set matrix C to make it work

Parameters
x_problemproblem on which vector is projected
y_problemproblem used to construct projection matrices
create_kspcreate ksp solver otherwise user need to set it up

◆ ConstrainMatrixCtx() [2/2]

ConstrainMatrixCtx::ConstrainMatrixCtx ( MoFEM::Interface m_field,
VecScatter  scatter,
bool  create_ksp = true,
bool  own_contrain_matrix = false 
)

Definition at line 49 of file ConstrainMatrixCtx.cpp.

52  : mField(m_field), INIT_DATA_CONSTRAINMATRIXCTX, sCatter(scatter) {
53  PetscLogEventRegister("ProjectionInit", 0, &MOFEM_EVENT_projInit);
54  PetscLogEventRegister("ProjectionQ", 0, &MOFEM_EVENT_projQ);
55  PetscLogEventRegister("ProjectionP", 0, &MOFEM_EVENT_projP);
56  PetscLogEventRegister("ProjectionR", 0, &MOFEM_EVENT_projR);
57  PetscLogEventRegister("ProjectionRT", 0, &MOFEM_EVENT_projRT);
58  PetscLogEventRegister("ProjectionCTC_QTKQ", 0, &MOFEM_EVENT_projCTC_QTKQ);
59 }
#define INIT_DATA_CONSTRAINMATRIXCTX
PetscLogEvent MOFEM_EVENT_projRT
MoFEM::Interface & mField
PetscLogEvent MOFEM_EVENT_projCTC_QTKQ
PetscLogEvent MOFEM_EVENT_projP
PetscLogEvent MOFEM_EVENT_projInit
PetscLogEvent MOFEM_EVENT_projQ
PetscLogEvent MOFEM_EVENT_projR

◆ ~ConstrainMatrixCtx()

virtual ConstrainMatrixCtx::~ConstrainMatrixCtx ( )
virtual

Definition at line 71 of file ConstrainMatrixCtx.hpp.

71  {
72  ierr = destroyQorP();
73  CHKERRABORT(mField.get_comm(), ierr);
74  ierr = destroyQTKQ();
75  CHKERRABORT(mField.get_comm(), ierr);
76  if (ownConstrainMatrix) {
77  ierr = MatDestroy(&C);
78  CHKERRABORT(mField.get_comm(), ierr);
79  }
80  };
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
MoFEMErrorCode destroyQTKQ()
destroy sub-matrices used for shell matrix QTKQ
MoFEMErrorCode destroyQorP()
destroy sub-matrices used for shell matrices P, Q, R, RT
virtual MPI_Comm & get_comm() const =0

Member Function Documentation

◆ destroyQorP()

MoFEMErrorCode ConstrainMatrixCtx::destroyQorP ( )

destroy sub-matrices used for shell matrices P, Q, R, RT

Definition at line 110 of file ConstrainMatrixCtx.cpp.

110  {
112  if (initQorP)
114  CHKERR MatDestroy(&CT);
115  CHKERR MatDestroy(&CCT);
116  if (createKSP) {
117  CHKERR KSPDestroy(&kSP);
118  }
119  CHKERR VecDestroy(&X);
120  CHKERR VecDestroy(&Cx);
121  CHKERR VecDestroy(&CCTm1_Cx);
122  CHKERR VecDestroy(&CT_CCTm1_Cx);
123  if (createScatter) {
124  CHKERR VecScatterDestroy(&sCatter);
125  }
126  initQorP = true;
128 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:440
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:339
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:409
#define CHKERR
Inline error check.
Definition: definitions.h:528

◆ destroyQTKQ()

MoFEMErrorCode ConstrainMatrixCtx::destroyQTKQ ( )

destroy sub-matrices used for shell matrix QTKQ

Definition at line 167 of file ConstrainMatrixCtx.cpp.

167  {
169  if (initQTKQ)
171  CHKERR MatDestroy(&CTC);
172  CHKERR VecDestroy(&Qx);
173  CHKERR VecDestroy(&KQx);
174  CHKERR VecDestroy(&CTCx);
175  initQTKQ = true;
177 }

◆ initializeQorP()

MoFEMErrorCode ConstrainMatrixCtx::initializeQorP ( Vec  x)

initialize vectors and matrices for Q and P shell matrices, scattering is set based on x_problem and y_problem

Parameters
xis a vector from problem x

Definition at line 61 of file ConstrainMatrixCtx.cpp.

61  {
63  if (initQorP) {
64  initQorP = false;
65 
66  PetscLogEventBegin(MOFEM_EVENT_projInit, 0, 0, 0, 0);
67  CHKERR MatTranspose(C, MAT_INITIAL_MATRIX, &CT);
68  // need to be calculated when C is changed
69  CHKERR MatTransposeMatMult(CT, CT, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &CCT);
70  if (createKSP) {
71  CHKERR KSPCreate(mField.get_comm(), &kSP);
72  // need to be recalculated when C is changed
73  CHKERR KSPSetOperators(kSP, CCT, CCT);
74  CHKERR KSPSetFromOptions(kSP);
75  CHKERR KSPSetInitialGuessKnoll(kSP, PETSC_TRUE);
76  CHKERR KSPGetTolerances(kSP, &rTol, &absTol, &dTol, &maxIts);
77  CHKERR KSPSetUp(kSP);
78  if (cancelKSPMonitor) {
79  CHKERR KSPMonitorCancel(kSP);
80  }
81  }
82 #if PETSC_VERSION_GE(3, 5, 3)
83  CHKERR MatCreateVecs(C, &X, PETSC_NULL);
84  CHKERR MatCreateVecs(C, PETSC_NULL, &Cx);
85  CHKERR MatCreateVecs(CCT, PETSC_NULL, &CCTm1_Cx);
86 #else
87  CHKERR MatGetVecs(C, &X, PETSC_NULL);
88  CHKERR MatGetVecs(C, PETSC_NULL, &Cx);
89  CHKERR MatGetVecs(CCT, PETSC_NULL, &CCTm1_Cx);
90 #endif
91  CHKERR VecDuplicate(X, &CT_CCTm1_Cx);
92  if (createScatter) {
93  CHKERR mField.getInterface<VecManager>()->vecScatterCreate(
94  x, xProblem, ROW, X, yProblem, COL, &sCatter);
95  }
96  PetscLogEventEnd(MOFEM_EVENT_projInit, 0, 0, 0, 0);
97  }
99 }
@ COL
Definition: definitions.h:116
@ ROW
Definition: definitions.h:116
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:33

◆ initializeQTKQ()

MoFEMErrorCode ConstrainMatrixCtx::initializeQTKQ ( )

initialize vectors and matrices for CTC+QTKQ shell matrices, scattering is set based on x_problem and y_problem

Definition at line 130 of file ConstrainMatrixCtx.cpp.

130  {
132  if (initQTKQ) {
133  initQTKQ = false;
134  PetscLogEventBegin(MOFEM_EVENT_projInit, 0, 0, 0, 0);
135  // need to be recalculated when C is changed
136  CHKERR MatTransposeMatMult(C, C, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &CTC);
137  if (debug) {
138  // MatView(CCT,PETSC_VIEWER_DRAW_WORLD);
139  int m, n;
140  MatGetSize(CCT, &m, &n);
141  PetscPrintf(mField.get_comm(), "CCT size (%d,%d)\n", m, n);
142  // std::string wait;
143  // std::cin >> wait;
144  }
145 #if PETSC_VERSION_GE(3, 5, 3)
146  CHKERR MatCreateVecs(K, &Qx, PETSC_NULL);
147  CHKERR MatCreateVecs(K, PETSC_NULL, &KQx);
148  CHKERR MatCreateVecs(CTC, PETSC_NULL, &CTCx);
149 #else
150  CHKERR MatGetVecs(K, &Qx, PETSC_NULL);
151  CHKERR MatGetVecs(K, PETSC_NULL, &KQx);
152  CHKERR MatGetVecs(CTC, PETSC_NULL, &CTCx);
153 #endif
154  PetscLogEventEnd(MOFEM_EVENT_projInit, 0, 0, 0, 0);
155  }
157 }
static Index< 'm', 3 > m
static Index< 'n', 3 > n
static const bool debug

◆ recalculateCTandCCT()

MoFEMErrorCode ConstrainMatrixCtx::recalculateCTandCCT ( )

re-calculate CT and CCT if C matrix has been changed since initialization

Definition at line 101 of file ConstrainMatrixCtx.cpp.

101  {
103  if (initQorP)
105  CHKERR MatTranspose(C, MAT_REUSE_MATRIX, &CT);
106  CHKERR MatTransposeMatMult(CT, CT, MAT_REUSE_MATRIX, PETSC_DEFAULT, &CCT);
108 }

◆ recalculateCTC()

MoFEMErrorCode ConstrainMatrixCtx::recalculateCTC ( )

re-calculate CTC matrix has been changed since initialization

Definition at line 159 of file ConstrainMatrixCtx.cpp.

159  {
161  if (initQTKQ)
163  CHKERR MatTransposeMatMult(C, C, MAT_REUSE_MATRIX, PETSC_DEFAULT, &CTC);
165 }

Friends And Related Function Documentation

◆ ConstrainMatrixDestroyOpPorQ

MoFEMErrorCode ConstrainMatrixDestroyOpPorQ ( )
friend

◆ ConstrainMatrixDestroyOpQTKQ

MoFEMErrorCode ConstrainMatrixDestroyOpQTKQ ( )
friend

◆ ConstrainMatrixMultOpCTC_QTKQ

MoFEMErrorCode ConstrainMatrixMultOpCTC_QTKQ ( Mat  CTC_QTKQ,
Vec  x,
Vec  f 
)
friend

Multiplication operator for RT = (CCT)^-TC.

Mat CTC_QTKQ; //for problem
projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&CTC_QTKQ);
MatShellSetOperation(CTC_QTKQ,MATOP_MULT,(void(*)(void))ConstrainMatrixMultOpCTC_QTKQ);
MatShellSetOperation(CTC_QTKQ,MATOP_DESTROY,(void(*)(void))ConstrainMatrixDestroyOpQTKQ);
static Index< 'M', 3 > M
structure for projection matrices
friend MoFEMErrorCode ConstrainMatrixDestroyOpQTKQ()
friend MoFEMErrorCode ConstrainMatrixMultOpCTC_QTKQ(Mat CTC_QTKQ, Vec x, Vec f)
Multiplication operator for RT = (CCT)^-TC.

Definition at line 277 of file ConstrainMatrixCtx.cpp.

277  {
279  void *void_ctx;
280  CHKERR MatShellGetContext(CTC_QTKQ, &void_ctx);
281  ConstrainMatrixCtx *ctx = (ConstrainMatrixCtx *)void_ctx;
282  PetscLogEventBegin(ctx->MOFEM_EVENT_projCTC_QTKQ, 0, 0, 0, 0);
283  Mat Q;
284  int M, N, m, n;
285  CHKERR MatGetSize(ctx->K, &M, &N);
286  CHKERR MatGetLocalSize(ctx->K, &m, &n);
287  CHKERR MatCreateShell(ctx->mField.get_comm(), m, n, M, N, ctx, &Q);
288  CHKERR MatShellSetOperation(Q, MATOP_MULT,
289  (void (*)(void))ProjectionMatrixMultOpQ);
290  CHKERR ctx->initializeQTKQ();
291  CHKERR MatMult(Q, x, ctx->Qx);
292  CHKERR MatMult(ctx->K, ctx->Qx, ctx->KQx);
293  CHKERR MatMult(Q, ctx->KQx, f);
294  CHKERR VecScatterBegin(ctx->sCatter, x, ctx->X, INSERT_VALUES,
295  SCATTER_FORWARD);
296  CHKERR VecScatterEnd(ctx->sCatter, x, ctx->X, INSERT_VALUES, SCATTER_FORWARD);
297  CHKERR MatMult(ctx->CTC, ctx->X, ctx->CTCx);
298  CHKERR VecScatterBegin(ctx->sCatter, ctx->CTCx, f, ADD_VALUES,
299  SCATTER_REVERSE);
300  CHKERR VecScatterEnd(ctx->sCatter, ctx->CTCx, f, ADD_VALUES, SCATTER_REVERSE);
301  CHKERR MatDestroy(&Q);
302  PetscLogEventEnd(ctx->MOFEM_EVENT_projCTC_QTKQ, 0, 0, 0, 0);
304 }
const int N
Definition: speed_test.cpp:3
MoFEMErrorCode initializeQTKQ()
initialize vectors and matrices for CTC+QTKQ shell matrices, scattering is set based on x_problem and...
friend MoFEMErrorCode ProjectionMatrixMultOpQ(Mat Q, Vec x, Vec f)
Multiplication operator for Q = I-CTC(CCT)^-1C.

◆ ConstrainMatrixMultOpP

MoFEMErrorCode ConstrainMatrixMultOpP ( Mat  P,
Vec  x,
Vec  f 
)
friend

Multiplication operator for P = CT(CCT)^-1C.

Mat P; //for problem
projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&P);
MatShellSetOperation(P,MATOP_MULT,(void(*)(void))ConstrainMatrixMultOpP);
friend MoFEMErrorCode ConstrainMatrixMultOpP(Mat P, Vec x, Vec f)
Multiplication operator for P = CT(CCT)^-1C.

Definition at line 213 of file ConstrainMatrixCtx.cpp.

213  {
215  void *void_ctx;
216  CHKERR MatShellGetContext(P, &void_ctx);
217  ConstrainMatrixCtx *ctx = (ConstrainMatrixCtx *)void_ctx;
218  PetscLogEventBegin(ctx->MOFEM_EVENT_projP, 0, 0, 0, 0);
219  CHKERR ctx->initializeQorP(x);
220  CHKERR VecScatterBegin(ctx->sCatter, x, ctx->X, INSERT_VALUES,
221  SCATTER_FORWARD);
222  CHKERR VecScatterEnd(ctx->sCatter, x, ctx->X, INSERT_VALUES, SCATTER_FORWARD);
223  CHKERR MatMult(ctx->C, ctx->X, ctx->Cx);
224  CHKERR KSPSolve(ctx->kSP, ctx->Cx, ctx->CCTm1_Cx);
225  CHKERR MatMult(ctx->CT, ctx->CCTm1_Cx, ctx->CT_CCTm1_Cx);
226  CHKERR VecZeroEntries(f);
227  CHKERR VecGhostUpdateBegin(f, INSERT_VALUES, SCATTER_FORWARD);
228  CHKERR VecGhostUpdateEnd(f, INSERT_VALUES, SCATTER_FORWARD);
229  CHKERR VecScatterBegin(ctx->sCatter, ctx->CT_CCTm1_Cx, f, INSERT_VALUES,
230  SCATTER_REVERSE);
231  CHKERR VecScatterEnd(ctx->sCatter, ctx->CT_CCTm1_Cx, f, INSERT_VALUES,
232  SCATTER_REVERSE);
233  PetscLogEventEnd(ctx->MOFEM_EVENT_projP, 0, 0, 0, 0);
235 }
MoFEMErrorCode initializeQorP(Vec x)
initialize vectors and matrices for Q and P shell matrices, scattering is set based on x_problem and ...

◆ ConstrainMatrixMultOpR

MoFEMErrorCode ConstrainMatrixMultOpR ( Mat  R,
Vec  x,
Vec  f 
)
friend

Multiplication operator for R = CT(CCT)^-1.

Mat R; //for problem
projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&R);
MatShellSetOperation(R,MATOP_MULT,(void(*)(void))ConstrainMatrixMultOpR);
friend MoFEMErrorCode ConstrainMatrixMultOpR(Mat R, Vec x, Vec f)
Multiplication operator for R = CT(CCT)^-1.

Definition at line 237 of file ConstrainMatrixCtx.cpp.

237  {
239  void *void_ctx;
240  CHKERR MatShellGetContext(R, &void_ctx);
241  ConstrainMatrixCtx *ctx = (ConstrainMatrixCtx *)void_ctx;
242  PetscLogEventBegin(ctx->MOFEM_EVENT_projR, 0, 0, 0, 0);
243  if (ctx->initQorP)
244  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE,
245  "you have to call first initQorP or use Q matrix");
246  CHKERR KSPSolve(ctx->kSP, x, ctx->CCTm1_Cx);
247  CHKERR MatMult(ctx->CT, ctx->CCTm1_Cx, ctx->CT_CCTm1_Cx);
248  CHKERR VecZeroEntries(f);
249  CHKERR VecGhostUpdateBegin(f, INSERT_VALUES, SCATTER_FORWARD);
250  CHKERR VecGhostUpdateEnd(f, INSERT_VALUES, SCATTER_FORWARD);
251  CHKERR VecScatterBegin(ctx->sCatter, ctx->CT_CCTm1_Cx, f, INSERT_VALUES,
252  SCATTER_REVERSE);
253  CHKERR VecScatterEnd(ctx->sCatter, ctx->CT_CCTm1_Cx, f, INSERT_VALUES,
254  SCATTER_REVERSE);
255  PetscLogEventEnd(ctx->MOFEM_EVENT_projR, 0, 0, 0, 0);
257 }
@ MOFEM_IMPOSIBLE_CASE
Definition: definitions.h:46

◆ ConstrainMatrixMultOpRT

MoFEMErrorCode ConstrainMatrixMultOpRT ( Mat  RT,
Vec  x,
Vec  f 
)
friend

Multiplication operator for RT = (CCT)^-TC.

Mat RT; //for problem
projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&RT);
MatShellSetOperation(RT,MATOP_MULT,(void(*)(void))ConstrainMatrixMultOpRT);
friend MoFEMErrorCode ConstrainMatrixMultOpRT(Mat RT, Vec x, Vec f)
Multiplication operator for RT = (CCT)^-TC.

Definition at line 259 of file ConstrainMatrixCtx.cpp.

259  {
261  void *void_ctx;
262  CHKERR MatShellGetContext(RT, &void_ctx);
263  ConstrainMatrixCtx *ctx = (ConstrainMatrixCtx *)void_ctx;
264  PetscLogEventBegin(ctx->MOFEM_EVENT_projRT, 0, 0, 0, 0);
265  if (ctx->initQorP)
266  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE,
267  "you have to call first initQorP or use Q matrix");
268  CHKERR VecScatterBegin(ctx->sCatter, x, ctx->X, INSERT_VALUES,
269  SCATTER_FORWARD);
270  CHKERR VecScatterEnd(ctx->sCatter, x, ctx->X, INSERT_VALUES, SCATTER_FORWARD);
271  CHKERR MatMult(ctx->C, ctx->X, ctx->Cx);
272  CHKERR KSPSolve(ctx->kSP, ctx->Cx, f);
273  PetscLogEventEnd(ctx->MOFEM_EVENT_projRT, 0, 0, 0, 0);
275 }

◆ ProjectionMatrixMultOpQ

MoFEMErrorCode ProjectionMatrixMultOpQ ( Mat  Q,
Vec  x,
Vec  f 
)
friend

Multiplication operator for Q = I-CTC(CCT)^-1C.

Mat Q; //for problem
projection_matrix_ctx(m_field,problem_name,contrains_problem_name);
CHKERR MatCreateShell(PETSC_COMM_WORLD,m,m,M,M,&projection_matrix_ctx,&Q);
MatShellSetOperation(Q,MATOP_MULT,(void(*)(void))ProjectionMatrixMultOpQ);
MatShellSetOperation(Q,MATOP_DESTROY,(void(*)(void))ConstrainMatrixDestroyOpPorQ);
friend MoFEMErrorCode ConstrainMatrixDestroyOpPorQ()

Definition at line 179 of file ConstrainMatrixCtx.cpp.

179  {
181  void *void_ctx;
182  CHKERR MatShellGetContext(Q, &void_ctx);
183  ConstrainMatrixCtx *ctx = (ConstrainMatrixCtx *)void_ctx;
184  PetscLogEventBegin(ctx->MOFEM_EVENT_projQ, 0, 0, 0, 0);
185  CHKERR ctx->initializeQorP(x);
186  CHKERR VecCopy(x, f);
187  CHKERR VecScatterBegin(ctx->sCatter, x, ctx->X, INSERT_VALUES,
188  SCATTER_FORWARD);
189  CHKERR VecScatterEnd(ctx->sCatter, x, ctx->X, INSERT_VALUES, SCATTER_FORWARD);
190  if (debug) {
191  // CHKERR VecView(ctx->X,PETSC_VIEWER_STDOUT_WORLD);
192  CHKERR VecScatterBegin(ctx->sCatter, ctx->X, f, INSERT_VALUES,
193  SCATTER_REVERSE);
194  CHKERR VecScatterEnd(ctx->sCatter, ctx->X, f, INSERT_VALUES,
195  SCATTER_REVERSE);
196  PetscBool flg;
197  CHKERR VecEqual(x, f, &flg);
198  if (flg == PETSC_FALSE)
199  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "scatter is not working");
200  }
201  CHKERR MatMult(ctx->C, ctx->X, ctx->Cx);
202  CHKERR KSPSolve(ctx->kSP, ctx->Cx, ctx->CCTm1_Cx);
203  CHKERR MatMult(ctx->CT, ctx->CCTm1_Cx, ctx->CT_CCTm1_Cx);
204  CHKERR VecScale(ctx->CT_CCTm1_Cx, -1);
205  CHKERR VecScatterBegin(ctx->sCatter, ctx->CT_CCTm1_Cx, f, ADD_VALUES,
206  SCATTER_REVERSE);
207  CHKERR VecScatterEnd(ctx->sCatter, ctx->CT_CCTm1_Cx, f, ADD_VALUES,
208  SCATTER_REVERSE);
209  PetscLogEventEnd(ctx->MOFEM_EVENT_projQ, 0, 0, 0, 0);
211 }

Member Data Documentation

◆ absTol

PetscReal ConstrainMatrixCtx::absTol

Definition at line 82 of file ConstrainMatrixCtx.hpp.

◆ C

Mat ConstrainMatrixCtx::C

Definition at line 33 of file ConstrainMatrixCtx.hpp.

◆ cancelKSPMonitor

bool ConstrainMatrixCtx::cancelKSPMonitor

Definition at line 39 of file ConstrainMatrixCtx.hpp.

◆ CCT

Mat ConstrainMatrixCtx::CCT

Definition at line 33 of file ConstrainMatrixCtx.hpp.

◆ CCTm1_Cx

Vec ConstrainMatrixCtx::CCTm1_Cx

Definition at line 34 of file ConstrainMatrixCtx.hpp.

◆ createKSP

bool ConstrainMatrixCtx::createKSP

Definition at line 37 of file ConstrainMatrixCtx.hpp.

◆ createScatter

bool ConstrainMatrixCtx::createScatter

Definition at line 38 of file ConstrainMatrixCtx.hpp.

◆ CT

Mat ConstrainMatrixCtx::CT

Definition at line 33 of file ConstrainMatrixCtx.hpp.

◆ CT_CCTm1_Cx

Vec ConstrainMatrixCtx::CT_CCTm1_Cx

Definition at line 34 of file ConstrainMatrixCtx.hpp.

◆ CTC

Mat ConstrainMatrixCtx::CTC

Definition at line 33 of file ConstrainMatrixCtx.hpp.

◆ CTCx

Vec ConstrainMatrixCtx::CTCx

Definition at line 34 of file ConstrainMatrixCtx.hpp.

◆ Cx

Vec ConstrainMatrixCtx::Cx

Definition at line 34 of file ConstrainMatrixCtx.hpp.

◆ dTol

PetscReal ConstrainMatrixCtx::dTol

Definition at line 82 of file ConstrainMatrixCtx.hpp.

◆ initQorP

bool ConstrainMatrixCtx::initQorP

Definition at line 36 of file ConstrainMatrixCtx.hpp.

◆ initQTKQ

bool ConstrainMatrixCtx::initQTKQ

Definition at line 36 of file ConstrainMatrixCtx.hpp.

◆ K

Mat ConstrainMatrixCtx::K

Definition at line 33 of file ConstrainMatrixCtx.hpp.

◆ KQx

Vec ConstrainMatrixCtx::KQx

Definition at line 35 of file ConstrainMatrixCtx.hpp.

◆ kSP

KSP ConstrainMatrixCtx::kSP

Definition at line 32 of file ConstrainMatrixCtx.hpp.

◆ maxIts

PetscInt ConstrainMatrixCtx::maxIts

Definition at line 83 of file ConstrainMatrixCtx.hpp.

◆ mField

MoFEM::Interface& ConstrainMatrixCtx::mField

Definition at line 30 of file ConstrainMatrixCtx.hpp.

◆ MOFEM_EVENT_projCTC_QTKQ

PetscLogEvent ConstrainMatrixCtx::MOFEM_EVENT_projCTC_QTKQ

Definition at line 53 of file ConstrainMatrixCtx.hpp.

◆ MOFEM_EVENT_projInit

PetscLogEvent ConstrainMatrixCtx::MOFEM_EVENT_projInit

Definition at line 48 of file ConstrainMatrixCtx.hpp.

◆ MOFEM_EVENT_projP

PetscLogEvent ConstrainMatrixCtx::MOFEM_EVENT_projP

Definition at line 50 of file ConstrainMatrixCtx.hpp.

◆ MOFEM_EVENT_projQ

PetscLogEvent ConstrainMatrixCtx::MOFEM_EVENT_projQ

Definition at line 49 of file ConstrainMatrixCtx.hpp.

◆ MOFEM_EVENT_projR

PetscLogEvent ConstrainMatrixCtx::MOFEM_EVENT_projR

Definition at line 51 of file ConstrainMatrixCtx.hpp.

◆ MOFEM_EVENT_projRT

PetscLogEvent ConstrainMatrixCtx::MOFEM_EVENT_projRT

Definition at line 52 of file ConstrainMatrixCtx.hpp.

◆ ownConstrainMatrix

bool ConstrainMatrixCtx::ownConstrainMatrix

Definition at line 40 of file ConstrainMatrixCtx.hpp.

◆ Qx

Vec ConstrainMatrixCtx::Qx

Definition at line 35 of file ConstrainMatrixCtx.hpp.

◆ rTol

PetscReal ConstrainMatrixCtx::rTol

Definition at line 82 of file ConstrainMatrixCtx.hpp.

◆ sCatter

VecScatter ConstrainMatrixCtx::sCatter

Definition at line 45 of file ConstrainMatrixCtx.hpp.

◆ X

Vec ConstrainMatrixCtx::X

Definition at line 35 of file ConstrainMatrixCtx.hpp.

◆ xProblem

string ConstrainMatrixCtx::xProblem

Definition at line 46 of file ConstrainMatrixCtx.hpp.

◆ yProblem

string ConstrainMatrixCtx::yProblem

Definition at line 46 of file ConstrainMatrixCtx.hpp.


The documentation for this struct was generated from the following files: