v0.14.0
Public Member Functions | Private Member Functions | Static Private Member Functions | Private Attributes | List of all members
Poisson2DLagrangeMultiplier Struct Reference
Collaboration diagram for Poisson2DLagrangeMultiplier:
[legend]

Public Member Functions

 Poisson2DLagrangeMultiplier (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProgram ()
 

Private Member Functions

MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode setIntegrationRules ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 

Static Private Member Functions

static double sourceTermFunction (const double x, const double y, const double z)
 
static double boundaryFunction (const double x, const double y, const double z)
 

Private Attributes

MoFEM::InterfacemField
 
SimplesimpleInterface
 
std::string domainField
 
std::string boundaryField
 
int oRder
 
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
 
Range boundaryEntitiesForFieldsplit
 

Detailed Description

Definition at line 16 of file poisson_2d_lagrange_multiplier.cpp.

Constructor & Destructor Documentation

◆ Poisson2DLagrangeMultiplier()

Poisson2DLagrangeMultiplier::Poisson2DLagrangeMultiplier ( MoFEM::Interface m_field)

Definition at line 62 of file poisson_2d_lagrange_multiplier.cpp.

64  : domainField("U"), boundaryField("L"), mField(m_field) {}

Member Function Documentation

◆ assembleSystem()

MoFEMErrorCode Poisson2DLagrangeMultiplier::assembleSystem ( )
private

Definition at line 135 of file poisson_2d_lagrange_multiplier.cpp.

135  {
137 
138  auto pipeline_mng = mField.getInterface<PipelineManager>();
139  auto det_ptr = boost::make_shared<VectorDouble>();
140  auto jac_ptr = boost::make_shared<MatrixDouble>();
141  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
142 
143  { // Push operators to the Pipeline that is responsible for calculating LHS of
144  // domain elements
145  pipeline_mng->getOpDomainLhsPipeline().push_back(
146  new OpCalculateHOJac<2>(jac_ptr));
147  pipeline_mng->getOpDomainLhsPipeline().push_back(
148  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
149  pipeline_mng->getOpDomainLhsPipeline().push_back(
150  new OpSetHOInvJacToScalarBases<2>(H1, inv_jac_ptr));
151  pipeline_mng->getOpDomainLhsPipeline().push_back(
152  new OpSetHOWeightsOnFace());
153  pipeline_mng->getOpDomainLhsPipeline().push_back(
155  }
156 
157  { // Push operators to the Pipeline that is responsible for calculating RHS of
158  // domain elements
159  pipeline_mng->getOpDomainRhsPipeline().push_back(
161  }
162 
163  { // Push operators to the Pipeline that is responsible for calculating LHS of
164  // boundary elements (Lagrange multiplier)
165  pipeline_mng->getOpBoundaryLhsPipeline().push_back(
167  pipeline_mng->getOpBoundaryLhsPipeline().push_back(
169  }
170 
171  { // Push operators to the Pipeline that is responsible for calculating RHS of
172  // boundary elements (Lagrange multiplier)
173  pipeline_mng->getOpBoundaryRhsPipeline().push_back(
175  }
176 
178 }

◆ boundaryCondition()

MoFEMErrorCode Poisson2DLagrangeMultiplier::boundaryCondition ( )
private

Definition at line 94 of file poisson_2d_lagrange_multiplier.cpp.

94  {
96 
97  // Get boundary edges marked in block named "BOUNDARY_CONDITION"
98  auto get_ents_on_mesh = [&]() {
99  Range boundary_entities;
101  std::string entity_name = it->getName();
102  if (entity_name.compare(0, 18, "BOUNDARY_CONDITION") == 0) {
103  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 1,
104  boundary_entities, true);
105  }
106  }
107  // Add vertices to boundary entities
108  Range boundary_vertices;
109  CHKERR mField.get_moab().get_connectivity(boundary_entities,
110  boundary_vertices, true);
111  boundary_entities.merge(boundary_vertices);
112 
113  // Store entities for fieldsplit (block) solver
114  boundaryEntitiesForFieldsplit = boundary_entities;
115 
116  return boundary_entities;
117  };
118 
119  auto mark_boundary_dofs = [&](Range &&skin_edges) {
120  auto problem_manager = mField.getInterface<ProblemsManager>();
121  auto marker_ptr = boost::make_shared<std::vector<unsigned char>>();
122  problem_manager->markDofs(simpleInterface->getProblemName(), ROW,
123  skin_edges, *marker_ptr);
124  return marker_ptr;
125  };
126 
127  // Get global local vector of marked DOFs. Is global, since is set for all
128  // DOFs on processor. Is local since only DOFs on processor are in the
129  // vector. To access DOFs use local indices.
130  boundaryMarker = mark_boundary_dofs(get_ents_on_mesh());
131 
133 }

◆ boundaryFunction()

static double Poisson2DLagrangeMultiplier::boundaryFunction ( const double  x,
const double  y,
const double  z 
)
inlinestaticprivate

Definition at line 40 of file poisson_2d_lagrange_multiplier.cpp.

41  {
42  return sin(x * 10.) * cos(y * 10.);
43  // return 0;
44  }

◆ outputResults()

MoFEMErrorCode Poisson2DLagrangeMultiplier::outputResults ( )
private

Definition at line 250 of file poisson_2d_lagrange_multiplier.cpp.

250  {
252 
253  auto pipeline_mng = mField.getInterface<PipelineManager>();
254  pipeline_mng->getDomainLhsFE().reset();
255  pipeline_mng->getBoundaryLhsFE().reset();
256 
257  auto d_ptr = boost::make_shared<VectorDouble>();
258  auto l_ptr = boost::make_shared<VectorDouble>();
259 
261 
262  auto post_proc_domain_fe = boost::make_shared<PostProcFaceEle>(mField);
263  post_proc_domain_fe->getOpPtrVector().push_back(
265  post_proc_domain_fe->getOpPtrVector().push_back(
266  new OpPPMap(post_proc_domain_fe->getPostProcMesh(),
267  post_proc_domain_fe->getMapGaussPts(), {{domainField, d_ptr}},
268  {}, {}, {}));
269  pipeline_mng->getDomainRhsFE() = post_proc_domain_fe;
270 
271  auto post_proc_boundary_fe = boost::make_shared<PostProcEdgeEle>(mField);
272  post_proc_boundary_fe->getOpPtrVector().push_back(
274  post_proc_boundary_fe->getOpPtrVector().push_back(
275  new OpPPMap(post_proc_boundary_fe->getPostProcMesh(),
276  post_proc_boundary_fe->getMapGaussPts(),
277  {{boundaryField, l_ptr}}, {}, {}, {}));
278  pipeline_mng->getBoundaryRhsFE() = post_proc_boundary_fe;
279 
280  CHKERR pipeline_mng->loopFiniteElements();
281  CHKERR post_proc_domain_fe->writeFile("out_result_domain.h5m");
282  CHKERR post_proc_boundary_fe->writeFile("out_result_boundary.h5m");
283 
285 }

◆ readMesh()

MoFEMErrorCode Poisson2DLagrangeMultiplier::readMesh ( )
private

◆ runProgram()

MoFEMErrorCode Poisson2DLagrangeMultiplier::runProgram ( )

Definition at line 287 of file poisson_2d_lagrange_multiplier.cpp.

287  {
289 
290  readMesh();
291  setupProblem();
293  assembleSystem();
295  solveSystem();
296  outputResults();
297 
299 }

◆ setIntegrationRules()

MoFEMErrorCode Poisson2DLagrangeMultiplier::setIntegrationRules ( )
private

Definition at line 180 of file poisson_2d_lagrange_multiplier.cpp.

180  {
182 
183  auto pipeline_mng = mField.getInterface<PipelineManager>();
184 
185  auto domain_rule_lhs = [](int, int, int p) -> int { return 2 * (p - 1); };
186  auto domain_rule_rhs = [](int, int, int p) -> int { return 2 * (p - 1); };
187  CHKERR pipeline_mng->setDomainLhsIntegrationRule(domain_rule_lhs);
188  CHKERR pipeline_mng->setDomainRhsIntegrationRule(domain_rule_rhs);
189 
190  auto boundary_rule_lhs = [](int, int, int p) -> int { return 2 * p; };
191  auto boundary_rule_rhs = [](int, int, int p) -> int { return 2 * p; };
192  CHKERR pipeline_mng->setBoundaryLhsIntegrationRule(boundary_rule_lhs);
193  CHKERR pipeline_mng->setBoundaryLhsIntegrationRule(boundary_rule_rhs);
194 
196 }

◆ setupProblem()

MoFEMErrorCode Poisson2DLagrangeMultiplier::setupProblem ( )
private

◆ solveSystem()

MoFEMErrorCode Poisson2DLagrangeMultiplier::solveSystem ( )
private

Definition at line 198 of file poisson_2d_lagrange_multiplier.cpp.

198  {
200 
201  auto pipeline_mng = mField.getInterface<PipelineManager>();
202 
203  auto ksp_solver = pipeline_mng->createKSP();
204  CHKERR KSPSetFromOptions(ksp_solver);
205 
206  // Create RHS and solution vectors
207  auto dm = simpleInterface->getDM();
208  auto F = createDMVector(dm);
209  auto D = vectorDuplicate(F);
210 
211  // Setup fieldsplit (block) solver - optional: yes/no
212  if (1) {
213  PC pc;
214  CHKERR KSPGetPC(ksp_solver, &pc);
215  PetscBool is_pcfs = PETSC_FALSE;
216  PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
217 
218  // Set up FIELDSPLIT, only when user set -pc_type fieldsplit
219  // Identify the index for boundary entities, remaining will be for domain
220  // Then split the fields for boundary and domain for solving
221  if (is_pcfs == PETSC_TRUE) {
222  IS is_domain, is_boundary;
223  cerr << "Running FIELDSPLIT..." << endl;
224  const MoFEM::Problem *problem_ptr;
225  CHKERR DMMoFEMGetProblemPtr(dm, &problem_ptr);
226  CHKERR mField.getInterface<ISManager>()->isCreateProblemFieldAndRank(
227  problem_ptr->getName(), ROW, domainField, 0, 1, &is_boundary,
229  // CHKERR ISView(is_boundary, PETSC_VIEWER_STDOUT_SELF);
230 
231  CHKERR PCFieldSplitSetIS(pc, NULL, is_boundary);
232 
233  CHKERR ISDestroy(&is_boundary);
234  }
235  }
236 
237  CHKERR KSPSetUp(ksp_solver);
238 
239  // Solve the system
240  CHKERR KSPSolve(ksp_solver, F, D);
241 
242  // Scatter result data on the mesh
243  CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
244  CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
245  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
246 
248 }

◆ sourceTermFunction()

static double Poisson2DLagrangeMultiplier::sourceTermFunction ( const double  x,
const double  y,
const double  z 
)
inlinestaticprivate

Definition at line 34 of file poisson_2d_lagrange_multiplier.cpp.

35  {
36  return 200 * sin(x * 10.) * cos(y * 10.);
37  // return 1;
38  }

Member Data Documentation

◆ boundaryEntitiesForFieldsplit

Range Poisson2DLagrangeMultiplier::boundaryEntitiesForFieldsplit
private

Definition at line 59 of file poisson_2d_lagrange_multiplier.cpp.

◆ boundaryField

std::string Poisson2DLagrangeMultiplier::boundaryField
private

Definition at line 52 of file poisson_2d_lagrange_multiplier.cpp.

◆ boundaryMarker

boost::shared_ptr<std::vector<unsigned char> > Poisson2DLagrangeMultiplier::boundaryMarker
private

Definition at line 56 of file poisson_2d_lagrange_multiplier.cpp.

◆ domainField

std::string Poisson2DLagrangeMultiplier::domainField
private

Definition at line 51 of file poisson_2d_lagrange_multiplier.cpp.

◆ mField

MoFEM::Interface& Poisson2DLagrangeMultiplier::mField
private

Definition at line 47 of file poisson_2d_lagrange_multiplier.cpp.

◆ oRder

int Poisson2DLagrangeMultiplier::oRder
private

Definition at line 53 of file poisson_2d_lagrange_multiplier.cpp.

◆ simpleInterface

Simple* Poisson2DLagrangeMultiplier::simpleInterface
private

Definition at line 48 of file poisson_2d_lagrange_multiplier.cpp.


The documentation for this struct was generated from the following file:
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
Poisson2DLagrangeMultiplier::sourceTermFunction
static double sourceTermFunction(const double x, const double y, const double z)
Definition: poisson_2d_lagrange_multiplier.cpp:34
H1
@ H1
continuous field
Definition: definitions.h:85
MoFEM::ProblemsManager
Problem manager is used to build and partition problems.
Definition: ProblemsManager.hpp:21
Poisson2DLagrangeMultiplier::assembleSystem
MoFEMErrorCode assembleSystem()
Definition: poisson_2d_lagrange_multiplier.cpp:135
Poisson2DLagrangeMultiplier::boundaryField
std::string boundaryField
Definition: poisson_2d_lagrange_multiplier.cpp:52
Poisson2DLagrangeMultiplier::solveSystem
MoFEMErrorCode solveSystem()
Definition: poisson_2d_lagrange_multiplier.cpp:198
MoFEM::OpSetHOInvJacToScalarBases< 2 >
Definition: HODataOperators.hpp:78
MoFEM::Simple::loadFile
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
Definition: Simple.cpp:194
MoFEM::PipelineManager
PipelineManager interface.
Definition: PipelineManager.hpp:24
MoFEM::DMoFEMMeshToLocalVector
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:523
MoFEM::OpCalculateHOJac< 2 >
Definition: HODataOperators.hpp:273
ROW
@ ROW
Definition: definitions.h:136
MoFEM::Simple::getOptions
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:180
MoFEM::Simple::getDM
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition: Simple.cpp:747
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::createDMVector
auto createDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:1099
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
Poisson2DLagrangeMultiplier::readMesh
MoFEMErrorCode readMesh()
Definition: poisson_2d_lagrange_multiplier.cpp:66
Poisson2DLagrangeMultiplier::outputResults
MoFEMErrorCode outputResults()
Definition: poisson_2d_lagrange_multiplier.cpp:250
Poisson2DLagrangeMultiplier::setIntegrationRules
MoFEMErrorCode setIntegrationRules()
Definition: poisson_2d_lagrange_multiplier.cpp:180
MoFEM::ISManager
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:23
Poisson2DLagrangeMultiplierOperators::OpBoundaryLhsC
Definition: poisson_2d_lagrange_multiplier.hpp:171
OpPPMap
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
Definition: photon_diffusion.cpp:29
Poisson2DLagrangeMultiplier::boundaryEntitiesForFieldsplit
Range boundaryEntitiesForFieldsplit
Definition: poisson_2d_lagrange_multiplier.cpp:59
Poisson2DLagrangeMultiplier::boundaryFunction
static double boundaryFunction(const double x, const double y, const double z)
Definition: poisson_2d_lagrange_multiplier.cpp:40
Poisson2DLagrangeMultiplier::boundaryCondition
MoFEMErrorCode boundaryCondition()
Definition: poisson_2d_lagrange_multiplier.cpp:94
MoFEM::OpCalculateScalarFieldValues
Get value at integration points for scalar field.
Definition: UserDataOperators.hpp:82
MoFEM::Simple::addDomainField
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
MoFEM::PipelineManager::createKSP
SmartPetscObj< KSP > createKSP(SmartPetscObj< DM > dm=nullptr)
Create KSP (linear) solver.
Definition: PipelineManager.cpp:49
Poisson2DLagrangeMultiplier::simpleInterface
Simple * simpleInterface
Definition: poisson_2d_lagrange_multiplier.cpp:48
MoFEM::OpSetHOWeightsOnFace
Modify integration weights on face to take in account higher-order geometry.
Definition: HODataOperators.hpp:122
MoFEM::ProblemsManager::markDofs
MoFEMErrorCode markDofs(const std::string problem_name, RowColData rc, const enum MarkOP op, const Range ents, std::vector< unsigned char > &marker) const
Create vector with marked indices.
Definition: ProblemsManager.cpp:3548
Poisson2DLagrangeMultiplier::setupProblem
MoFEMErrorCode setupProblem()
Definition: poisson_2d_lagrange_multiplier.cpp:76
Poisson2DLagrangeMultiplier::domainField
std::string domainField
Definition: poisson_2d_lagrange_multiplier.cpp:51
Poisson2DLagrangeMultiplierOperators::OpDomainRhsF
Definition: poisson_2d_lagrange_multiplier.hpp:106
AINSWORTH_BERNSTEIN_BEZIER_BASE
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:64
MoFEM::Simple::setFieldOrder
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:545
Range
MoFEM::PipelineManager::getDomainLhsFE
boost::shared_ptr< FEMethod > & getDomainLhsFE()
Definition: PipelineManager.hpp:401
MoFEM::vectorDuplicate
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
Definition: PetscSmartObj.hpp:221
_IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
Definition: MeshsetsManager.hpp:71
BLOCKSET
@ BLOCKSET
Definition: definitions.h:161
MoFEM::DMMoFEMGetProblemPtr
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMoFEM.cpp:426
MoFEM::OpInvertMatrix
Definition: UserDataOperators.hpp:3249
MoFEM::Simple::addBoundaryField
MoFEMErrorCode addBoundaryField(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 boundary.
Definition: Simple.cpp:354
ReactionDiffusionEquation::D
const double D
diffusivity
Definition: reaction_diffusion.cpp:20
MoFEM::Problem
keeps basic data about problem
Definition: ProblemsMultiIndices.hpp:54
Poisson2DLagrangeMultiplier::oRder
int oRder
Definition: poisson_2d_lagrange_multiplier.cpp:53
MoFEM::Simple::setUp
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:683
MoFEM::Problem::getName
auto getName() const
Definition: ProblemsMultiIndices.hpp:372
Poisson2DLagrangeMultiplier::mField
MoFEM::Interface & mField
Definition: poisson_2d_lagrange_multiplier.cpp:47
MoFEM::Simple::getProblemName
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:389
convert.int
int
Definition: convert.py:64
MoFEM::PetscOptionsGetInt
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
Definition: DeprecatedPetsc.hpp:142
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
Poisson2DLagrangeMultiplierOperators::OpBoundaryRhsG
Definition: poisson_2d_lagrange_multiplier.hpp:238
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
Poisson2DLagrangeMultiplierOperators::OpDomainLhsK
Definition: poisson_2d_lagrange_multiplier.hpp:25
Poisson2DLagrangeMultiplier::boundaryMarker
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
Definition: poisson_2d_lagrange_multiplier.cpp:56
F
@ F
Definition: free_surface.cpp:394
MoFEM::OpPostProcMapInMoab
Post post-proc data at points from hash maps.
Definition: PostProcBrokenMeshInMoabBase.hpp:698