v0.14.0
Loading...
Searching...
No Matches
Functions | Variables
magnetostatic.cpp File Reference
#include <BasicFiniteElements.hpp>
#include <MagneticElement.hpp>

Go to the source code of this file.

Functions

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

Variables

static char help []
 

Function Documentation

◆ main()

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

Definition at line 19 of file magnetostatic.cpp.

19 {
20
21 const string default_options = "-ksp_type fgmres \n"
22 "-pc_type lu \n"
23 "-pc_factor_mat_solver_type mumps \n"
24 "-mat_mumps_icntl_20 0 \n"
25 "-mat_mumps_icntl_13 1 \n"
26 "-ksp_monitor \n"
27 "-mat_mumps_icntl_24 1 \n"
28 "-mat_mumps_icntl_13 1 \n";
29
30 string param_file = "param_file.petsc";
31 if (!static_cast<bool>(ifstream(param_file))) {
32 std::ofstream file(param_file.c_str(), std::ios::ate);
33 if (file.is_open()) {
34 file << default_options;
35 file.close();
36 }
37 }
38
39 MoFEM::Core::Initialize(&argc, &argv, param_file.c_str(), help);
40
41 try {
42
43 moab::Core mb_instance;
44 moab::Interface &moab = mb_instance;
45
46 ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
47 auto moab_comm_wrap =
48 boost::make_shared<WrapMPIComm>(PETSC_COMM_WORLD, false);
49 if (pcomm == NULL)
50 pcomm =
51 new ParallelComm(&moab, moab_comm_wrap->get_comm());
52
53 // Read parameters from line command
54 PetscBool flg_file;
55 char mesh_file_name[255];
56 PetscInt order = 2;
57 PetscBool is_partitioned = PETSC_FALSE;
58 PetscBool regression_test = PETSC_FALSE;
59
60 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "Magnetostatics options",
61 "none");
62 CHKERR PetscOptionsString("-my_file", "mesh file name", "", "mesh.h5m",
63 mesh_file_name, 255, &flg_file);
64 CHKERR PetscOptionsInt("-my_order", "default approximation order", "",
65 order, &order, PETSC_NULL);
66 CHKERR PetscOptionsBool(
67 "-regression_test",
68 "if set norm of solution vector is check agains expected value ",
69 "",
70 PETSC_FALSE, &regression_test, PETSC_NULL);
71
72 ierr = PetscOptionsEnd();
74
75 // Read mesh to MOAB
76 const char *option;
77 option = "PARALLEL=READ_PART;"
78 "PARALLEL_RESOLVE_SHARED_ENTS;"
79 "PARTITION=PARALLEL_PARTITION;";
80 CHKERR moab.load_file(mesh_file_name, 0, option);
81
82 // Create mofem interface
83 MoFEM::Core core(moab);
84 MoFEM::Interface &m_field = core;
85
86 MagneticElement magnetic(m_field);
87 magnetic.blockData.oRder = order;
88 CHKERR magnetic.getNaturalBc();
89 CHKERR magnetic.getEssentialBc();
90 CHKERR magnetic.createFields();
91 CHKERR magnetic.createElements();
92 CHKERR magnetic.createProblem();
93 CHKERR magnetic.solveProblem(regression_test == PETSC_TRUE);
94 CHKERR magnetic.postProcessResults();
95 CHKERR magnetic.destroyProblem();
96
97 // write solution to file (can be used by lorentz_force example)
98 CHKERR moab.write_file("solution.h5m", "MOAB", "PARALLEL=WRITE_PART");
99 }
101
103 return 0;
104}
const std::string default_options
std::string param_file
#define CATCH_ERRORS
Catch errors.
#define MYPCOMM_INDEX
default communicator number PCOMM
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
#define CHKERR
Inline error check.
constexpr int order
static char help[]
char mesh_file_name[255]
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Implementation of magneto-static problem (basic Implementation)
Core (interface) class.
Definition Core.hpp:82
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition Core.cpp:72
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition Core.cpp:112
Deprecated interface functions.

Variable Documentation

◆ help

char help[]
static
Initial value:
= "-my_file mesh file name\n"
"-my_order default approximation order\n"
"-my_is_partitioned set if mesh is partitioned\n"
"-regression_test check solution vector against expected value\n"
"\n\n"

Definition at line 13 of file magnetostatic.cpp.