2 #include <BasicFiniteElements.hpp>
8 static char help[] =
"...\n\n";
19 , cOmm(m_field.get_comm())
20 , rAnk(m_field.get_comm_rank())
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>(
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());
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);
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);
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);
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);
77 MoFEMErrorCode extract_bd_ents(std::string ESSENTIAL, std::string NATURAL);
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);
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);
96 boost::shared_ptr<VectorDouble> &mass_ptr);
101 boost::shared_ptr<VolEle> &initial_ele,
131 boost::shared_ptr<PostProcVolumeOnRefinedMesh>
post_proc;
159 boost::shared_ptr<ForcesAndSourcesCore>
null;
164 CHKERR m_field.getInterface(simple_interface);
165 CHKERR simple_interface->getOptions();
166 CHKERR simple_interface->loadFile();
172 CHKERR simple_interface->addDomainField(
"u",
L2,
186 CHKERR simple_interface->setFieldOrder(
"u",
order - 1);
187 CHKERR simple_interface->setFieldOrder(
"f",
order);
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);
197 std::string natural) {
200 string name = it->getName();
201 if (name.compare(0, 14, natural) == 0) {
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);
218 block_id,
BLOCKSET, 3, surface,
true);
224 std::string mass_field, boost::shared_ptr<VectorDouble> &mass_ptr) {
226 vol_ele_slow_rhs->getOpPtrVector().push_back(
233 boost::shared_ptr<PreviousData> &data_v,
234 boost::shared_ptr<PreviousData> &data_w,
235 boost::shared_ptr<PreviousData> &data_s) {
238 vol_ele_slow_rhs->getOpPtrVector().push_back(
244 natural_bdry_ele_slow_rhs->getOpPtrVector().push_back(
250 std::string mass_field,
251 boost::shared_ptr<VectorDouble> &mass_ptr,
252 boost::shared_ptr<VectorDouble> &mass_dot_ptr) {
256 vol_ele_stiff_rhs->getOpPtrVector().push_back(
258 vol_ele_stiff_rhs->getOpPtrVector().push_back(
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(
272 vol_ele_stiff_rhs->getOpPtrVector().push_back(
278 std::string mass_field,
279 boost::shared_ptr<VectorDouble> &mass_ptr) {
281 vol_ele_stiff_lhs->getOpPtrVector().push_back(
291 vol_ele_stiff_lhs->getOpPtrVector().push_back(
294 vol_ele_stiff_lhs->getOpPtrVector().push_back(
297 vol_ele_stiff_lhs->getOpPtrVector().push_back(
300 vol_ele_stiff_lhs->getOpPtrVector().push_back(
303 vol_ele_stiff_lhs->getOpPtrVector().push_back(
306 vol_ele_stiff_lhs->getOpPtrVector().push_back(
309 vol_ele_stiff_lhs->getOpPtrVector().push_back(
312 vol_ele_stiff_lhs->getOpPtrVector().push_back(
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;
323 vol_ele_stiff_rhs->getRuleHook = vol_rule;
325 vol_ele_stiff_lhs->getRuleHook = vol_rule;
331 boost::shared_ptr<VolEle> &initial_ele,
double init_val) {
333 initial_ele->getOpPtrVector().push_back(
341 "SimpleProblem",
"f", essential_bdry_ents);
347 CHKERR TSSetType(ts, TSARKIMEX);
348 CHKERR TSARKIMEXSetType(ts, TSARKIMEXA2);
351 vol_ele_stiff_lhs,
null,
null);
354 vol_ele_stiff_rhs,
null,
null);
357 vol_ele_slow_rhs,
null,
null);
359 natural_bdry_ele_slow_rhs,
null,
null);
365 post_proc->addFieldValuesPostProc(mass_field);
372 monitor_ptr,
null,
null);
385 CHKERR TSSetDuration(ts, PETSC_DEFAULT, ftime);
386 CHKERR TSSetSolution(ts, X);
387 CHKERR TSSetFromOptions(ts);
399 CHKERR simple_interface->setUp();
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);
409 CHKERR push_slow_rhs(data_u, data_v, data_w, data_s);
414 CHKERR extract_bd_ents(
"ESSENTIAL",
"NATURAL");
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);
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(
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(
444 CHKERR push_stiff_rhs(data_u, data_v, data_w, data_s);
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);
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(
471 CHKERR set_integration_rule();
472 dm = simple_interface->getDM();
475 boost::shared_ptr<VolEle> initial_mass_ele(
new VolEle(m_field));
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);
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);
486 dm, simple_interface->getDomainFEName(), initial_mass_ele);
492 post_proc->generateReferenceElementMesh();
494 CHKERR post_proc_fields(
"U");
495 CHKERR post_proc_fields(
"V");
496 CHKERR post_proc_fields(
"W");
497 CHKERR post_proc_fields(
"S");
499 monitor_ptr = boost::shared_ptr<Monitor>(
500 new Monitor(cOmm, rAnk, dm, post_proc));
506 int main(
int argc,
char *argv[]) {
524 DMType dm_name =
"DMMOFEM";
MoFEM::FaceElementForcesAndSourcesCore FaceEle
#define CATCH_ERRORS
Catch errors.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
@ L2
field with C-1 continuity
@ HDIV
field with continuous normal traction
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
int main(int argc, char *argv[])
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
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
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
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
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.
implementation of Data Operators for Forces and Sources
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.
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
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
OpCalculateScalarFieldValuesDot OpCalculateScalarValuesDot
DeprecatedCoreInterface Interface
MoFEMErrorCode post_proc_fields(std::string mass_field)
MoFEMErrorCode extract_bd_ents(std::string ESSENTIAL, std::string NATURAL)
boost::shared_ptr< VectorDouble > divs_w_ptr
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)
boost::shared_ptr< VectorDouble > values_w_ptr
MoFEMErrorCode setup_system()
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
ElectroPhysioProblem(MoFEM::Core &core, const int order)
boost::shared_ptr< FaceEle > natural_bdry_ele_slow_rhs
Range essential_bdry_ents
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
Simple * simple_interface
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()
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
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.
VolumeElementForcesAndSourcesCore VolEle