v0.13.2
Loading...
Searching...
No Matches
electro_physio.cpp
Go to the documentation of this file.
1#include <stdlib.h>
3#include <EP_Operators.hpp>
4
5using namespace MoFEM;
6using namespace ElecPhys;
7
8static char help[] = "...\n\n";
9
10const double init_u = 0;
11const double init_v = 0;
12const double init_w = 1;
13const double init_s = 0.02155;
14
16public:
18 : m_field(core), order(order)
19 , cOmm(m_field.get_comm())
20 , rAnk(m_field.get_comm_rank())
21 {
22 vol_ele_slow_rhs = boost::shared_ptr<VolEle>(new VolEle(m_field));
24 boost::shared_ptr<FaceEle>(new FaceEle(m_field));
25 vol_ele_stiff_rhs = boost::shared_ptr<VolEle>(new VolEle(m_field));
26 vol_ele_stiff_lhs = boost::shared_ptr<VolEle>(new VolEle(m_field));
27 post_proc = boost::shared_ptr<PostProcVolumeOnRefinedMesh>(
29
30 data_u = boost::shared_ptr<PreviousData>(new PreviousData());
31 data_v = boost::shared_ptr<PreviousData>(new PreviousData());
32 data_w = boost::shared_ptr<PreviousData>(new PreviousData());
33 data_s = boost::shared_ptr<PreviousData>(new PreviousData());
34
36 boost::shared_ptr<MatrixDouble>(data_u, &data_u->flux_values);
38 boost::shared_ptr<VectorDouble>(data_u, &data_u->flux_divs);
40 boost::shared_ptr<VectorDouble>(data_u, &data_u->mass_values);
42 boost::shared_ptr<VectorDouble>(data_u, &data_u->mass_dots);
43
45 boost::shared_ptr<MatrixDouble>(data_v, &data_v->flux_values);
47 boost::shared_ptr<VectorDouble>(data_v, &data_v->flux_divs);
49 boost::shared_ptr<VectorDouble>(data_v, &data_v->mass_values);
51 boost::shared_ptr<VectorDouble>(data_v, &data_v->mass_dots);
52
53
55 boost::shared_ptr<MatrixDouble>(data_w, &data_w->flux_values);
57 boost::shared_ptr<VectorDouble>(data_w, &data_w->flux_divs);
59 boost::shared_ptr<VectorDouble>(data_w, &data_w->mass_values);
61 boost::shared_ptr<VectorDouble>(data_w, &data_w->mass_dots);
62
64 boost::shared_ptr<MatrixDouble>(data_s, &data_s->flux_values);
66 boost::shared_ptr<VectorDouble>(data_s, &data_s->flux_divs);
68 boost::shared_ptr<VectorDouble>(data_s, &data_s->mass_values);
70 boost::shared_ptr<VectorDouble>(data_s, &data_s->mass_dots);
71 }
73
74private:
77 MoFEMErrorCode extract_bd_ents(std::string ESSENTIAL, std::string NATURAL);
78
80
81 MoFEMErrorCode update_slow_rhs(std::string mass_field,
82 boost::shared_ptr<VectorDouble> &value_ptr);
83 MoFEMErrorCode push_slow_rhs(boost::shared_ptr<PreviousData> &data_u,
84 boost::shared_ptr<PreviousData> &data_v,
85 boost::shared_ptr<PreviousData> &data_w,
86 boost::shared_ptr<PreviousData> &data_s);
87
88 MoFEMErrorCode update_stiff_rhs(std::string mass_field,
89 boost::shared_ptr<VectorDouble> &mass_ptr,
90 boost::shared_ptr<VectorDouble> &mass_dot_ptr);
91 MoFEMErrorCode push_stiff_rhs(boost::shared_ptr<PreviousData> &data_u,
92 boost::shared_ptr<PreviousData> &data_v,
93 boost::shared_ptr<PreviousData> &data_w,
94 boost::shared_ptr<PreviousData> &data_s);
95 MoFEMErrorCode update_stiff_lhs(std::string mass_fiedl,
96 boost::shared_ptr<VectorDouble> &mass_ptr);
98
100 MoFEMErrorCode apply_IC(std::string mass_field, Range &surface,
101 boost::shared_ptr<VolEle> &initial_ele,
102 double init_val);
105 MoFEMErrorCode post_proc_fields(std::string mass_field);
108
113
116
117 Range inner_surface_u; // nb_species times
121
122 MPI_Comm cOmm;
123 const int rAnk;
124
125 int order;
126
127 boost::shared_ptr<VolEle> vol_ele_slow_rhs;
128 boost::shared_ptr<VolEle> vol_ele_stiff_rhs;
129 boost::shared_ptr<VolEle> vol_ele_stiff_lhs;
130 boost::shared_ptr<FaceEle> natural_bdry_ele_slow_rhs;
131 boost::shared_ptr<PostProcVolumeOnRefinedMesh> post_proc;
132 boost::shared_ptr<Monitor> monitor_ptr;
133
134 boost::shared_ptr<PreviousData> data_u; // nb_species times
135 boost::shared_ptr<PreviousData> data_v;
136 boost::shared_ptr<PreviousData> data_w;
137 boost::shared_ptr<PreviousData> data_s;
138
139 boost::shared_ptr<MatrixDouble> flux_u_ptr; // nb_species times
140 boost::shared_ptr<MatrixDouble> flux_v_ptr;
141 boost::shared_ptr<MatrixDouble> flux_w_ptr;
142 boost::shared_ptr<MatrixDouble> flux_s_ptr;
143
144 boost::shared_ptr<VectorDouble> divs_u_ptr; // nb_species times
145 boost::shared_ptr<VectorDouble> divs_v_ptr;
146 boost::shared_ptr<VectorDouble> divs_w_ptr;
147 boost::shared_ptr<VectorDouble> divs_s_ptr;
148
149 boost::shared_ptr<VectorDouble> values_u_ptr; // nb_species times
150 boost::shared_ptr<VectorDouble> values_v_ptr;
151 boost::shared_ptr<VectorDouble> values_w_ptr;
152 boost::shared_ptr<VectorDouble> values_s_ptr;
153
154 boost::shared_ptr<VectorDouble> dots_u_ptr;
155 boost::shared_ptr<VectorDouble> dots_v_ptr;
156 boost::shared_ptr<VectorDouble> dots_w_ptr;
157 boost::shared_ptr<VectorDouble> dots_s_ptr;
158
159 boost::shared_ptr<ForcesAndSourcesCore> null;
160};
161
168}
169
174
176 1);
177
179 1);
180
185
188
193
195}
197 std::string natural) {
200 string name = it->getName();
201 if (name.compare(0, 14, natural) == 0) {
202
203 CHKERR it->getMeshsetIdEntitiesByDimension(m_field.get_moab(), 2,
204 natural_bdry_ents, true);
205 } else if (name.compare(0, 14, essential) == 0) {
206 CHKERR it->getMeshsetIdEntitiesByDimension(m_field.get_moab(), 2,
207 essential_bdry_ents, true);
208 }
209 }
211}
213 Range &surface) {
216 BLOCKSET)) {
217 CHKERR m_field.getInterface<MeshsetsManager>()->getEntitiesByDimension(
218 block_id, BLOCKSET, 3, surface, true);
219 }
221}
222
224 std::string mass_field, boost::shared_ptr<VectorDouble> &mass_ptr) {
226 vol_ele_slow_rhs->getOpPtrVector().push_back(
227 new OpCalculateScalarFieldValues(mass_field, mass_ptr));
229}
230
232ElectroPhysioProblem::push_slow_rhs(boost::shared_ptr<PreviousData> &data_u,
233 boost::shared_ptr<PreviousData> &data_v,
234 boost::shared_ptr<PreviousData> &data_w,
235 boost::shared_ptr<PreviousData> &data_s) {
237
238 vol_ele_slow_rhs->getOpPtrVector().push_back(
240
241 // natural_bdry_ele_slow_rhs->getOpPtrVector().push_back(
242 // new OpAssembleNaturalBCRhsTau(natural_bdry_ents));
243
244 natural_bdry_ele_slow_rhs->getOpPtrVector().push_back(
247}
248
250 std::string mass_field,
251 boost::shared_ptr<VectorDouble> &mass_ptr,
252 boost::shared_ptr<VectorDouble> &mass_dot_ptr) {
253
255
256 vol_ele_stiff_rhs->getOpPtrVector().push_back(
257 new OpCalculateScalarFieldValues(mass_field, mass_ptr));
258 vol_ele_stiff_rhs->getOpPtrVector().push_back(
259 new OpCalculateScalarValuesDot(mass_field, mass_dot_ptr));
261}
262
264ElectroPhysioProblem::push_stiff_rhs(boost::shared_ptr<PreviousData> &data_u,
265 boost::shared_ptr<PreviousData> &data_v,
266 boost::shared_ptr<PreviousData> &data_w,
267 boost::shared_ptr<PreviousData> &data_s) {
269 vol_ele_stiff_rhs->getOpPtrVector().push_back(
271
272 vol_ele_stiff_rhs->getOpPtrVector().push_back(
275}
276
278 std::string mass_field,
279 boost::shared_ptr<VectorDouble> &mass_ptr) {
281 vol_ele_stiff_lhs->getOpPtrVector().push_back(
282 new OpCalculateScalarFieldValues(mass_field, mass_ptr));
283
284
286}
287
291 vol_ele_stiff_lhs->getOpPtrVector().push_back(
293
294 vol_ele_stiff_lhs->getOpPtrVector().push_back(
295 new OpAssembleLhsVV());
296
297 vol_ele_stiff_lhs->getOpPtrVector().push_back(
299
300 vol_ele_stiff_lhs->getOpPtrVector().push_back(
301 new OpAssembleLhsVTau());
302
303 vol_ele_stiff_lhs->getOpPtrVector().push_back(
305
306 vol_ele_stiff_lhs->getOpPtrVector().push_back(
308
309 vol_ele_stiff_lhs->getOpPtrVector().push_back(
311
312 vol_ele_stiff_lhs->getOpPtrVector().push_back(
315}
316
319 auto vol_rule = [](int, int, int p) -> int { return 2 * p + 1; };
320 vol_ele_slow_rhs->getRuleHook = vol_rule;
321 natural_bdry_ele_slow_rhs->getRuleHook = vol_rule;
322
323 vol_ele_stiff_rhs->getRuleHook = vol_rule;
324
325 vol_ele_stiff_lhs->getRuleHook = vol_rule;
327}
328
331 boost::shared_ptr<VolEle> &initial_ele, double init_val) {
333 initial_ele->getOpPtrVector().push_back(
334 new OpInitialMass(mass_field, surface, init_val));
336}
337
340 CHKERR m_field.getInterface<ProblemsManager>()->removeDofsOnEntities(
341 "SimpleProblem", "f", essential_bdry_ents);
342
344}
347 CHKERR TSSetType(ts, TSARKIMEX);
348 CHKERR TSARKIMEXSetType(ts, TSARKIMEXA2);
349
352
355
361}
362
365 post_proc->addFieldValuesPostProc(mass_field);
367}
368
374}
375
378 // Create solution vector
381 CHKERR DMoFEMMeshToLocalVector(dm, X, INSERT_VALUES, SCATTER_FORWARD);
382 // Solve problem
383 double ftime = 1;
384 CHKERR TSSetDM(ts, dm);
385 CHKERR TSSetDuration(ts, PETSC_DEFAULT, ftime);
386 CHKERR TSSetSolution(ts, X);
387 CHKERR TSSetFromOptions(ts);
388 CHKERR TSSolve(ts, X);
390}
391
394
395 CHKERR setup_system(); // only once
396 CHKERR add_fe(); // nb_species times
397
398
400
401
402
403
408
410
411 // natural_bdry_ele_slow_rhs->getOpPtrVector().push_back(
412 // new OpSetContrariantPiolaTransformOnEdge());
413
414 CHKERR extract_bd_ents("ESSENTIAL", "NATURAL"); // nb_species times
415
416
417
418
419
420
425
426 vol_ele_stiff_rhs->getOpPtrVector().push_back(
428 vol_ele_stiff_rhs->getOpPtrVector().push_back(
430 vol_ele_stiff_rhs->getOpPtrVector().push_back(
432 vol_ele_stiff_rhs->getOpPtrVector().push_back(
434
435 vol_ele_stiff_rhs->getOpPtrVector().push_back(
437 vol_ele_stiff_rhs->getOpPtrVector().push_back(
439 vol_ele_stiff_rhs->getOpPtrVector().push_back(
441 vol_ele_stiff_rhs->getOpPtrVector().push_back(
443
445
446
447
448
449
454
455 vol_ele_stiff_lhs->getOpPtrVector().push_back(
457 vol_ele_stiff_lhs->getOpPtrVector().push_back(
459 vol_ele_stiff_lhs->getOpPtrVector().push_back(
461 vol_ele_stiff_lhs->getOpPtrVector().push_back(
463
464
465
467
468
469
470
474
475 boost::shared_ptr<VolEle> initial_mass_ele(new VolEle(m_field));
476
477 CHKERR apply_IC("u", inner_surface_u, initial_mass_ele, init_u);
478 CHKERR apply_IC("U", inner_surface_u, initial_mass_ele, init_u); // nb_species times
479 CHKERR apply_IC("V", inner_surface_u, initial_mass_ele, init_v);
480 CHKERR apply_IC("W", inner_surface_u, initial_mass_ele, init_w);
481 CHKERR apply_IC("S", inner_surface_u, initial_mass_ele, init_s);
482
483
484
486 dm, simple_interface->getDomainFEName(), initial_mass_ele);
487
488 CHKERR apply_BC();
489
490
491 CHKERR loop_fe(); // only once
492 post_proc->generateReferenceElementMesh(); // only once
493
498
499 monitor_ptr = boost::shared_ptr<Monitor>(
500 new Monitor(cOmm, rAnk, dm, post_proc));
502 CHKERR solve();
504 }
505
506int main(int argc, char *argv[]) {
507 const char param_file[] = "param_file.petsc";
509 // for(int nsteps = 0; nsteps < 20; ++nsteps){
510 // if(((nsteps + 0) % 4) == 0){
511 // cout << "U" << nsteps << endl;
512 // } else if (((nsteps + 3) % 4) == 0) {
513 // cout << "V" << nsteps << endl;
514 // } else if (((nsteps + 2) % 4) == 0) {
515 // cout << "W" << nsteps << endl;
516 // } else if (((nsteps + 1) % 4) == 0) {
517 // cout << "S" << nsteps << endl;
518 // }
519 // }
520 try {
521 moab::Core mb_instance;
522 moab::Interface &moab = mb_instance;
523 MoFEM::Core core(moab);
524 DMType dm_name = "DMMOFEM";
525 CHKERR DMRegister_MoFEM(dm_name);
526
527 int order = 1;
528 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
529 ElectroPhysioProblem electro_physio_problem(core, order + 1);
530 CHKERR electro_physio_problem.run_analysis();
531 }
534 return 0;
535}
static Index< 'p', 3 > p
std::string param_file
int main()
Definition: adol-c_atom.cpp:46
MoFEM::FaceElementForcesAndSourcesCore FaceEle
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:66
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ BLOCKSET
Definition: definitions.h:148
#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
const double init_w
static char help[]
const double init_v
const double init_u
const double init_s
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set TS implicit function evaluation function
Definition: DMMoFEM.cpp:786
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMoFEM.cpp:509
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:47
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:572
PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMoFEM.cpp:839
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMoFEM.cpp:1153
bool checkMeshset(const int ms_id, const CubitBCType cubit_bc_type) const
check for CUBIT Id and CUBIT type
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
Definition: DMMoFEM.cpp:1042
PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS the right hand side function
Definition: DMMoFEM.cpp:868
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
OpCalculateScalarFieldValuesDot OpCalculateScalarValuesDot
auto createTS(MPI_Comm comm)
MoFEMErrorCode post_proc_fields(std::string mass_field)
MoFEMErrorCode solve()
MoFEMErrorCode extract_bd_ents(std::string ESSENTIAL, std::string NATURAL)
boost::shared_ptr< VectorDouble > divs_w_ptr
SmartPetscObj< DM > dm
MoFEMErrorCode output_result()
boost::shared_ptr< PreviousData > data_s
boost::shared_ptr< VectorDouble > values_u_ptr
MoFEMErrorCode run_analysis()
boost::shared_ptr< VectorDouble > divs_v_ptr
boost::shared_ptr< ForcesAndSourcesCore > null
MoFEMErrorCode apply_IC(std::string mass_field, Range &surface, boost::shared_ptr< VolEle > &initial_ele, double init_val)
boost::shared_ptr< VectorDouble > dots_s_ptr
MoFEMErrorCode update_stiff_rhs(std::string mass_field, boost::shared_ptr< VectorDouble > &mass_ptr, boost::shared_ptr< VectorDouble > &mass_dot_ptr)
MoFEMErrorCode add_fe()
boost::shared_ptr< VectorDouble > values_w_ptr
MoFEMErrorCode setup_system()
MoFEMErrorCode loop_fe()
MoFEMErrorCode extract_initial_ents(int block_id, Range &surface)
boost::shared_ptr< VectorDouble > dots_u_ptr
boost::shared_ptr< VectorDouble > divs_s_ptr
boost::shared_ptr< MatrixDouble > flux_u_ptr
SmartPetscObj< TS > ts
ElectroPhysioProblem(MoFEM::Core &core, const int order)
boost::shared_ptr< FaceEle > natural_bdry_ele_slow_rhs
boost::shared_ptr< VectorDouble > values_v_ptr
boost::shared_ptr< PreviousData > data_w
boost::shared_ptr< VolEle > vol_ele_stiff_lhs
MoFEMErrorCode update_stiff_lhs(std::string mass_fiedl, boost::shared_ptr< VectorDouble > &mass_ptr)
boost::shared_ptr< MatrixDouble > flux_w_ptr
MoFEMErrorCode apply_BC()
boost::shared_ptr< Monitor > monitor_ptr
MoFEMErrorCode push_slow_rhs(boost::shared_ptr< PreviousData > &data_u, boost::shared_ptr< PreviousData > &data_v, boost::shared_ptr< PreviousData > &data_w, boost::shared_ptr< PreviousData > &data_s)
MoFEMErrorCode set_integration_rule()
boost::shared_ptr< PreviousData > data_v
boost::shared_ptr< VectorDouble > dots_v_ptr
boost::shared_ptr< VectorDouble > divs_u_ptr
boost::shared_ptr< VolEle > vol_ele_slow_rhs
boost::shared_ptr< MatrixDouble > flux_s_ptr
boost::shared_ptr< VolEle > vol_ele_stiff_rhs
MoFEM::Interface & m_field
boost::shared_ptr< PostProcVolumeOnRefinedMesh > post_proc
MoFEMErrorCode push_stiff_rhs(boost::shared_ptr< PreviousData > &data_u, boost::shared_ptr< PreviousData > &data_v, boost::shared_ptr< PreviousData > &data_w, boost::shared_ptr< PreviousData > &data_s)
MoFEMErrorCode update_slow_rhs(std::string mass_field, boost::shared_ptr< VectorDouble > &value_ptr)
boost::shared_ptr< PreviousData > data_u
boost::shared_ptr< VectorDouble > dots_w_ptr
boost::shared_ptr< VectorDouble > values_s_ptr
boost::shared_ptr< MatrixDouble > flux_v_ptr
MoFEMErrorCode push_stiff_lhs()
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
Core (interface) class.
Definition: Core.hpp:82
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:72
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
Deprecated interface functions.
Interface for managing meshsets containing materials and boundary conditions.
Get vector field for H-div approximation.
Calculate divergence of vector field.
Get value at integration points for scalar field.
Problem manager is used to build and partition problems.
Simple interface for fast problem set-up.
Definition: Simple.hpp:27
const std::string getBoundaryFEName() const
Get the Boundary FE Name.
Definition: Simple.hpp:334
MoFEMErrorCode addDataField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on domain.
Definition: Simple.cpp:320
MoFEMErrorCode addDomainField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on domain.
Definition: Simple.cpp:264
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:180
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition: Simple.cpp:667
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:194
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:473
MoFEMErrorCode addBoundaryField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on boundary.
Definition: Simple.cpp:282
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:611
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition: Simple.hpp:327
intrusive_ptr for managing petsc objects
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Post processing.
VolumeElementForcesAndSourcesCore VolEle