v0.5.86
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); CHKERRQ(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); CHKERRQ_MOAB(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.query_interface(meshsets_manager_ptr); CHKERRQ(ierr);
398  ierr = meshsets_manager_ptr->setMeshsetFromFile(); CHKERRQ(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); CHKERRQ(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); CHKERRQ(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); CHKERRQ(ierr);
430  ierr = ufe.addFiniteElements("FLUXES","VALUES"); CHKERRQ(ierr);
431  // Set boundary conditions
432  ierr = ufe.addBoundaryElements(ref_level);
433  ierr = ufe.buildProblem(ref_level); CHKERRQ(ierr);
434  ierr = ufe.createMatrices(); CHKERRQ(ierr);
435  ierr = ufe.solveLinearProblem(); CHKERRQ(ierr);
436  ierr = ufe.calculateResidual(); CHKERRQ(ierr);
437  ierr = ufe.evaluateError(); CHKERRQ(ierr);
438  ierr = ufe.destroyMatrices(); CHKERRQ(ierr);
439  ierr = ufe.postProc("out_0.h5m"); CHKERRQ(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); CHKERRQ(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(); CHKERRQ(ierr);
449  ierr = ufe.refineMesh(ufe,nb_levels,order); CHKERRQ(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); CHKERRQ(ierr);
454  ierr = ufe.createMatrices(); CHKERRQ(ierr);
455  ierr = ufe.solveLinearProblem(); CHKERRQ(ierr);
456  ierr = ufe.calculateResidual(); CHKERRQ(ierr);
457  ierr = ufe.evaluateError(); CHKERRQ(ierr);
458  ierr = ufe.destroyMatrices(); CHKERRQ(ierr);
459  ierr = ufe.postProc(
460  static_cast<std::ostringstream&>
461  (std::ostringstream().seekp(0) << "out_" << nb_levels << ".h5m").str()
462  ); CHKERRQ(ierr);
463  }
464 
465  } catch (MoFEMException const &e) {
466  SETERRQ(PETSC_COMM_SELF,e.errorCode,e.errorMessage);
467  }
468 
469  ierr = PetscFinalize(); CHKERRQ(ierr);
470 
471  return 0;
472 
473 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:362
MoFEMErrorCodes errorCode
Definition: Common.hpp:105
static MoABErrorCode rval
Definition: Common.hpp:25
PetscErrorCode setMeshsetFromFile(const string file_name, const bool clean_file_options=true)
add blocksets reading config file
Interface for managing meshsets containing materials and boundary conditions.
Core (interface) classThis is the implementation of abstract MoFEM::Interface class. Similarly to the convention used in MoAB, we use small letters to name function of purely abstract classes. This is an exception used only here. For more details about naming functions see Coding practice.
Definition: Core.hpp:44
char errorMessage[255]
Definition: Common.hpp:106
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:48
PetscErrorCode query_interface(IFace *&ptr) const
Definition: Interface.hpp:47
CHKERRQ(ierr)
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
virtual DEPRECATED PetscErrorCode seed_ref_level_3D(const EntityHandle meshset, const BitRefLevel &bit, int verb=-1)=0
seed 2D entities in the meshset and their adjacencies (only TETs adjacencies) in a particular BitRefL...
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:305
InterfaceThis interface is used by user to:
Definition: Interface.hpp:38
Application of mix transport data structure.
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.