v0.6.20
Classes | Functions | Variables
thermal_unsteady.cpp File Reference

Example of thermal unsteady analyze. More...

#include <BasicFiniteElements.hpp>
#include <boost/program_options.hpp>

Go to the source code of this file.

Classes

struct  BlockOptionData
 
struct  MonitorPostProc
 

Functions

int main (int argc, char *argv[])
 

Variables

static char help []
 

Detailed Description

Example of thermal unsteady analyze.

TODO:

Todo:
Make it work in distributed meshes with multigird solver. At the moment it is not working efficient as can.

Definition in file thermal_unsteady.cpp.

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)

Definition at line 119 of file thermal_unsteady.cpp.

119  {
120 
121  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
122 
123  try {
124 
125  PetscBool flg = PETSC_TRUE;
126  char mesh_file_name[255];
127  CHKERR PetscOptionsGetString(PETSC_NULL, PETSC_NULL, "-my_file",
128  mesh_file_name, 255, &flg);
129  if(flg != PETSC_TRUE) {
130  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND,
131  "*** ERROR -my_file (MESH FILE NEEDED)");
132  }
133  const char *option;
134  option = "";
135 
136  char time_data_file_for_ground_surface[255];
137  PetscBool ground_temperature_analysis = PETSC_FALSE;
138  CHKERR PetscOptionsGetString(PETSC_NULL,PETSC_NULL,"-my_ground_analysis_data",
139  time_data_file_for_ground_surface,255,&ground_temperature_analysis);
140  if(ground_temperature_analysis) {
141 #ifndef WITH_ADOL_C
142  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_INSTALLED,
143  "*** ERROR to do ground thermal analysis MoFEM need to be compiled "
144  "with ADOL-C");
145 #endif // WITH_ADOL_C
146  }
147 
148  //create MoAB database
149  moab::Core mb_instance;
150  moab::Interface& moab = mb_instance;
151  ParallelComm* pcomm = ParallelComm::get_pcomm(&moab,MYPCOMM_INDEX);
152  if(pcomm == NULL) pcomm = new ParallelComm(&moab,PETSC_COMM_WORLD);
153  CHKERR moab.load_file(mesh_file_name, 0, option);
154  //create MoFEM database
155  MoFEM::Core core(moab);
156  MoFEM::Interface& m_field = core;
157 
158  DMType dm_name = "DMTHERMAL";
159  CHKERR DMRegister_MoFEM(dm_name);
160  // create dm instance
161  DM dm;
162  CHKERR DMCreate(PETSC_COMM_WORLD, &dm);
163  CHKERR DMSetType(dm, dm_name);
164 
165  //set entities bit level (this allow to set refinement levels for h-adaptivity)
166  //only one level is used in this example
167  BitRefLevel bit_level0;
168  bit_level0.set(0);
169  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(0, 3,
170  bit_level0);
171 
172  //Fields H1 space rank 1
173  CHKERR m_field.add_field("TEMP", H1, AINSWORTH_LEGENDRE_BASE, 1,
174  MB_TAG_SPARSE, MF_ZERO);
175  CHKERR m_field.add_field("TEMP_RATE", H1, AINSWORTH_LEGENDRE_BASE, 1,
176  MB_TAG_SPARSE, MF_ZERO);
177 
178  //Add field H1 space rank 3 to approximate geometry using hierarchical basis
179  //For 10 node tets, before use, geometry is projected on that field (see below)
180  CHKERR m_field.add_field(
181  "MESH_NODE_POSITIONS",H1,AINSWORTH_LEGENDRE_BASE,3,MB_TAG_SPARSE,MF_ZERO
182  );
183 
184  //meshset consisting all entities in mesh
185  EntityHandle root_set = moab.get_root_set();
186  //add entities to field (root_mesh, i.e. on all mesh etities fields are approx.)
187  CHKERR m_field.add_ents_to_field_by_type(root_set,MBTET,"TEMP");
188  CHKERR m_field.add_ents_to_field_by_type(root_set,MBTET,"TEMP_RATE");
189 
190  int order;
191  CHKERR PetscOptionsGetInt(PETSC_NULL, PETSC_NULL, "-my_order", &order, &flg);
192 
193  if (flg != PETSC_TRUE) {
194  order = 1;
195  }
196  // set app. order
197  // see Hierarchic Finite Element Bases on Unstructured Tetrahedral Meshes
198  // (Mark Ainsworth & Joe Coyle) for simplicity of example to all entities is
199  // applied the same order
200  CHKERR m_field.set_field_order(root_set,MBTET,"TEMP",order);
201  CHKERR m_field.set_field_order(root_set,MBTRI,"TEMP",order);
202  CHKERR m_field.set_field_order(root_set,MBEDGE,"TEMP",order);
203  CHKERR m_field.set_field_order(root_set,MBVERTEX,"TEMP",1);
204 
205  CHKERR m_field.set_field_order(root_set,MBTET,"TEMP_RATE",order);
206  CHKERR m_field.set_field_order(root_set,MBTRI,"TEMP_RATE",order);
207  CHKERR m_field.set_field_order(root_set,MBEDGE,"TEMP_RATE",order);
208  CHKERR m_field.set_field_order(root_set,MBVERTEX,"TEMP_RATE",1);
209 
210  //geometry approximation is set to 2nd oreder
211  CHKERR m_field.add_ents_to_field_by_type(root_set, MBTET,
212  "MESH_NODE_POSITIONS");
213  CHKERR m_field.set_field_order(0,MBTET,"MESH_NODE_POSITIONS",2);
214  CHKERR m_field.set_field_order(0,MBTRI,"MESH_NODE_POSITIONS",2);
215  CHKERR m_field.set_field_order(0,MBEDGE,"MESH_NODE_POSITIONS",2);
216  CHKERR m_field.set_field_order(0,MBVERTEX,"MESH_NODE_POSITIONS",1);
217 
218  // configure blocks by parsing config file
219  // it allow to set approximation order for each block independently
220  PetscBool block_config;
221  char block_config_file[255];
222  CHKERR PetscOptionsGetString(PETSC_NULL, PETSC_NULL, "-my_block_config",
223  block_config_file, 255, &block_config);
224  std::map<int,BlockOptionData> block_data;
225  bool solar_radiation = false;
226  if (block_config) {
227  try {
228  ifstream ini_file(block_config_file);
229  // std::cerr << block_config_file << std::endl;
230  po::variables_map vm;
231  po::options_description config_file_options;
233 
234  std::ostringstream str_order;
235  str_order << "block_" << it->getMeshsetId() << ".temperature_order";
236  config_file_options.add_options()(
237  str_order.str().c_str(),
238  po::value<int>(&block_data[it->getMeshsetId()].oRder)
239  ->default_value(order));
240 
241  std::ostringstream str_cond;
242  str_cond << "block_" << it->getMeshsetId() << ".heat_conductivity";
243  config_file_options.add_options()(
244  str_cond.str().c_str(),
245  po::value<double>(&block_data[it->getMeshsetId()].cOnductivity)
246  ->default_value(-1));
247 
248  std::ostringstream str_capa;
249  str_capa << "block_" << it->getMeshsetId() << ".heat_capacity";
250  config_file_options.add_options()(
251  str_capa.str().c_str(),
252  po::value<double>(&block_data[it->getMeshsetId()].cApacity)
253  ->default_value(-1));
254 
255  std::ostringstream str_init_temp;
256  str_init_temp << "block_" << it->getMeshsetId()
257  << ".initail_temperature";
258  config_file_options.add_options()(
259  str_init_temp.str().c_str(),
260  po::value<double>(&block_data[it->getMeshsetId()].initTemp)
261  ->default_value(0));
262  }
263  config_file_options.add_options()(
264  "climate_model.solar_radiation",
265  po::value<bool>(&solar_radiation)->default_value(false));
266 
267  po::parsed_options parsed =
268  parse_config_file(ini_file, config_file_options, true);
269  store(parsed,vm);
270  po::notify(vm);
271 
273  if (block_data[it->getMeshsetId()].oRder == -1)
274  continue;
275  if (block_data[it->getMeshsetId()].oRder == order)
276  continue;
277  PetscPrintf(PETSC_COMM_WORLD, "Set block %d oRder to %d\n",
278  it->getMeshsetId(), block_data[it->getMeshsetId()].oRder);
279  Range block_ents;
280  CHKERR moab.get_entities_by_handle(it->meshset, block_ents, true);
281  CHKERRG(rval);
282  Range ents_to_set_order;
283  CHKERR moab.get_adjacencies(block_ents, 3, false, ents_to_set_order,
284  moab::Interface::UNION);
285  ents_to_set_order = ents_to_set_order.subset_by_type(MBTET);
286  CHKERR moab.get_adjacencies(block_ents, 2, false, ents_to_set_order,
287  moab::Interface::UNION);
288  CHKERR moab.get_adjacencies(block_ents, 1, false, ents_to_set_order,
289  moab::Interface::UNION);
290  CHKERR m_field.set_field_order(ents_to_set_order, "TEMP",
291  block_data[it->getMeshsetId()].oRder);
292  CHKERR m_field.set_field_order(ents_to_set_order, "TEMP_RATE",
293  block_data[it->getMeshsetId()].oRder);
294  }
295  std::vector<std::string> additional_parameters;
296  additional_parameters =
297  collect_unrecognized(parsed.options, po::include_positional);
298  for (std::vector<std::string>::iterator vit =
299  additional_parameters.begin();
300  vit != additional_parameters.end(); vit++) {
301  CHKERR PetscPrintf(PETSC_COMM_WORLD,
302  "** WARRING Unrecognised option %s\n", vit->c_str());
303  }
304 
305  } catch (const std::exception& ex) {
306  std::ostringstream ss;
307  ss << ex.what() << std::endl;
308  SETERRQ(PETSC_COMM_SELF,MOFEM_STD_EXCEPTION_THROW,ss.str().c_str());
309  }
310  }
311 
312  // this default class to calculate thermal elements
313  ThermalElement thermal_elements(m_field);
314  CHKERR thermal_elements.addThermalElements("TEMP");
315  CHKERR thermal_elements.addThermalFluxElement("TEMP");
316  CHKERR thermal_elements.addThermalConvectionElement("TEMP");
317  CHKERR thermal_elements.addThermalRadiationElement("TEMP");
318  // add rate of temperature to data field of finite element
319  CHKERR m_field.modify_finite_element_add_field_data("THERMAL_FE",
320  "TEMP_RATE");
321  // and temperature element default element operators at integration (gauss)
322  // points
323  CHKERR thermal_elements.setTimeSteppingProblem("TEMP", "TEMP_RATE");
324 
325  // set block material data from option file
326  std::map<int, ThermalElement::BlockData>::iterator mit;
327  mit = thermal_elements.setOfBlocks.begin();
328  for (; mit != thermal_elements.setOfBlocks.end(); mit++) {
329  // std::cerr << mit->first << std::endl;
330  // std::cerr << block_data[mit->first].cOnductivity << " " <<
331  // block_data[mit->first].cApacity << std::endl;
332  if (block_data[mit->first].cOnductivity != -1) {
333  PetscPrintf(PETSC_COMM_WORLD, "Set block %d heat conductivity to %3.2e\n",
334  mit->first, block_data[mit->first].cOnductivity);
335  for (int dd = 0; dd < 3; dd++) {
336  mit->second.cOnductivity_mat(dd, dd) =
337  block_data[mit->first].cOnductivity;
338  }
339  }
340  if (block_data[mit->first].cApacity != -1) {
341  PetscPrintf(PETSC_COMM_WORLD, "Set block %d heat capacity to %3.2e\n",
342  mit->first, block_data[mit->first].cApacity);
343  mit->second.cApacity = block_data[mit->first].cApacity;
344  }
345  }
346 
347 #ifdef __GROUNDSURFACETEMERATURE_HPP
348  GroundSurfaceTemerature ground_surface(m_field);
349  CrudeClimateModel time_data(time_data_file_for_ground_surface);
350  GroundSurfaceTemerature::PreProcess exectuteGenericClimateModel(&time_data);
351  if (ground_temperature_analysis) {
352  CHKERR ground_surface.addSurfaces("TEMP");
353  CHKERR ground_surface.setOperators(&time_data, "TEMP");
354  }
355 #endif //__GROUNDSURFACETEMERATURE_HPP
356 
357  //build database, i.e. declare dofs, elements and adjacencies
358 
359  // build field
360  CHKERR m_field.build_fields();
361  // priject 10 node tet approximation of geometry on hierarchical basis
362  Projection10NodeCoordsOnField ent_method_material(m_field,
363  "MESH_NODE_POSITIONS");
364  CHKERR m_field.loop_dofs("MESH_NODE_POSITIONS", ent_method_material);
366  if (block_data[it->getMeshsetId()].initTemp != 0) {
367  Range block_ents;
368  CHKERR moab.get_entities_by_handle(it->meshset, block_ents, true);
369  CHKERRG(rval);
370  Range vertices;
371  CHKERR moab.get_connectivity(block_ents, vertices, true);
372  CHKERRG(rval);
373  CHKERR m_field.getInterface<FieldBlas>()->setField(
374  block_data[it->getMeshsetId()].initTemp, MBVERTEX, vertices, "TEMP");
375  }
376  }
377 
378  // build finite elemnts
379  CHKERR m_field.build_finite_elements();
380  // build adjacencies
381  CHKERR m_field.build_adjacencies(bit_level0);
382 
383  // delete old temperature recorded series
384  SeriesRecorder *recorder_ptr;
385  CHKERR m_field.getInterface(recorder_ptr);
386  if (recorder_ptr->check_series("THEMP_SERIES")) {
387  /*for(_IT_SERIES_STEPS_BY_NAME_FOR_LOOP_(recorder_ptr,"THEMP_SERIES",sit)) {
388  CHKERR
389  recorder_ptr->load_series_data("THEMP_SERIES",sit->get_step_number());
390  }*/
391  CHKERR recorder_ptr->delete_recorder_series("THEMP_SERIES");
392  }
393 
394  // set dm data structure which created mofem data structures
395  CHKERR DMMoFEMCreateMoFEM(dm, &m_field, dm_name, bit_level0);
396  CHKERR DMSetFromOptions(dm);
397  // add elements to dm
398  CHKERR DMMoFEMAddElement(dm, "THERMAL_FE");
399  CHKERR DMMoFEMAddElement(dm, "THERMAL_FLUX_FE");
400  CHKERR DMMoFEMAddElement(dm, "THERMAL_CONVECTION_FE");
401  CHKERR DMMoFEMAddElement(dm, "THERMAL_RADIATION_FE");
402 #ifdef __GROUNDSURFACETEMERATURE_HPP
403  if (ground_temperature_analysis) {
404  CHKERR DMMoFEMAddElement(dm, "GROUND_SURFACE_FE");
405  }
406 #endif //__GROUNDSURFACETEMERATURE_HPP
407 
408  CHKERR DMSetUp(dm);
409 
410  // create matrices
411  Vec T, F;
413  CHKERR VecDuplicate(T, &F);
414  Mat A;
416 
417  DirichletTemperatureBc dirichlet_bc(m_field, "TEMP", A, T, F);
418  ThermalElement::UpdateAndControl update_velocities(m_field, "TEMP",
419  "TEMP_RATE");
420  ThermalElement::TimeSeriesMonitor monitor(m_field, "THEMP_SERIES", "TEMP");
421  MonitorPostProc post_proc(m_field);
422 
423  // Initialize data with values save of on the field
424  CHKERR VecZeroEntries(T);
425  CHKERR DMoFEMMeshToLocalVector(dm, T, INSERT_VALUES, SCATTER_FORWARD);
426  CHKERR DMoFEMPreProcessFiniteElements(dm, &dirichlet_bc);
427  CHKERR DMoFEMMeshToGlobalVector(dm, T, INSERT_VALUES, SCATTER_REVERSE);
428 
429  // preprocess
430  CHKERR DMMoFEMTSSetIFunction(dm, DM_NO_ELEMENT, NULL, &update_velocities,
431  NULL);
432  CHKERR DMMoFEMTSSetIFunction(dm, DM_NO_ELEMENT, NULL, &dirichlet_bc, NULL);
433  CHKERR DMMoFEMTSSetIJacobian(dm, DM_NO_ELEMENT, NULL, &dirichlet_bc, NULL);
434 #ifdef __GROUNDSURFACETEMERATURE_HPP
436  &exectuteGenericClimateModel, NULL);
437  { // add preprocessor, calculating angle on which sun ray on the surface
438  if (solar_radiation) {
439  boost::ptr_vector<
441  hi_it;
442  it = ground_surface.preProcessShade.begin();
443  hi_it = ground_surface.preProcessShade.end();
444  for (; it != hi_it; it++) {
445  CHKERR DMMoFEMTSSetIFunction(dm, DM_NO_ELEMENT, NULL, &*it, NULL);
446  }
447  }
448  }
449 #endif //__GROUNDSURFACETEMERATURE_HPP
450 
451  // loops rhs
452  CHKERR DMMoFEMTSSetIFunction(dm, "THERMAL_FE", &thermal_elements.feRhs, NULL,
453  NULL);
454  CHKERR DMMoFEMTSSetIFunction(dm, "THERMAL_FLUX_FE", &thermal_elements.feFlux,
455  NULL, NULL);
456  CHKERR DMMoFEMTSSetIFunction(dm, "THERMAL_CONVECTION_FE",
457  &thermal_elements.feConvectionRhs, NULL, NULL);
458  CHKERR DMMoFEMTSSetIFunction(dm, "THERMAL_RADIATION_FE",
459  &thermal_elements.feRadiationRhs, NULL, NULL);
460 #ifdef __GROUNDSURFACETEMERATURE_HPP
461  if (ground_temperature_analysis) {
462  CHKERR DMMoFEMTSSetIFunction(dm, "GROUND_SURFACE_FE",
463  &ground_surface.getFeGroundSurfaceRhs(), NULL,
464  NULL);
465  }
466 #endif //__GROUNDSURFACETEMERATURE_HPP
467 
468  // loops lhs
469  CHKERR DMMoFEMTSSetIJacobian(dm, "THERMAL_FE", &thermal_elements.feLhs, NULL,
470  NULL);
471  CHKERR DMMoFEMTSSetIJacobian(dm, "THERMAL_CONVECTION_FE",
472  &thermal_elements.feConvectionLhs, NULL, NULL);
473  CHKERR DMMoFEMTSSetIJacobian(dm, "THERMAL_RADIATION_FE",
474  &thermal_elements.feRadiationLhs, NULL, NULL);
475 #ifdef __GROUNDSURFACETEMERATURE_HPP
476  if (ground_temperature_analysis) {
477  CHKERR DMMoFEMTSSetIJacobian(dm, "GROUND_SURFACE_FE",
478  &ground_surface.getFeGroundSurfaceLhs(), NULL,
479  NULL);
480  }
481 #endif //__GROUNDSURFACETEMERATURE_HPP
482 
483  //postprocess
484  CHKERR DMMoFEMTSSetIFunction(dm,DM_NO_ELEMENT,NULL,NULL,&dirichlet_bc);
485  CHKERR DMMoFEMTSSetIJacobian(dm,DM_NO_ELEMENT,NULL,NULL,&dirichlet_bc);
486 
487  TsCtx *ts_ctx;
488  DMMoFEMGetTsCtx(dm,&ts_ctx);
489  //add monitor operator
490  ts_ctx->get_postProcess_to_do_Monitor().push_back(&monitor);
491  ts_ctx->get_postProcess_to_do_Monitor().push_back(&post_proc);
492 
493  //create time solver
494  TS ts;
495  CHKERR TSCreate(PETSC_COMM_WORLD,&ts);
496  CHKERR TSSetType(ts,TSBEULER);
497 
498  CHKERR TSSetIFunction(ts,F,PETSC_NULL,PETSC_NULL);
499  CHKERR TSSetIJacobian(ts,A,A,PETSC_NULL,PETSC_NULL);
500  //add monitor to TS solver
501  CHKERR TSMonitorSet(ts,f_TSMonitorSet,ts_ctx,PETSC_NULL); // !!!
502 
503  CHKERR recorder_ptr->add_series_recorder("THEMP_SERIES");
504  //start to record
505  CHKERR recorder_ptr->initialize_series_recorder("THEMP_SERIES");
506 
507  double ftime = 1;
508  CHKERR TSSetDuration(ts,PETSC_DEFAULT,ftime);
509  CHKERR TSSetFromOptions(ts);
510  CHKERR TSSetDM(ts,dm);
511 
512  CHKERR TSSolve(ts,T);
513  CHKERR TSGetTime(ts,&ftime);
514 
515  //end recoder
516  CHKERR recorder_ptr->finalize_series_recorder("THEMP_SERIES");
517 
518  PetscInt steps,snesfails,rejects,nonlinits,linits;
519  CHKERR TSGetTimeStepNumber(ts,&steps);
520  CHKERR TSGetSNESFailures(ts,&snesfails);
521  CHKERR TSGetStepRejections(ts,&rejects);
522  CHKERR TSGetSNESIterations(ts,&nonlinits);
523  CHKERR TSGetKSPIterations(ts,&linits);
524 
525  PetscPrintf(PETSC_COMM_WORLD,
526  "steps %D (%D rejected, %D SNES fails), ftime %g, nonlinits %D, "
527  "linits %D\n",
528  steps, rejects, snesfails, ftime, nonlinits, linits);
529 
530  // save solution, if boundary conditions are defined you can use that file in
531  // mechanical problem to calculate thermal stresses
532  PetscBool is_partitioned = PETSC_FALSE;
533  CHKERR PetscOptionsGetBool(PETSC_NULL, PETSC_NULL, "-dm_is_partitioned",
534  &is_partitioned, PETSC_NULL);
535  if (is_partitioned) {
536  CHKERR moab.write_file("solution.h5m");
537  CHKERRG(rval);
538  } else {
539  if (pcomm->rank() == 0) {
540  CHKERR moab.write_file("solution.h5m");
541  CHKERRG(rval);
542  }
543  }
544 
545  CHKERR MatDestroy(&A);
546  CHKERR VecDestroy(&T);
547  CHKERR VecDestroy(&F);
548 
549  CHKERR TSDestroy(&ts);
550 
551  }
552  CATCH_ERRORS;
553 
554  return MoFEM::Core::Finalize();
555 }
virtual MoFEMErrorCode delete_recorder_series(const std::string &series_name)
#define DM_NO_ELEMENT
Definition: DMMoFEM.hpp:23
virtual MoFEMErrorCode finalize_series_recorder(const std::string &serie_name)
BasicMethodsSequence & get_postProcess_to_do_Monitor()
Definition: TsCtx.hpp:87
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:469
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:449
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:558
PetscErrorCode f_TSMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Definition: TsCtx.cpp:187
Core (interface) class.
Definition: Core.hpp:50
static char help[]
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:915
TS monitore it records temperature at time steps mofem_thermal_elem.
Projection of edge entities with one mid-node on hierarchical basis.
virtual MoFEMErrorCode build_finite_elements(int verb=-1)=0
Build finite elementsBuild finite element data structures. Have to be run before problem and adjacenc...
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM...
Definition: DMMMoFEM.cpp:895
PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const std::string &fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::FEMethod > pre_only, boost::shared_ptr< MoFEM::FEMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:787
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:51
virtual MoFEMErrorCode add_field(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_EXCL, int verb=-1)=0
Add field.
PetscErrorCode PetscOptionsGetString(PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Common.hpp:147
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
virtual MoFEMErrorCode initialize_series_recorder(const std::string &serie_name)
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
virtual bool check_series(const std::string &name) const
check if series is in database
Basic algebra on fields.
Definition: FieldBlas.hpp:33
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:481
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=-1)=0
Add entities to field meshsetThe lower dimension entities are added depending on the space type...
Managing BitRefLevels.
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=-1)=0
Set order approximation of the entities in the field.
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:91
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:141
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:28
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:492
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:150
#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...
#define CHKERR
Inline error check.
Definition: definitions.h:610
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::FEMethod *pre_only, MoFEM::FEMethod *post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:737
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:316
virtual MoFEMErrorCode build_fields(int verb=-1)=0
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, EntMethod &method, int lower_rank, int upper_rank, int verb=-1)=0
Make a loop over entities.
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition: DMMMoFEM.cpp:870
PetscErrorCode monitor(TS ts, PetscInt step, PetscReal crtime, Vec u, void *ctx)
virtual MoFEMErrorCode add_series_recorder(const std::string &series_name)
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:462
continuous field
Definition: definitions.h:176
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=-1)=0
build adjacencies
Implementation of ground surface temperature.
this calass is to control time stepping mofem_thermal_elem
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:61
structure grouping operators and data used for thermal problemsIn order to assemble matrices and righ...

Variable Documentation

◆ help

char help[]
static
Initial value:
=
"-my_file mesh file\n"
"-order set approx. order to all blocks\n"
"-my_block_config set block data\n"
"-my_ground_analysis_data data for crude climate model\n"
"\n"

Definition at line 44 of file thermal_unsteady.cpp.