v0.13.2
Loading...
Searching...
No Matches
SnesCtx.cpp
Go to the documentation of this file.
1
2
3namespace MoFEM {
4
14}
15
18 loops_to_do_Mat.clear();
19 loops_to_do_Rhs.clear();
20 preProcess_Mat.clear();
21 postProcess_Mat.clear();
22 preProcess_Rhs.clear();
23 postProcess_Rhs.clear();
25}
26
27PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx) {
28 SnesCtx *snes_ctx = (SnesCtx *)ctx;
29 // PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
31 PetscLogEventBegin(snes_ctx->MOFEM_EVENT_SnesRhs, 0, 0, 0, 0);
32 CHKERR VecGhostUpdateBegin(x, INSERT_VALUES, SCATTER_FORWARD);
33 CHKERR VecGhostUpdateEnd(x, INSERT_VALUES, SCATTER_FORWARD);
34 if (snes_ctx->vErify) {
35 // Verify finite elements, check for not a number
36 CHKERR VecAssemblyBegin(f);
37 CHKERR VecAssemblyEnd(f);
38 MPI_Comm comm = PetscObjectComm((PetscObject)f);
39 PetscSynchronizedPrintf(comm, "SNES Verify x\n");
40 const Problem *prb_ptr;
41 CHKERR snes_ctx->mField.get_problem(snes_ctx->problemName, &prb_ptr);
42 CHKERR snes_ctx->mField.getInterface<Tools>()->checkVectorForNotANumber(
43 prb_ptr, COL, x);
44 }
45 CHKERR snes_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
46 snes_ctx->problemName, COL, x, INSERT_VALUES, SCATTER_REVERSE);
47
48 auto zero_ghost_vec = [](Vec g) {
50 Vec l;
51 CHKERR VecGhostGetLocalForm(g, &l);
52 double *a;
53 CHKERR VecGetArray(l, &a);
54 int s;
55 CHKERR VecGetLocalSize(l, &s);
56 for (int i = 0; i != s; ++i)
57 a[i] = 0;
58 CHKERR VecRestoreArray(l, &a);
59 CHKERR VecGhostRestoreLocalForm(g, &l);
61 };
62 CHKERR zero_ghost_vec(f);
63
64 snes_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
65 auto cache_ptr = boost::make_shared<CacheTuple>();
67 cache_ptr);
68
69 auto set = [&](auto &fe) {
70 fe.snes = snes;
71 fe.snes_x = x;
72 fe.snes_f = f;
74 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
76 fe.cacheWeakPtr = cache_ptr;
77 };
78
79 auto unset = [&](auto &fe) {
80 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
81 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
82 fe.data_ctx = PetscData::CtxSetNone;
83 };
84
85 for (auto &bit : snes_ctx->preProcess_Rhs) {
86 bit->vecAssembleSwitch = boost::move(snes_ctx->vecAssembleSwitch);
87 set(*bit);
89 snes_ctx->problemName, *bit);
90 unset(*bit);
91 snes_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
92 }
93
94 for (auto &lit : snes_ctx->loops_to_do_Rhs) {
95 lit.second->vecAssembleSwitch = boost::move(snes_ctx->vecAssembleSwitch);
96 set(*lit.second);
98 snes_ctx->problemName, lit.first, *lit.second, nullptr, snes_ctx->bH,
99 cache_ptr);
100 unset(*lit.second);
101 if (snes_ctx->vErify) {
102 // Verify finite elements, check for not a number
103 CHKERR VecAssemblyBegin(f);
104 CHKERR VecAssemblyEnd(f);
105 MPI_Comm comm = PetscObjectComm((PetscObject)f);
106 PetscSynchronizedPrintf(comm, "SNES Verify f FE < %s >\n",
107 lit.first.c_str());
108 const Problem *prb_ptr;
109 CHKERR snes_ctx->mField.get_problem(snes_ctx->problemName, &prb_ptr);
110 CHKERR snes_ctx->mField.getInterface<Tools>()->checkVectorForNotANumber(
111 prb_ptr, ROW, f);
112 }
113
114 snes_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
115 }
116
117 for (auto &bit : snes_ctx->postProcess_Rhs) {
118 bit->vecAssembleSwitch = boost::move(snes_ctx->vecAssembleSwitch);
119 set(*bit);
121 snes_ctx->problemName, *bit);
122 unset(*bit);
123 snes_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
124 }
125
126 if (snes_ctx->vecAssembleSwitch) {
127 CHKERR VecGhostUpdateBegin(f, ADD_VALUES, SCATTER_REVERSE);
128 CHKERR VecGhostUpdateEnd(f, ADD_VALUES, SCATTER_REVERSE);
129 CHKERR VecAssemblyBegin(f);
130 CHKERR VecAssemblyEnd(f);
131 }
132 PetscLogEventEnd(snes_ctx->MOFEM_EVENT_SnesRhs, 0, 0, 0, 0);
134}
135
136PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx) {
137 SnesCtx *snes_ctx = (SnesCtx *)ctx;
138 // PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
140 PetscLogEventBegin(snes_ctx->MOFEM_EVENT_SnesMat, 0, 0, 0, 0);
141 if (snes_ctx->zeroPreCondMatrixB)
142 CHKERR MatZeroEntries(B);
143
144 snes_ctx->matAssembleSwitch = boost::movelib::make_unique<bool>(true);
145 auto cache_ptr = boost::make_shared<CacheTuple>();
147 cache_ptr);
148
149 auto set = [&](auto &fe) {
150 fe.snes = snes;
151 fe.snes_x = x;
152 fe.snes_A = A;
153 fe.snes_B = B;
155 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
157 fe.cacheWeakPtr = cache_ptr;
158 };
159
160 auto unset = [&](auto &fe) {
161 fe.snes_ctx = SnesMethod::CTX_SNESNONE;
162 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
163 fe.data_ctx = PetscData::CtxSetNone;
164 };
165
166
167 CHKERR VecGhostUpdateBegin(x, INSERT_VALUES, SCATTER_FORWARD);
168 CHKERR VecGhostUpdateEnd(x, INSERT_VALUES, SCATTER_FORWARD);
169 CHKERR snes_ctx->mField.getInterface<VecManager>()->setLocalGhostVector(
170 snes_ctx->problemName, COL, x, INSERT_VALUES, SCATTER_REVERSE);
171 for (auto &bit : snes_ctx->preProcess_Mat) {
172 bit->matAssembleSwitch = boost::move(snes_ctx->matAssembleSwitch);
173 set(*bit);
175 snes_ctx->problemName, *bit);
176 unset(*bit);
177 snes_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
178 }
179
180
181
182 for (auto &lit : snes_ctx->loops_to_do_Mat) {
183 lit.second->matAssembleSwitch = boost::move(snes_ctx->matAssembleSwitch);
184 set(*lit.second);
186 snes_ctx->problemName, lit.first, *(lit.second), nullptr, snes_ctx->bH,
187 cache_ptr);
188 unset(*lit.second);
189 snes_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
190 }
191
192 for (auto &bit : snes_ctx->postProcess_Mat) {
193 bit->matAssembleSwitch = boost::move(snes_ctx->matAssembleSwitch);
194 set(*bit);
196 snes_ctx->problemName, *bit);
197 unset(*bit);
198 snes_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
199 }
200
201 if (*snes_ctx->matAssembleSwitch) {
202 CHKERR MatAssemblyBegin(B, snes_ctx->typeOfAssembly);
203 CHKERR MatAssemblyEnd(B, snes_ctx->typeOfAssembly);
204 }
205 PetscLogEventEnd(snes_ctx->MOFEM_EVENT_SnesMat, 0, 0, 0, 0);
207}
208
209MoFEMErrorCode SnesMoFEMSetAssemblyType(SNES snes, MatAssemblyType type) {
210 SnesCtx *snes_ctx;
211 // PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
213 CHKERR SNESGetApplicationContext(snes, &snes_ctx);
214 snes_ctx->typeOfAssembly = type;
216}
217
219 SnesCtx *snes_ctx;
221 CHKERR SNESGetApplicationContext(snes, &snes_ctx);
222 snes_ctx->bH = bh;
224}
225
226MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm,
227 SnesCtx *snes_ctx) {
229 auto &m_field = snes_ctx->mField;
230 auto problem_ptr = m_field.get_problem(snes_ctx->problemName);
231 auto fields_ptr = m_field.get_fields();
232 auto dofs = problem_ptr->numeredRowDofsPtr;
233
234 std::vector<double> lnorms(fields_ptr->size(), 0),
235 norms(fields_ptr->size(), 0);
236
237 Vec res;
238 CHKERR SNESGetFunction(snes, &res, NULL, NULL);
239
240 const double *r;
241 CHKERR VecGetArrayRead(res, &r);
242 {
243 int f = 0;
244 for (auto fi : *fields_ptr) {
245 const auto lo_uid = FieldEntity::getLoBitNumberUId(fi->bitNumber);
246 const auto hi_uid = FieldEntity::getHiBitNumberUId(fi->bitNumber);
247 const auto hi = dofs->get<Unique_mi_tag>().upper_bound(hi_uid);
248 for (auto lo = dofs->get<Unique_mi_tag>().lower_bound(lo_uid); lo != hi;
249 ++lo) {
250 const DofIdx loc_idx = (*lo)->getPetscLocalDofIdx();
251 if (loc_idx >= 0 && loc_idx < problem_ptr->nbLocDofsRow) {
252 lnorms[f] += PetscRealPart(PetscSqr(r[loc_idx]));
253 }
254 }
255 ++f;
256 }
257 }
258 CHKERR VecRestoreArrayRead(res, &r);
259
260 MPIU_Allreduce(&*lnorms.begin(), &*norms.begin(), lnorms.size(), MPIU_REAL,
261 MPIU_SUM, PetscObjectComm((PetscObject)snes));
262
263 std::stringstream s;
264 int tl;
265 CHKERR PetscObjectGetTabLevel((PetscObject)snes, &tl);
266 for (auto t = 0; t != tl; ++t)
267 s << " ";
268 s << its << " Function norm " << boost::format("%14.12e") % (double)fgnorm
269 << " [";
270 {
271 int f = 0;
272 for (auto fi : *fields_ptr) {
273 if (f > 0)
274 s << ", ";
275 s << boost::format("%14.12e") % (double)PetscSqrtReal(norms[f]);
276 ++f;
277 }
278 s << "]";
279 }
280
281 MOFEM_LOG("SNES_WORLD", Sev::inform) << s.str();
282
284}
285
286} // namespace MoFEM
constexpr double a
@ COL
Definition: definitions.h:123
@ ROW
Definition: definitions.h:123
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:97
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
virtual MoFEMErrorCode problem_basic_method_postProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
virtual MoFEMErrorCode loop_finite_elements(const std::string problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
auto bit
set bit
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
int DofIdx
Index of DOF.
Definition: Types.hpp:18
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx)
This is MoFEM implementation for the left hand side (tangent matrix) evaluation in SNES solver.
Definition: SnesCtx.cpp:136
MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *snes_ctx)
Sens monitor printing residual field by field.
Definition: SnesCtx.cpp:226
PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx)
This is MoFEM implementation for the right hand side (residual vector) evaluation in SNES solver.
Definition: SnesCtx.cpp:27
MoFEMErrorCode SnesMoFEMSetAssemblyType(SNES snes, MatAssemblyType type)
Set assembly type at the end of SnesMat.
Definition: SnesCtx.cpp:209
MoFEMErrorCode SnesMoFEMSetBehavior(SNES snes, MoFEMTypes bh)
Set behavior if finite element in sequence does not exist.
Definition: SnesCtx.cpp:218
constexpr AssemblyType A
constexpr double t
plate stiffness
Definition: plate.cpp:59
constexpr double g
virtual MoFEMErrorCode cache_problem_entities(const std::string prb_name, CacheTupleWeakPtr cache_ptr)=0
Cache variables.
virtual MoFEMErrorCode problem_basic_method_preProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:37
static constexpr Switches CtxSetX
Definition: LoopMethods.hpp:39
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:35
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:36
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:38
keeps basic data about problem
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredRowDofsPtr
store DOFs on rows for this problem
Interface for nonlinear (SNES) solver.
Definition: SnesCtx.hpp:13
PetscLogEvent MOFEM_EVENT_SnesMat
Log events to assemble tangent matrix.
Definition: SnesCtx.hpp:128
BasicMethodsSequence preProcess_Rhs
Sequence of methods run before residual is assembled.
Definition: SnesCtx.hpp:39
BasicMethodsSequence postProcess_Mat
Definition: SnesCtx.hpp:36
MoFEMTypes bH
If set to MF_EXIST check if element exist, default MF_EXIST.
Definition: SnesCtx.hpp:20
BasicMethodsSequence postProcess_Rhs
Sequence of methods run after residual is assembled.
Definition: SnesCtx.hpp:41
MoFEMErrorCode copyLoops(const SnesCtx &snes_ctx)
Copy sequences from other SNES contex.
Definition: SnesCtx.cpp:5
MatAssemblyType typeOfAssembly
type of assembly at the end
Definition: SnesCtx.hpp:23
bool zeroPreCondMatrixB
Definition: SnesCtx.hpp:21
MoFEMErrorCode clearLoops()
Clear loops.
Definition: SnesCtx.cpp:16
BasicMethodsSequence preProcess_Mat
Definition: SnesCtx.hpp:34
bool vErify
If true verify vector.
Definition: SnesCtx.hpp:24
FEMethodsSequence loops_to_do_Mat
Definition: SnesCtx.hpp:30
std::string problemName
problem name
Definition: SnesCtx.hpp:18
FEMethodsSequence loops_to_do_Rhs
Definition: SnesCtx.hpp:32
boost::movelib::unique_ptr< bool > matAssembleSwitch
Definition: SnesCtx.hpp:126
MoFEM::Interface & mField
database Interface
Definition: SnesCtx.hpp:15
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Definition: SnesCtx.hpp:125
PetscLogEvent MOFEM_EVENT_SnesRhs
Log events to assemble residual.
Definition: SnesCtx.hpp:127
Auxiliary tools.
Definition: Tools.hpp:19
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:23