v0.13.0
electro_physio.cpp
Go to the documentation of this file.
1 #include <stdlib.h>
2 #include <BasicFiniteElements.hpp>
3 #include <EP_Operators.hpp>
4 
5 using namespace MoFEM;
6 using namespace ElecPhys;
7 
8 static char help[] = "...\n\n";
9 
10 const double init_u = 0;
11 const double init_v = 0;
12 const double init_w = 1;
13 const double init_s = 0.02155;
14 
16 public:
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));
23  natural_bdry_ele_slow_rhs =
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>(
28  new PostProcVolumeOnRefinedMesh(m_field));
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 
35  flux_u_ptr =
36  boost::shared_ptr<MatrixDouble>(data_u, &data_u->flux_values);
37  divs_u_ptr =
38  boost::shared_ptr<VectorDouble>(data_u, &data_u->flux_divs);
39  values_u_ptr =
40  boost::shared_ptr<VectorDouble>(data_u, &data_u->mass_values);
41  dots_u_ptr =
42  boost::shared_ptr<VectorDouble>(data_u, &data_u->mass_dots);
43 
44  flux_v_ptr =
45  boost::shared_ptr<MatrixDouble>(data_v, &data_v->flux_values);
46  divs_v_ptr =
47  boost::shared_ptr<VectorDouble>(data_v, &data_v->flux_divs);
48  values_v_ptr =
49  boost::shared_ptr<VectorDouble>(data_v, &data_v->mass_values);
50  dots_v_ptr =
51  boost::shared_ptr<VectorDouble>(data_v, &data_v->mass_dots);
52 
53 
54  flux_w_ptr =
55  boost::shared_ptr<MatrixDouble>(data_w, &data_w->flux_values);
56  divs_w_ptr =
57  boost::shared_ptr<VectorDouble>(data_w, &data_w->flux_divs);
58  values_w_ptr =
59  boost::shared_ptr<VectorDouble>(data_w, &data_w->mass_values);
60  dots_w_ptr =
61  boost::shared_ptr<VectorDouble>(data_w, &data_w->mass_dots);
62 
63  flux_s_ptr =
64  boost::shared_ptr<MatrixDouble>(data_s, &data_s->flux_values);
65  divs_s_ptr =
66  boost::shared_ptr<VectorDouble>(data_s, &data_s->flux_divs);
67  values_s_ptr =
68  boost::shared_ptr<VectorDouble>(data_s, &data_s->mass_values);
69  dots_s_ptr =
70  boost::shared_ptr<VectorDouble>(data_s, &data_s->mass_dots);
71  }
72  MoFEMErrorCode run_analysis();
73 
74 private:
75  MoFEMErrorCode setup_system();
76  MoFEMErrorCode add_fe();
77  MoFEMErrorCode extract_bd_ents(std::string ESSENTIAL, std::string NATURAL);
78 
79  MoFEMErrorCode extract_initial_ents(int block_id, Range &surface);
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);
97  MoFEMErrorCode push_stiff_lhs();
98 
99  MoFEMErrorCode set_integration_rule();
100  MoFEMErrorCode apply_IC(std::string mass_field, Range &surface,
101  boost::shared_ptr<VolEle> &initial_ele,
102  double init_val);
103  MoFEMErrorCode apply_BC();
104  MoFEMErrorCode loop_fe();
105  MoFEMErrorCode post_proc_fields(std::string mass_field);
106  MoFEMErrorCode output_result();
107  MoFEMErrorCode solve();
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 
164  CHKERR m_field.getInterface(simple_interface);
165  CHKERR simple_interface->getOptions();
166  CHKERR simple_interface->loadFile();
168 }
169 
172  CHKERR simple_interface->addDomainField("u", L2,
174 
175  CHKERR simple_interface->addDomainField("f", HDIV, DEMKOWICZ_JACOBI_BASE,
176  1);
177 
178  CHKERR simple_interface->addBoundaryField("f", HDIV, DEMKOWICZ_JACOBI_BASE,
179  1);
180 
181  CHKERR simple_interface->addDataField("U", L2, AINSWORTH_LEGENDRE_BASE, 1);
182  CHKERR simple_interface->addDataField("V", L2, AINSWORTH_LEGENDRE_BASE, 1);
183  CHKERR simple_interface->addDataField("W", L2, AINSWORTH_LEGENDRE_BASE, 1);
184  CHKERR simple_interface->addDataField("S", L2, AINSWORTH_LEGENDRE_BASE, 1);
185 
186  CHKERR simple_interface->setFieldOrder("u", order - 1);
187  CHKERR simple_interface->setFieldOrder("f", order);
188 
189  CHKERR simple_interface->setFieldOrder("U", order - 1);
190  CHKERR simple_interface->setFieldOrder("V", order - 1);
191  CHKERR simple_interface->setFieldOrder("W", order - 1);
192  CHKERR simple_interface->setFieldOrder("S", order - 1);
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) {
215  if (m_field.getInterface<MeshsetsManager>()->checkMeshset(block_id,
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 
232 ElectroPhysioProblem::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(
239  new OpAssembleSlowRhsV(data_u, data_v, data_w, data_s));
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(
245  new OpEssentialBC(essential_bdry_ents));
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 
264 ElectroPhysioProblem::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(
270  new OpAssembleStiffRhsTau<3>(data_u, data_v, data_w, data_s));
271 
272  vol_ele_stiff_rhs->getOpPtrVector().push_back(
273  new OpAssembleStiffRhsV<3>(data_u, data_v, data_w, data_s));
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(
292  new OpAssembleLhsTauTau<3>());
293 
294  vol_ele_stiff_lhs->getOpPtrVector().push_back(
295  new OpAssembleLhsVV());
296 
297  vol_ele_stiff_lhs->getOpPtrVector().push_back(
298  new OpAssembleLhsTauV<3>());
299 
300  vol_ele_stiff_lhs->getOpPtrVector().push_back(
301  new OpAssembleLhsVTau());
302 
303  vol_ele_stiff_lhs->getOpPtrVector().push_back(
304  new OpDamp_dofs_to_field_data("U"));
305 
306  vol_ele_stiff_lhs->getOpPtrVector().push_back(
307  new OpDamp_dofs_to_field_data("V"));
308 
309  vol_ele_stiff_lhs->getOpPtrVector().push_back(
310  new OpDamp_dofs_to_field_data("W"));
311 
312  vol_ele_stiff_lhs->getOpPtrVector().push_back(
313  new OpDamp_dofs_to_field_data("S"));
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 
330 ElectroPhysioProblem::apply_IC(std::string mass_field, Range &surface,
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 
350  CHKERR DMMoFEMTSSetIJacobian(dm, simple_interface->getDomainFEName(),
351  vol_ele_stiff_lhs, null, null);
352 
353  CHKERR DMMoFEMTSSetIFunction(dm, simple_interface->getDomainFEName(),
354  vol_ele_stiff_rhs, null, null);
355 
356  CHKERR DMMoFEMTSSetRHSFunction(dm, simple_interface->getDomainFEName(),
357  vol_ele_slow_rhs, null, null);
358  CHKERR DMMoFEMTSSetRHSFunction(dm, simple_interface->getBoundaryFEName(),
359  natural_bdry_ele_slow_rhs, null, null);
361 }
362 
365  post_proc->addFieldValuesPostProc(mass_field);
367 }
368 
371  CHKERR DMMoFEMTSSetMonitor(dm, ts, simple_interface->getDomainFEName(),
372  monitor_ptr, null, null);
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 
399  CHKERR simple_interface->setUp();
400 
401 
402 
403 
404  CHKERR update_slow_rhs("U", values_u_ptr);
405  CHKERR update_slow_rhs("V", values_v_ptr);
406  CHKERR update_slow_rhs("W", values_w_ptr);
407  CHKERR update_slow_rhs("S", values_s_ptr);
408 
409  CHKERR push_slow_rhs(data_u, data_v, data_w, data_s);
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 
421  CHKERR update_stiff_rhs("U", values_u_ptr, dots_u_ptr);
422  CHKERR update_stiff_rhs("V", values_v_ptr, dots_v_ptr);
423  CHKERR update_stiff_rhs("W", values_w_ptr, dots_w_ptr);
424  CHKERR update_stiff_rhs("S", values_s_ptr, dots_s_ptr);
425 
426  vol_ele_stiff_rhs->getOpPtrVector().push_back(
427  new OpCalculateHVecVectorField<3>("f", flux_u_ptr));
428  vol_ele_stiff_rhs->getOpPtrVector().push_back(
429  new OpCalculateHVecVectorField<3>("f", flux_v_ptr));
430  vol_ele_stiff_rhs->getOpPtrVector().push_back(
431  new OpCalculateHVecVectorField<3>("f", flux_w_ptr));
432  vol_ele_stiff_rhs->getOpPtrVector().push_back(
433  new OpCalculateHVecVectorField<3>("f", flux_s_ptr));
434 
435  vol_ele_stiff_rhs->getOpPtrVector().push_back(
436  new OpCalculateHdivVectorDivergence<3, 3>("f", divs_u_ptr));
437  vol_ele_stiff_rhs->getOpPtrVector().push_back(
438  new OpCalculateHdivVectorDivergence<3, 3>("f", divs_v_ptr));
439  vol_ele_stiff_rhs->getOpPtrVector().push_back(
440  new OpCalculateHdivVectorDivergence<3, 3>("f", divs_w_ptr));
441  vol_ele_stiff_rhs->getOpPtrVector().push_back(
442  new OpCalculateHdivVectorDivergence<3, 3>("f", divs_s_ptr));
443 
444  CHKERR push_stiff_rhs(data_u, data_v, data_w, data_s);
445 
446 
447 
448 
449 
450  CHKERR update_stiff_lhs("U", values_u_ptr);
451  CHKERR update_stiff_lhs("V", values_v_ptr);
452  CHKERR update_stiff_lhs("W", values_w_ptr);
453  CHKERR update_stiff_lhs("S", values_s_ptr);
454 
455  vol_ele_stiff_lhs->getOpPtrVector().push_back(
456  new OpCalculateHVecVectorField<3>("f", flux_u_ptr));
457  vol_ele_stiff_lhs->getOpPtrVector().push_back(
458  new OpCalculateHVecVectorField<3>("f", flux_v_ptr));
459  vol_ele_stiff_lhs->getOpPtrVector().push_back(
460  new OpCalculateHVecVectorField<3>("f", flux_w_ptr));
461  vol_ele_stiff_lhs->getOpPtrVector().push_back(
462  new OpCalculateHVecVectorField<3>("f", flux_s_ptr));
463 
464 
465 
466  CHKERR push_stiff_lhs();
467 
468 
469 
470 
471  CHKERR set_integration_rule();
472  dm = simple_interface->getDM();
473  ts = createTS(m_field.get_comm());
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 
494  CHKERR post_proc_fields("U");
495  CHKERR post_proc_fields("V");
496  CHKERR post_proc_fields("W");
497  CHKERR post_proc_fields("S");
498 
499  monitor_ptr = boost::shared_ptr<Monitor>(
500  new Monitor(cOmm, rAnk, dm, post_proc));
501  CHKERR output_result();
502  CHKERR solve();
504  }
505 
506 int main(int argc, char *argv[]) {
507  const char param_file[] = "param_file.petsc";
508  MoFEM::Core::Initialize(&argc, &argv, param_file, help);
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  }
532  CATCH_ERRORS;
534  return 0;
535 }
static Index< 'p', 3 > p
std::string param_file
MoFEM::FaceElementForcesAndSourcesCore FaceEle
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:79
@ L2
field with C-1 continuity
Definition: definitions.h:101
@ HDIV
field with continuous normal traction
Definition: definitions.h:100
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ BLOCKSET
Definition: definitions.h:161
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
const double init_w
int main(int argc, char *argv[])
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: DMMMoFEM.cpp:758
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:481
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:544
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: DMMMoFEM.cpp:811
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMMoFEM.cpp:1105
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:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
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: DMMMoFEM.cpp:994
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: DMMMoFEM.cpp:840
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
auto createTS
OpCalculateScalarFieldValuesDot OpCalculateScalarValuesDot
CoreTmp< 0 > Core
Definition: Core.hpp:1096
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
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
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()
Core (interface) class.
Definition: Core.hpp:92
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:85
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:125
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:33
Post processing.
VolumeElementForcesAndSourcesCore VolEle