h_adaptive_transport.cpp File Reference

Example implementation of transport problem using mixed formulation. More...

#include <MoFEM.hpp>
#include <MethodForForceScaling.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, BcFluxData >	BcFluxMap

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

Variables
static char	help []

Detailed Description

Example implementation of transport problem using mixed 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 37 of file h_adaptive_transport.cpp.

Function Documentation

main()

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

Definition at line 357 of file h_adaptive_transport.cpp.

                                 {
 
  const string default_options = "-ksp_type fgmres \n"
                                 "-pc_type lu \n"
                                 "-pc_factor_mat_solver_type mumps \n"
                                 "-mat_mumps_icntl_20 0 \n"
                                 "-ksp_monitor\n";
 
  string param_file = "param_file.petsc";
  if (!static_cast<bool>(ifstream(param_file))) {
    std::ofstream file(param_file.c_str(), std::ios::ate);
    if (file.is_open()) {
      file << default_options;
      file.close();
    }
  }
 
  MoFEM::Core::Initialize(&argc, &argv, param_file.c_str(), help);
 
  try {
 
    moab::Core mb_instance;
    moab::Interface &moab = mb_instance;
    int rank;
    MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
 
    // get file name form command line
    PetscBool flg = PETSC_TRUE;
    char mesh_file_name[255];
    CHKERR PetscOptionsGetString(PETSC_NULL, PETSC_NULL, "-my_file",
                                 mesh_file_name, 255, &flg);
    if (flg != PETSC_TRUE) {
      SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
              "*** ERROR -my_file (MESH FILE NEEDED)");
    }
 
    const char *option;
    option = "";
    CHKERR moab.load_file(mesh_file_name, 0, option);
 
    // Create mofem interface
    MoFEM::Core core(moab);
    MoFEM::Interface &m_field = core;
 
    // Add meshsets with material and boundary conditions
    MeshsetsManager *meshsets_manager_ptr;
    CHKERR m_field.getInterface(meshsets_manager_ptr);
    CHKERR meshsets_manager_ptr->setMeshsetFromFile();
 
    PetscPrintf(PETSC_COMM_WORLD,
                "Read meshsets add added meshsets for bc.cfg\n");
    for (_IT_CUBITMESHSETS_FOR_LOOP_(m_field, it)) {
      PetscPrintf(PETSC_COMM_WORLD, "%s",
                  static_cast<std::ostringstream &>(
                      std::ostringstream().seekp(0) << *it << endl)
                      .str()
                      .c_str());
      cerr << *it << endl;
    }
 
    // set entities bit level
    BitRefLevel ref_level;
    ref_level.set(0);
    CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
        0, 3, ref_level);
 
    // set app. order
    // see Hierarchic Finite Element Bases on Unstructured Tetrahedral Meshes
    // (Mark Ainsworth & Joe Coyle)
    PetscInt order;
    CHKERR PetscOptionsGetInt(PETSC_NULL, PETSC_NULL, "-my_order", &order,
                              &flg);
    if (flg != PETSC_TRUE) {
      order = 0;
    }
 
    // finite elements
 
    BcFluxMap bc_flux_map;
    MyTransport ufe(m_field, bc_flux_map);
 
    // Initially calculate problem on coarse mesh
 
    CHKERR ufe.addFields("VALUES", "FLUXES", order);
    CHKERR ufe.addFiniteElements("FLUXES", "VALUES");
    // Set boundary conditions
    CHKERR ufe.addBoundaryElements(ref_level);
    CHKERR ufe.buildProblem(ref_level);
    CHKERR ufe.createMatrices();
    CHKERR ufe.solveLinearProblem();
    CHKERR ufe.calculateResidual();
    CHKERR ufe.evaluateError();
    CHKERR ufe.destroyMatrices();
    CHKERR ufe.postProc("out_0.h5m");
 
    int nb_levels = 5; // default number of refinement levels
    // get number of refinement levels form command line
    CHKERR PetscOptionsGetInt(PETSC_NULL, PETSC_NULL, "-nb_levels", &nb_levels,
                              PETSC_NULL);
 
    // refine mesh, solve problem and do it again until number of refinement
    // levels are exceeded.
    for (int ll = 1; ll != nb_levels; ll++) {
      const int nb_levels = ll;
      CHKERR ufe.squashBits();
      CHKERR ufe.refineMesh(ufe, nb_levels, order);
      ref_level = BitRefLevel().set(nb_levels);
      bc_flux_map.clear();
      CHKERR ufe.addBoundaryElements(ref_level);
      CHKERR ufe.buildProblem(ref_level);
      CHKERR ufe.createMatrices();
      CHKERR ufe.solveLinearProblem();
      CHKERR ufe.calculateResidual();
      CHKERR ufe.evaluateError();
      CHKERR ufe.destroyMatrices();
      CHKERR ufe.postProc(
          static_cast<std::ostringstream &>(std::ostringstream().seekp(0)
                                            << "out_" << nb_levels << ".h5m")
              .str());
    }
  }
  CATCH_ERRORS;
 
  MoFEM::Core::Finalize();
 
  return 0;
}

Variable Documentation

help

char help[]

static

Initial value:

= "-my_file input file"
                     "-my_order order of approximation"
                     "-nb_levels number of refinements levels"
                     "\n\n"

Definition at line 23 of file h_adaptive_transport.cpp.