v0.6.9
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 364 of file h_adaptive_transport.cpp.

364  {
365 
366  PetscInitialize(&argc,&argv,(char *)0,help);
367 
368  try {
369 
370  moab::Core mb_instance;
371  moab::Interface& moab = mb_instance;
372  int rank;
373  MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
374 
375  // get file name form command line
376  PetscBool flg = PETSC_TRUE;
377  char mesh_file_name[255];
378  ierr = PetscOptionsGetString(PETSC_NULL,PETSC_NULL,"-my_file",mesh_file_name,255,&flg); CHKERRG(ierr);
379  if(flg != PETSC_TRUE) {
380  SETERRQ(PETSC_COMM_SELF,MOFEM_INVALID_DATA,"*** ERROR -my_file (MESH FILE NEEDED)");
381  }
382 
383  // create MOAB communicator
384  ParallelComm* pcomm = ParallelComm::get_pcomm(&moab,MYPCOMM_INDEX);
385  if(pcomm == NULL) pcomm = new ParallelComm(&moab,PETSC_COMM_WORLD);
386 
387  const char *option;
388  option = "";
389  rval = moab.load_file(mesh_file_name, 0, option); CHKERRG(rval);
390 
391  //Create mofem interface
392  MoFEM::Core core(moab);
393  MoFEM::Interface& m_field = core;
394 
395  // Add meshsets with material and boundary conditions
396  MeshsetsManager *meshsets_manager_ptr;
397  ierr = m_field.getInterface(meshsets_manager_ptr); CHKERRG(ierr);
398  ierr = meshsets_manager_ptr->setMeshsetFromFile(); CHKERRG(ierr);
399 
400  PetscPrintf(PETSC_COMM_WORLD,"Read meshsets add added meshsets for bc.cfg\n");
401  for(_IT_CUBITMESHSETS_FOR_LOOP_(m_field,it)) {
402  PetscPrintf(
403  PETSC_COMM_WORLD,
404  "%s",static_cast<std::ostringstream&>(std::ostringstream().seekp(0) << *it << endl).str().c_str()
405  );
406  cerr << *it << endl;
407  }
408 
409  //set entities bit level
410  BitRefLevel ref_level;
411  ref_level.set(0);
412  ierr = m_field.seed_ref_level_3D(0,ref_level); CHKERRG(ierr);
413 
414  //set app. order
415  //see Hierarchic Finite Element Bases on Unstructured Tetrahedral Meshes (Mark Ainsworth & Joe Coyle)
416  PetscInt order;
417  ierr = PetscOptionsGetInt(PETSC_NULL,PETSC_NULL,"-my_order",&order,&flg); CHKERRG(ierr);
418  if(flg != PETSC_TRUE) {
419  order = 0;
420  }
421 
422  //finite elements
423 
424  BcFluxMap bc_flux_map;
425  MyTransport ufe(m_field,bc_flux_map);
426 
427  // Initially calculate problem on coarse mesh
428 
429  ierr = ufe.addFields("VALUES","FLUXES",order); CHKERRG(ierr);
430  ierr = ufe.addFiniteElements("FLUXES","VALUES"); CHKERRG(ierr);
431  // Set boundary conditions
432  ierr = ufe.addBoundaryElements(ref_level);
433  ierr = ufe.buildProblem(ref_level); CHKERRG(ierr);
434  ierr = ufe.createMatrices(); CHKERRG(ierr);
435  ierr = ufe.solveLinearProblem(); CHKERRG(ierr);
436  ierr = ufe.calculateResidual(); CHKERRG(ierr);
437  ierr = ufe.evaluateError(); CHKERRG(ierr);
438  ierr = ufe.destroyMatrices(); CHKERRG(ierr);
439  ierr = ufe.postProc("out_0.h5m"); CHKERRG(ierr);
440 
441  int nb_levels = 5; // default number of refinement levels
442  // get number of refinement levels form command line
443  ierr = PetscOptionsGetInt(PETSC_NULL,PETSC_NULL,"-nb_levels",&nb_levels,PETSC_NULL); CHKERRG(ierr);
444 
445  // refine mesh, solve problem and do it again until number of refinement levels are exceeded.
446  for(int ll = 1;ll!=nb_levels;ll++) {
447  const int nb_levels = ll;
448  ierr = ufe.squashBits(); CHKERRG(ierr);
449  ierr = ufe.refineMesh(ufe,nb_levels,order); CHKERRG(ierr);
450  ref_level = BitRefLevel().set(nb_levels);
451  bc_flux_map.clear();
452  ierr = ufe.addBoundaryElements(ref_level);
453  ierr = ufe.buildProblem(ref_level); CHKERRG(ierr);
454  ierr = ufe.createMatrices(); CHKERRG(ierr);
455  ierr = ufe.solveLinearProblem(); CHKERRG(ierr);
456  ierr = ufe.calculateResidual(); CHKERRG(ierr);
457  ierr = ufe.evaluateError(); CHKERRG(ierr);
458  ierr = ufe.destroyMatrices(); CHKERRG(ierr);
459  ierr = ufe.postProc(
460  static_cast<std::ostringstream&>
461  (std::ostringstream().seekp(0) << "out_" << nb_levels << ".h5m").str()
462  ); CHKERRG(ierr);
463  }
464 
465  } catch (MoFEMException const &e) {
466  SETERRQ(PETSC_COMM_SELF,e.errorCode,e.errorMessage);
467  }
468 
469  ierr = PetscFinalize(); CHKERRG(ierr);
470 
471  return 0;
472 
473 }
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:511
Interface for managing meshsets containing materials and boundary conditions.
Exception to catch.
Definition: Common.hpp:26
Core (interface) class.
Definition: Core.hpp:50
char errorMessage[255]
Definition: Common.hpp:28
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
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
const int errorCode
Definition: Common.hpp:27
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:315
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
Application of mix transport data structure.
DEPRECATED MoFEMErrorCode seed_ref_level_3D(const EntityHandle meshset, const BitRefLevel &bit, int verb=-1)
seed 2D entities in the meshset and their adjacencies (only TETs adjacencies) in a particular BitRefL...
static char help[]
map< int, BcFluxData > BcFluxMap

Variable Documentation

◆ help

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

Definition at line 32 of file h_adaptive_transport.cpp.