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 15 of file poisson_2d_lagrange_multiplier.cpp.

Constructor & Destructor Documentation

◆ Poisson2DLagrangeMultiplier()

Poisson2DLagrangeMultiplier::Poisson2DLagrangeMultiplier ( MoFEM::Interface m_field)

Definition at line 61 of file poisson_2d_lagrange_multiplier.cpp.

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

Member Function Documentation

◆ assembleSystem()

MoFEMErrorCode Poisson2DLagrangeMultiplier::assembleSystem ( )
private

Definition at line 134 of file poisson_2d_lagrange_multiplier.cpp.

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

◆ boundaryCondition()

MoFEMErrorCode Poisson2DLagrangeMultiplier::boundaryCondition ( )
private

Definition at line 93 of file poisson_2d_lagrange_multiplier.cpp.

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

◆ boundaryFunction()

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

Definition at line 39 of file poisson_2d_lagrange_multiplier.cpp.

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

◆ outputResults()

MoFEMErrorCode Poisson2DLagrangeMultiplier::outputResults ( )
private

Definition at line 249 of file poisson_2d_lagrange_multiplier.cpp.

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

◆ readMesh()

MoFEMErrorCode Poisson2DLagrangeMultiplier::readMesh ( )
private

◆ runProgram()

MoFEMErrorCode Poisson2DLagrangeMultiplier::runProgram ( )

Definition at line 286 of file poisson_2d_lagrange_multiplier.cpp.

286  {
288 
289  readMesh();
290  setupProblem();
292  assembleSystem();
294  solveSystem();
295  outputResults();
296 
298 }

◆ setIntegrationRules()

MoFEMErrorCode Poisson2DLagrangeMultiplier::setIntegrationRules ( )
private

Definition at line 179 of file poisson_2d_lagrange_multiplier.cpp.

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

◆ setupProblem()

MoFEMErrorCode Poisson2DLagrangeMultiplier::setupProblem ( )
private

◆ solveSystem()

MoFEMErrorCode Poisson2DLagrangeMultiplier::solveSystem ( )
private

Definition at line 197 of file poisson_2d_lagrange_multiplier.cpp.

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

◆ sourceTermFunction()

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

Definition at line 33 of file poisson_2d_lagrange_multiplier.cpp.

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

Member Data Documentation

◆ boundaryEntitiesForFieldsplit

Range Poisson2DLagrangeMultiplier::boundaryEntitiesForFieldsplit
private

Definition at line 58 of file poisson_2d_lagrange_multiplier.cpp.

◆ boundaryField

std::string Poisson2DLagrangeMultiplier::boundaryField
private

Definition at line 51 of file poisson_2d_lagrange_multiplier.cpp.

◆ boundaryMarker

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

Definition at line 55 of file poisson_2d_lagrange_multiplier.cpp.

◆ domainField

std::string Poisson2DLagrangeMultiplier::domainField
private

Definition at line 50 of file poisson_2d_lagrange_multiplier.cpp.

◆ mField

MoFEM::Interface& Poisson2DLagrangeMultiplier::mField
private

Definition at line 46 of file poisson_2d_lagrange_multiplier.cpp.

◆ oRder

int Poisson2DLagrangeMultiplier::oRder
private

Definition at line 52 of file poisson_2d_lagrange_multiplier.cpp.

◆ simpleInterface

Simple* Poisson2DLagrangeMultiplier::simpleInterface
private

Definition at line 47 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:33
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:134
Poisson2DLagrangeMultiplier::boundaryField
std::string boundaryField
Definition: poisson_2d_lagrange_multiplier.cpp:51
Poisson2DLagrangeMultiplier::solveSystem
MoFEMErrorCode solveSystem()
Definition: poisson_2d_lagrange_multiplier.cpp:197
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:693
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:65
Poisson2DLagrangeMultiplier::outputResults
MoFEMErrorCode outputResults()
Definition: poisson_2d_lagrange_multiplier.cpp:249
Poisson2DLagrangeMultiplier::setIntegrationRules
MoFEMErrorCode setIntegrationRules()
Definition: poisson_2d_lagrange_multiplier.cpp:179
MoFEM::ISManager
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:23
Poisson2DLagrangeMultiplierOperators::OpBoundaryLhsC
Definition: poisson_2d_lagrange_multiplier.hpp:169
OpPPMap
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
Definition: photon_diffusion.cpp:29
Poisson2DLagrangeMultiplier::boundaryEntitiesForFieldsplit
Range boundaryEntitiesForFieldsplit
Definition: poisson_2d_lagrange_multiplier.cpp:58
Poisson2DLagrangeMultiplier::boundaryFunction
static double boundaryFunction(const double x, const double y, const double z)
Definition: poisson_2d_lagrange_multiplier.cpp:39
Poisson2DLagrangeMultiplier::boundaryCondition
MoFEMErrorCode boundaryCondition()
Definition: poisson_2d_lagrange_multiplier.cpp:93
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:47
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:3200
Poisson2DLagrangeMultiplier::setupProblem
MoFEMErrorCode setupProblem()
Definition: poisson_2d_lagrange_multiplier.cpp:75
Poisson2DLagrangeMultiplier::domainField
std::string domainField
Definition: poisson_2d_lagrange_multiplier.cpp:50
Poisson2DLagrangeMultiplierOperators::OpDomainRhsF
Definition: poisson_2d_lagrange_multiplier.hpp:104
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:491
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:300
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:52
MoFEM::Simple::setUp
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:629
MoFEM::Problem::getName
auto getName() const
Definition: ProblemsMultiIndices.hpp:372
Poisson2DLagrangeMultiplier::mField
MoFEM::Interface & mField
Definition: poisson_2d_lagrange_multiplier.cpp:46
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:236
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:23
Poisson2DLagrangeMultiplier::boundaryMarker
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
Definition: poisson_2d_lagrange_multiplier.cpp:55
F
@ F
Definition: free_surface.cpp:394
MoFEM::OpPostProcMapInMoab
Post post-proc data at points from hash maps.
Definition: PostProcBrokenMeshInMoabBase.hpp:698