v0.14.0
Loading...
Searching...
No Matches
poisson_2d_dis_galerkin.cpp
Go to the documentation of this file.
1/**
2 * \file poisson_2d_dis_galerkin.cpp
3 * \example poisson_2d_dis_galerkin.cpp
4 *
5 * Example of implementation for discontinuous Galerkin.
6 */
7
9
10constexpr int BASE_DIM = 1;
11constexpr int FIELD_DIM = 1;
12constexpr int SPACE_DIM = 2;
13
15
18
24using FaceSideOp = FaceSideEle::UserDataOperator;
25
27
28static double penalty = 1e6;
29static double phi =
30 -1; // 1 - symmetric Nitsche, 0 - nonsymmetric, -1 antisymmetrica
31static double nitsche = 1;
32
34
39
40PetscBool is_test = PETSC_FALSE;
41
42auto u_exact = [](const double x, const double y, const double) {
43 if (is_test)
44 return x * x * y * y;
45 else
46 return cos(2 * x * M_PI) * cos(2 * y * M_PI);
47};
48
49auto u_grad_exact = [](const double x, const double y, const double) {
50 if (is_test)
51 return FTensor::Tensor1<double, 2>{2 * x * y * y, 2 * x * x * y};
52 else
54
55 -2 * M_PI * cos(2 * M_PI * y) * sin(2 * M_PI * x),
56 -2 * M_PI * cos(2 * M_PI * x) * sin(2 * M_PI * y)
57
58 };
59};
60
61auto source = [](const double x, const double y, const double) {
62 if (is_test)
63 return -(2 * x * x + 2 * y * y);
64 else
65 return 8 * M_PI * M_PI * cos(2 * x * M_PI) * cos(2 * y * M_PI);
66};
67
68using namespace MoFEM;
70
71static char help[] = "...\n\n";
72
74public:
76
77 // Declaration of the main function to run analysis
79
80private:
81 // Declaration of other main functions called in runProgram()
90
91 // MoFEM interfaces
94
95 // Field name and approximation order
96 std::string domainField;
97 int oRder;
98};
99
101 : domainField("U"), mField(m_field), oRder(4) {}
102
103//! [Read mesh]
106
109
110 // Only L2 field is set in this example. Two lines bellow forces simple
111 // interface to creat lower dimension (edge) elements, despite that fact that
112 // there is no field spanning on such elements. We need them for DG method.
115
117
119}
120//! [Read mesh]
121
122//! [Setup problem]
125
126 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &oRder, PETSC_NULL);
127 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-penalty", &penalty,
128 PETSC_NULL);
129 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-phi", &phi, PETSC_NULL);
130 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-nitsche", &nitsche,
131 PETSC_NULL);
132 PetscOptionsGetBool(PETSC_NULL, "", "-is_test", &is_test, PETSC_NULL);
133
134 MOFEM_LOG("WORLD", Sev::inform) << "Set order: " << oRder;
135 MOFEM_LOG("WORLD", Sev::inform) << "Set penalty: " << penalty;
136 MOFEM_LOG("WORLD", Sev::inform) << "Set phi: " << phi;
137 MOFEM_LOG("WORLD", Sev::inform) << "Set nitche: " << nitsche;
138 MOFEM_LOG("WORLD", Sev::inform) << "Set test: " << (is_test == PETSC_TRUE);
139
144
145 // This is only for debigging and experimentation, to see boundary and edge
146 // elements.
147 auto save_shared = [&](auto meshset, std::string prefix) {
149 auto file_name =
150 prefix + "_" +
151 boost::lexical_cast<std::string>(mField.get_comm_rank()) + ".vtk";
152 CHKERR mField.get_moab().write_file(file_name.c_str(), "VTK", "", &meshset,
153 1);
155 };
156
157 // CHKERR save_shared(simpleInterface->getBoundaryMeshSet(), "bdy");
158 // CHKERR save_shared(simpleInterface->getSkeletonMeshSet(), "skel");
159
161}
162//! [Setup problem]
163
164//! [Boundary condition]
168}
169
170//! [Boundary condition]
171
172//! [Assemble system]
175
176 auto pipeline_mng = mField.getInterface<PipelineManager>();
177
178 auto add_base_ops = [&](auto &pipeline) {
179 auto det_ptr = boost::make_shared<VectorDouble>();
180 auto jac_ptr = boost::make_shared<MatrixDouble>();
181 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
182 pipeline.push_back(new OpCalculateHOJacForFace(jac_ptr));
183 pipeline.push_back(new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
184 pipeline.push_back(new OpSetInvJacL2ForFace(inv_jac_ptr));
185 };
186
187 add_base_ops(pipeline_mng->getOpDomainLhsPipeline());
188 pipeline_mng->getOpDomainLhsPipeline().push_back(new OpDomainGradGrad(
190 [](const double, const double, const double) { return 1; }));
191 pipeline_mng->getOpDomainRhsPipeline().push_back(
193
194 // Push operators to the Pipeline for Skeleton
195 auto side_fe_ptr = boost::make_shared<FaceSideEle>(mField);
196 add_base_ops(side_fe_ptr->getOpPtrVector());
197 side_fe_ptr->getOpPtrVector().push_back(
199
200 // Push operators to the Pipeline for Skeleton
201 pipeline_mng->getOpSkeletonLhsPipeline().push_back(
202 new OpL2LhsPenalty(side_fe_ptr));
203
204 // Push operators to the Pipeline for Boundary
205 pipeline_mng->getOpBoundaryLhsPipeline().push_back(
206 new OpL2LhsPenalty(side_fe_ptr));
207 pipeline_mng->getOpBoundaryRhsPipeline().push_back(
208 new OpL2BoundaryRhs(side_fe_ptr, u_exact));
209
211}
212//! [Assemble system]
213
214//! [Set integration rules]
217
218 auto rule_lhs = [](int, int, int p) -> int { return 2 * p; };
219 auto rule_rhs = [](int, int, int p) -> int { return 2 * p; };
220 auto rule_2 = [this](int, int, int) { return 2 * oRder; };
221
222 auto pipeline_mng = mField.getInterface<PipelineManager>();
223 CHKERR pipeline_mng->setDomainLhsIntegrationRule(rule_lhs);
224 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule_rhs);
225
226 CHKERR pipeline_mng->setSkeletonLhsIntegrationRule(rule_2);
227 CHKERR pipeline_mng->setSkeletonRhsIntegrationRule(rule_2);
228 CHKERR pipeline_mng->setBoundaryLhsIntegrationRule(rule_2);
229 CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(rule_2);
230
232}
233//! [Set integration rules]
234
235//! [Solve system]
238
239 auto pipeline_mng = mField.getInterface<PipelineManager>();
240
241 auto ksp_solver = pipeline_mng->createKSP();
242 CHKERR KSPSetFromOptions(ksp_solver);
243 CHKERR KSPSetUp(ksp_solver);
244
245 // Create RHS and solution vectors
246 auto dm = simpleInterface->getDM();
247 auto F = createDMVector(dm);
248 auto D = vectorDuplicate(F);
249
250 CHKERR KSPSolve(ksp_solver, F, D);
251
252 // Scatter result data on the mesh
253 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
254 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
255 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
256
258}
259//! [Solve system]
260
263
264 auto pipeline_mng = mField.getInterface<PipelineManager>();
265 pipeline_mng->getDomainRhsFE().reset();
266 pipeline_mng->getDomainLhsFE().reset();
267 pipeline_mng->getSkeletonRhsFE().reset();
268 pipeline_mng->getSkeletonLhsFE().reset();
269 pipeline_mng->getBoundaryRhsFE().reset();
270 pipeline_mng->getBoundaryLhsFE().reset();
271
272 auto rule = [](int, int, int p) -> int { return 2 * p; };
273 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule);
274 auto rule_2 = [this](int, int, int) { return 2 * oRder; };
275 CHKERR pipeline_mng->setSkeletonRhsIntegrationRule(rule_2);
276
277 auto add_base_ops = [&](auto &pipeline) {
278 auto det_ptr = boost::make_shared<VectorDouble>();
279 auto jac_ptr = boost::make_shared<MatrixDouble>();
280 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
281 pipeline.push_back(new OpCalculateHOJacForFace(jac_ptr));
282 pipeline.push_back(new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
283 pipeline.push_back(new OpSetInvJacL2ForFace(inv_jac_ptr));
284 };
285
286 auto u_vals_ptr = boost::make_shared<VectorDouble>();
287 auto u_grad_ptr = boost::make_shared<MatrixDouble>();
288
289 add_base_ops(pipeline_mng->getOpDomainRhsPipeline());
290 pipeline_mng->getOpDomainRhsPipeline().push_back(
292 pipeline_mng->getOpDomainRhsPipeline().push_back(
294
295 enum NORMS { L2 = 0, SEMI_NORM, H1, LAST_NORM };
296 std::array<int, LAST_NORM> error_indices;
297 for (int i = 0; i != LAST_NORM; ++i)
298 error_indices[i] = i;
299 auto l2_vec = createVectorMPI(
300 mField.get_comm(), (!mField.get_comm_rank()) ? LAST_NORM : 0, LAST_NORM);
301
302 auto error_op = new DomainEleOp(domainField, DomainEleOp::OPROW);
303 error_op->doWorkRhsHook = [&](DataOperator *op_ptr, int side, EntityType type,
306 auto o = static_cast<DomainEleOp *>(op_ptr);
307
309
310 if (const size_t nb_dofs = data.getIndices().size()) {
311
312 const int nb_integration_pts = o->getGaussPts().size2();
313 auto t_w = o->getFTensor0IntegrationWeight();
314 auto t_val = getFTensor0FromVec(*u_vals_ptr);
315 auto t_grad = getFTensor1FromMat<2>(*u_grad_ptr);
316 auto t_coords = o->getFTensor1CoordsAtGaussPts();
317
318 std::array<double, LAST_NORM> error;
319 std::fill(error.begin(), error.end(), 0);
320
321 for (int gg = 0; gg != nb_integration_pts; ++gg) {
322
323 const double alpha = t_w * o->getMeasure();
324 const double diff =
325 t_val - u_exact(t_coords(0), t_coords(1), t_coords(2));
326
327 auto t_exact_grad = u_grad_exact(t_coords(0), t_coords(1), t_coords(2));
328
329 const double diff_grad =
330 (t_grad(i) - t_exact_grad(i)) * (t_grad(i) - t_exact_grad(i));
331
332 error[L2] += alpha * pow(diff, 2);
333 error[SEMI_NORM] += alpha * diff_grad;
334
335 ++t_w;
336 ++t_val;
337 ++t_grad;
338 ++t_coords;
339 }
340
341 error[H1] = error[L2] + error[SEMI_NORM];
342
343 CHKERR VecSetValues(l2_vec, LAST_NORM, error_indices.data(), error.data(),
344 ADD_VALUES);
345 }
346
348 };
349
350 auto side_fe_ptr = boost::make_shared<FaceSideEle>(mField);
351 add_base_ops(side_fe_ptr->getOpPtrVector());
352 side_fe_ptr->getOpPtrVector().push_back(
354 std::array<VectorDouble, 2> side_vals;
355 std::array<double, 2> area_map;
356
357 auto side_vals_op = new DomainEleOp(domainField, DomainEleOp::OPROW);
358 side_vals_op->doWorkRhsHook = [&](DataOperator *op_ptr, int side,
359 EntityType type,
362 auto o = static_cast<FaceSideOp *>(op_ptr);
363
364 const auto nb_in_loop = o->getFEMethod()->nInTheLoop;
365 area_map[nb_in_loop] = o->getMeasure();
366 side_vals[nb_in_loop] = *u_vals_ptr;
367 if (!nb_in_loop) {
368 area_map[1] = 0;
369 side_vals[1].clear();
370 }
371
373 };
374 side_fe_ptr->getOpPtrVector().push_back(side_vals_op);
375
376 auto do_work_rhs_error = [&](DataOperator *op_ptr, int side, EntityType type,
379 auto o = static_cast<BoundaryEleOp *>(op_ptr);
380
381 CHKERR o->loopSideFaces("dFE", *side_fe_ptr);
382 const auto in_the_loop = side_fe_ptr->nInTheLoop;
383
384#ifndef NDEBUG
385 const std::array<std::string, 2> ele_type_name = {"BOUNDARY", "SKELETON"};
386 MOFEM_LOG("SELF", Sev::noisy)
387 << "do_work_rhs_error in_the_loop " << ele_type_name[in_the_loop];
388#endif
389
390 const double s = o->getMeasure() / (area_map[0] + area_map[1]);
391 const double p = penalty * s;
392
393 constexpr std::array<int, 2> sign_array{1, -1};
394
395 std::array<double, LAST_NORM> error;
396 std::fill(error.begin(), error.end(), 0);
397
398 const int nb_integration_pts = o->getGaussPts().size2();
399
400 if (!in_the_loop) {
401 side_vals[1].resize(nb_integration_pts, false);
402 auto t_coords = o->getFTensor1CoordsAtGaussPts();
403 auto t_val_m = getFTensor0FromVec(side_vals[1]);
404 for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
405 t_val_m = u_exact(t_coords(0), t_coords(1), t_coords(2));
406 ++t_coords;
407 ++t_val_m;
408 }
409 }
410
411 auto t_val_p = getFTensor0FromVec(side_vals[0]);
412 auto t_val_m = getFTensor0FromVec(side_vals[1]);
413 auto t_w = o->getFTensor0IntegrationWeight();
414
415 for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
416
417 const double alpha = o->getMeasure() * t_w;
418 const auto diff = t_val_p - t_val_m;
419 error[SEMI_NORM] += alpha * p * diff * diff;
420
421 ++t_w;
422 ++t_val_p;
423 ++t_val_m;
424 }
425
426
427 error[H1] = error[SEMI_NORM];
428 CHKERR VecSetValues(l2_vec, LAST_NORM, error_indices.data(), error.data(),
429 ADD_VALUES);
430
432 };
433
434 auto skeleton_error_op = new BoundaryEleOp(NOSPACE, BoundaryEleOp::OPSPACE);
435 skeleton_error_op->doWorkRhsHook = do_work_rhs_error;
436 auto boundary_error_op = new BoundaryEleOp(NOSPACE, BoundaryEleOp::OPSPACE);
437 boundary_error_op->doWorkRhsHook = do_work_rhs_error;
438
439 pipeline_mng->getOpDomainRhsPipeline().push_back(error_op);
440 pipeline_mng->getOpSkeletonLhsPipeline().push_back(skeleton_error_op);
441 pipeline_mng->getOpBoundaryRhsPipeline().push_back(boundary_error_op);
442
443 CHKERR pipeline_mng->loopFiniteElements();
444
445 CHKERR VecAssemblyBegin(l2_vec);
446 CHKERR VecAssemblyEnd(l2_vec);
447
448 if (mField.get_comm_rank() == 0) {
449 const double *array;
450 CHKERR VecGetArrayRead(l2_vec, &array);
451 MOFEM_LOG_C("SELF", Sev::inform, "Error Norm L2 %6.4e",
452 std::sqrt(array[L2]));
453 MOFEM_LOG_C("SELF", Sev::inform, "Error Norm Energetic %6.4e",
454 std::sqrt(array[SEMI_NORM]));
455 MOFEM_LOG_C("SELF", Sev::inform, "Error Norm H1 %6.4e",
456 std::sqrt(array[H1]));
457
458 if(is_test) {
459 constexpr double eps = 1e-12;
460 if (std::sqrt(array[H1]) > eps)
461 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID, "Error is too big");
462 }
463
464 CHKERR VecRestoreArrayRead(l2_vec, &array);
465 const MoFEM::Problem *problem_ptr;
467 MOFEM_LOG_C("SELF", Sev::inform, "Nb. DOFs %d",
468 problem_ptr->getNbDofsRow());
469 }
470
471
472
474}
475
476//! [Output results]
479
480 auto pipeline_mng = mField.getInterface<PipelineManager>();
481 pipeline_mng->getDomainLhsFE().reset();
482 pipeline_mng->getSkeletonRhsFE().reset();
483 pipeline_mng->getSkeletonLhsFE().reset();
484 pipeline_mng->getBoundaryRhsFE().reset();
485 pipeline_mng->getBoundaryLhsFE().reset();
486
487 auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
488
489 auto u_ptr = boost::make_shared<VectorDouble>();
490 post_proc_fe->getOpPtrVector().push_back(
492
494
495 post_proc_fe->getOpPtrVector().push_back(
496
497 new OpPPMap(
498
499 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
500
501 {{"U", u_ptr}},
502
503 {},
504
505 {},
506
507 {})
508
509 );
510
511 pipeline_mng->getDomainRhsFE() = post_proc_fe;
512 CHKERR pipeline_mng->loopFiniteElements();
513 CHKERR post_proc_fe->writeFile("out_result.h5m");
514
516}
517//! [Output results]
518
519//! [Run program]
522
531
533}
534//! [Run program]
535
536//! [Main]
537int main(int argc, char *argv[]) {
538
539 // Initialisation of MoFEM/PETSc and MOAB data structures
540 const char param_file[] = "param_file.petsc";
542
543 // Error handling
544 try {
545 // Register MoFEM discrete manager in PETSc
546 DMType dm_name = "DMMOFEM";
547 CHKERR DMRegister_MoFEM(dm_name);
548
549 // Create MOAB instance
550 moab::Core mb_instance; // mesh database
551 moab::Interface &moab = mb_instance; // mesh database interface
552
553 // Create MoFEM instance
554 MoFEM::Core core(moab); // finite element database
555 MoFEM::Interface &m_field = core; // finite element interface
556
557 // Run the main analysis
558 Poisson2DiscontGalerkin poisson_problem(m_field);
559 CHKERR poisson_problem.runProgram();
560 }
562
563 // Finish work: cleaning memory, getting statistics, etc.
565
566 return 0;
567}
568//! [Main]
static Index< 'o', 3 > o
static Index< 'p', 3 > p
std::string param_file
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
int main()
Definition: adol-c_atom.cpp:46
static const double eps
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, FIELD_DIM > OpDomainSource
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:62
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ H1
continuous field
Definition: definitions.h:85
@ NOSPACE
Definition: definitions.h:83
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
ElementsAndOps< SPACE_DIM >::FaceSideEle FaceSideEle
@ F
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMoFEM.cpp:509
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMoFEM.cpp:412
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:47
auto createDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:1003
SmartPetscObj< KSP > createKSP(SmartPetscObj< DM > dm=nullptr)
Create KSP (linear) solver.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
FTensor::Index< 'i', SPACE_DIM > i
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradGrad< 1, 1, SPACE_DIM > OpDomainGradGrad
Definition: helmholtz.cpp:25
double D
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
PetscErrorCode PetscOptionsGetScalar(PetscOptions *, const char pre[], const char name[], PetscScalar *dval, PetscBool *set)
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
auto createVectorMPI(MPI_Comm comm, PetscInt n, PetscInt N)
Create MPI Vector.
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
auto u_grad_exact
static char help[]
constexpr int SPACE_DIM
static double nitsche
constexpr int FIELD_DIM
constexpr int BASE_DIM
DomainEle::UserDataOperator DomainEleOp
static double penalty
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradGrad< BASE_DIM, FIELD_DIM, SPACE_DIM > OpDomainGradGrad
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< BASE_DIM, FIELD_DIM > OpDomainSource
static double phi
BoundaryEle::UserDataOperator BoundaryEleOp
PetscBool is_test
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
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
base operator to do operations at Gauss Pt. level
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
@ OPROW
operator doWork function is executed on FE rows
@ OPSPACE
operator do Work is execute on space data
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Get value at integration points for scalar field.
Post post-proc data at points from hash maps.
PipelineManager interface.
boost::shared_ptr< FEMethod > & getDomainRhsFE()
boost::shared_ptr< FEMethod > & getDomainLhsFE()
keeps basic data about problem
DofIdx getNbDofsRow() const
Simple interface for fast problem set-up.
Definition: Simple.hpp:27
bool & getAddBoundaryFE()
Get the addSkeletonFE.
Definition: Simple.hpp:422
MoFEMErrorCode addDomainField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on domain.
Definition: Simple.cpp:264
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:180
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition: Simple.cpp:667
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:194
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:473
bool & getAddSkeletonFE()
Get the addSkeletonFE.
Definition: Simple.hpp:411
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:611
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
MoFEMErrorCode boundaryCondition()
[Setup problem]
MoFEMErrorCode checkResults()
[Solve system]
MoFEMErrorCode setupProblem()
[Read mesh]
MoFEMErrorCode setIntegrationRules()
[Assemble system]
MoFEMErrorCode assembleSystem()
[Boundary condition]
Poisson2DiscontGalerkin(MoFEM::Interface &m_field)
MoFEMErrorCode solveSystem()
[Set integration rules]
MoFEMErrorCode readMesh()
[Read mesh]
MoFEMErrorCode outputResults()
[Output results]
MoFEMErrorCode runProgram()
[Output results]
Operator tp collect data from elements on the side of Edge/Face.
Operator to evaluate Dirichlet boundary conditions using DG.