v0.8.14
Classes | Typedefs | Functions | Variables
h_adaptive_transport.cpp File Reference

Example implementation of transport problem using ultra-week formulation. More...

#include <BasicFiniteElements.hpp>
#include <MixTransportElement.hpp>

Go to the source code of this file.

Classes

struct  BcFluxData
 
struct  MyTransport
 Application of mix transport data structure. More...
 

Typedefs

typedef map< int, BcFluxDataBcFluxMap
 

Functions

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

Variables

static char help [] = "...\n\n"
 

Detailed Description

Example implementation of transport problem using ultra-week formulation.

Todo:
Should be implemented and tested problem from this article Demkowicz, Leszek, and Jayadeep Gopalakrishnan. "Analysis of the DPG method for the Poisson equation." SIAM Journal on Numerical Analysis 49.5 (2011): 1788-1809.

Definition in file h_adaptive_transport.cpp.

Typedef Documentation

◆ BcFluxMap

typedef map<int, BcFluxData> BcFluxMap

Definition at line 43 of file h_adaptive_transport.cpp.

Function Documentation

◆ main()

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

Definition at line 374 of file h_adaptive_transport.cpp.

374  {
375 
376  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
377 
378  try {
379 
380  moab::Core mb_instance;
381  moab::Interface &moab = mb_instance;
382  int rank;
383  MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
384 
385  // get file name form command line
386  PetscBool flg = PETSC_TRUE;
387  char mesh_file_name[255];
388  CHKERR PetscOptionsGetString(PETSC_NULL, PETSC_NULL, "-my_file",
389  mesh_file_name, 255, &flg);
390  if (flg != PETSC_TRUE) {
391  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
392  "*** ERROR -my_file (MESH FILE NEEDED)");
393  }
394 
395  // create MOAB communicator
396  ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
397  if (pcomm == NULL)
398  pcomm = new ParallelComm(&moab, PETSC_COMM_WORLD);
399 
400  const char *option;
401  option = "";
402  CHKERR moab.load_file(mesh_file_name, 0, option);
403 
404  // Create mofem interface
405  MoFEM::Core core(moab);
406  MoFEM::Interface &m_field = core;
407 
408  // Add meshsets with material and boundary conditions
409  MeshsetsManager *meshsets_manager_ptr;
410  CHKERR m_field.getInterface(meshsets_manager_ptr);
411  CHKERR meshsets_manager_ptr->setMeshsetFromFile();
412 
413  PetscPrintf(PETSC_COMM_WORLD,
414  "Read meshsets add added meshsets for bc.cfg\n");
415  for (_IT_CUBITMESHSETS_FOR_LOOP_(m_field, it)) {
416  PetscPrintf(PETSC_COMM_WORLD, "%s",
417  static_cast<std::ostringstream &>(
418  std::ostringstream().seekp(0) << *it << endl)
419  .str()
420  .c_str());
421  cerr << *it << endl;
422  }
423 
424  // set entities bit level
425  BitRefLevel ref_level;
426  ref_level.set(0);
427  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
428  0, 3, ref_level);
429 
430  // set app. order
431  // see Hierarchic Finite Element Bases on Unstructured Tetrahedral Meshes
432  // (Mark Ainsworth & Joe Coyle)
433  PetscInt order;
434  CHKERR PetscOptionsGetInt(PETSC_NULL, PETSC_NULL, "-my_order", &order,
435  &flg);
436  if (flg != PETSC_TRUE) {
437  order = 0;
438  }
439 
440  // finite elements
441 
442  BcFluxMap bc_flux_map;
443  MyTransport ufe(m_field, bc_flux_map);
444 
445  // Initially calculate problem on coarse mesh
446 
447  CHKERR ufe.addFields("VALUES", "FLUXES", order);
448  CHKERR ufe.addFiniteElements("FLUXES", "VALUES");
449  // Set boundary conditions
450  CHKERR ufe.addBoundaryElements(ref_level);
451  CHKERR ufe.buildProblem(ref_level);
452  CHKERR ufe.createMatrices();
453  CHKERR ufe.solveLinearProblem();
454  CHKERR ufe.calculateResidual();
455  CHKERR ufe.evaluateError();
456  CHKERR ufe.destroyMatrices();
457  CHKERR ufe.postProc("out_0.h5m");
458 
459  int nb_levels = 5; // default number of refinement levels
460  // get number of refinement levels form command line
461  CHKERR PetscOptionsGetInt(PETSC_NULL, PETSC_NULL, "-nb_levels", &nb_levels,
462  PETSC_NULL);
463 
464  // refine mesh, solve problem and do it again until number of refinement
465  // levels are exceeded.
466  for (int ll = 1; ll != nb_levels; ll++) {
467  const int nb_levels = ll;
468  CHKERR ufe.squashBits();
469  CHKERR ufe.refineMesh(ufe, nb_levels, order);
470  ref_level = BitRefLevel().set(nb_levels);
471  bc_flux_map.clear();
472  CHKERR ufe.addBoundaryElements(ref_level);
473  CHKERR ufe.buildProblem(ref_level);
474  CHKERR ufe.createMatrices();
475  CHKERR ufe.solveLinearProblem();
476  CHKERR ufe.calculateResidual();
477  CHKERR ufe.evaluateError();
478  CHKERR ufe.destroyMatrices();
479  CHKERR ufe.postProc(
480  static_cast<std::ostringstream &>(std::ostringstream().seekp(0)
481  << "out_" << nb_levels << ".h5m")
482  .str());
483  }
484  }
485  CATCH_ERRORS;
486 
488 
489  return 0;
490 }
Interface for managing meshsets containing materials and boundary conditions.
Core (interface) class.
Definition: Core.hpp:50
char mesh_file_name[255]
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
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.
MoFEMErrorCode setMeshsetFromFile(const string file_name, const bool clean_file_options=true)
add blocksets reading config file
Managing BitRefLevels.
#define CHKERR
Inline error check.
Definition: definitions.h:578
map< int, BcFluxData > BcFluxMap
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:282
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:429
Application of mix transport data structure.
static char help[]
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:61

Variable Documentation

◆ help

char help[] = "...\n\n"
static

Definition at line 32 of file h_adaptive_transport.cpp.