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

Public Member Functions

 Poisson (moab::Core &mb_instance, MoFEM::Core &core, int order)
 
MoFEMErrorCode runAnalysis ()
 

Private Member Functions

MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupSystem ()
 
MoFEMErrorCode setIntegrationRule ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode applyBc ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode setOutput ()
 

Private Attributes

boost::shared_ptr< VolElevolPipelineLhs
 
boost::shared_ptr< VolElevolPipelineRhs
 
boost::shared_ptr< VolElepostFaceVolume
 
MoFEM::InterfacemField
 
moab::Interface & moab
 
SimplesimpleInterface
 
DM dm
 
KSP ksp
 
MatrixDouble invJac
 
MPI_Comm mpiComm
 
const int mpiRank
 
int order
 
std::string domainField
 

Detailed Description

Definition at line 10 of file poisson_prob.cpp.

Constructor & Destructor Documentation

◆ Poisson()

Poisson::Poisson ( moab::Core &  mb_instance,
MoFEM::Core core,
int  order 
)

Definition at line 46 of file poisson_prob.cpp.

47  : domainField("U"), moab(mb_instance), mField(core),
49  order(order) {
50  volPipelineLhs = boost::shared_ptr<VolEle>(new VolEle(mField));
51  volPipelineRhs = boost::shared_ptr<VolEle>(new VolEle(mField));
52 }
const int mpiRank
MoFEM::VolumeElementForcesAndSourcesCore VolEle
MoFEM::Interface & mField
MPI_Comm mpiComm
std::string domainField
virtual int get_comm_rank() const =0
boost::shared_ptr< VolEle > volPipelineLhs
moab::Interface & moab
boost::shared_ptr< VolEle > volPipelineRhs
virtual MPI_Comm & get_comm() const =0

Member Function Documentation

◆ applyBc()

MoFEMErrorCode Poisson::applyBc ( )
private

Definition at line 168 of file poisson_prob.cpp.

168  {
170  Range bdry_entities;
172  string name = it->getName();
173  if (name.compare(0, 14, "ESSENTIAL") == 0) {
174  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 1,
175  bdry_entities, true);
176  }
177  }
178 
179  cerr << bdry_entities;
180  Range bdr_verts;
181  CHKERR mField.get_moab().get_connectivity(bdry_entities, bdr_verts, true);
182  bdry_entities.merge(bdr_verts);
183 
184  CHKERR mField.getInterface<ProblemsManager>()->removeDofsOnEntities(
185  simpleInterface->getProblemName(), domainField, bdry_entities);
186 
187 
189 }
Problem manager is used to build and partition problems.
virtual moab::Interface & get_moab()=0
MoFEM::Interface & mField
Simple * simpleInterface
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::string domainField
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:302
#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.
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ assembleSystem()

MoFEMErrorCode Poisson::assembleSystem ( )
private

Definition at line 84 of file poisson_prob.cpp.

84  {
86 
87  volPipelineLhs->getOpPtrVector().push_back(new OpCalculateInvJacForFace(invJac));
88  volPipelineLhs->getOpPtrVector().push_back(new OpSetInvJacH1ForFace(invJac));
89 
90  volPipelineLhs->getOpPtrVector().push_back(
92  volPipelineRhs->getOpPtrVector().push_back(new OpRhsU(domainField));
93 
94  // dm = simpleInterface->getDM();
96  boost::shared_ptr<VolEle> null;
97 
99 
103 }
MatrixDouble invJac
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition: Simple.cpp:706
PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set compute operator for KSP solver via sub-matrix and IS.
Definition: DMMMoFEM.cpp:556
Simple * simpleInterface
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set KSP operators and push mofem finite element methods.
Definition: DMMMoFEM.cpp:597
std::string domainField
Calculate inverse of jacobian for face element.
boost::shared_ptr< VolEle > volPipelineLhs
#define CHKERR
Inline error check.
Definition: definitions.h:602
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition: Simple.hpp:281
Transform local reference derivatives of shape functions to global derivatives.
boost::shared_ptr< VolEle > volPipelineRhs
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
MoFEMErrorCode setIntegrationRule()

◆ readMesh()

MoFEMErrorCode Poisson::readMesh ( )
private

Definition at line 54 of file poisson_prob.cpp.

54  {
56 
60 
62 }
MoFEM::Interface & mField
Simple * simpleInterface
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:212
#define CHKERR
Inline error check.
Definition: definitions.h:602
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:198
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ runAnalysis()

MoFEMErrorCode Poisson::runAnalysis ( )

Definition at line 191 of file poisson_prob.cpp.

191  {
193 
194  CHKERR readMesh();
196  CHKERR applyBc();
197 
199 
200  CHKERR setOutput();
201 
203 
205 }
MoFEMErrorCode setupSystem()
MoFEMErrorCode setOutput()
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
MoFEMErrorCode solveSystem()
MoFEMErrorCode applyBc()
#define CHKERR
Inline error check.
Definition: definitions.h:602
MoFEMErrorCode readMesh()
MoFEMErrorCode assembleSystem()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ setIntegrationRule()

MoFEMErrorCode Poisson::setIntegrationRule ( )
private

Definition at line 75 of file poisson_prob.cpp.

75  {
77  auto vol_ruleLhs = [](int, int, int p) -> int { return 2 * (p - 1); };
78  auto vol_ruleRhs = [](int, int, int p) -> int { return 2 * (p - 1); };
79  volPipelineLhs->getRuleHook = vol_ruleLhs;
80  volPipelineRhs->getRuleHook = vol_ruleRhs;
82 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
boost::shared_ptr< VolEle > volPipelineLhs
boost::shared_ptr< VolEle > volPipelineRhs
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ setOutput()

MoFEMErrorCode Poisson::setOutput ( )
private

Definition at line 149 of file poisson_prob.cpp.

149  {
151 
153  boost::shared_ptr<VolEle>(new PostProcFaceOnRefinedMesh(mField));
154 
155  CHKERR
156  boost::static_pointer_cast<PostProcFaceOnRefinedMesh>(postFaceVolume)
157  ->generateReferenceElementMesh();
158  CHKERR
159  boost::static_pointer_cast<PostProcFaceOnRefinedMesh>(postFaceVolume)
160  ->addFieldValuesPostProc(domainField);
161  // CHKERR
162  // boost::static_pointer_cast<PostProcFaceOnRefinedMesh>(postFaceVolume)
163  // ->addFieldValuesGradientPostProc(domainField);
164 
166 }
MoFEM::Interface & mField
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::string domainField
#define CHKERR
Inline error check.
Definition: definitions.h:602
Postprocess on face.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
boost::shared_ptr< VolEle > postFaceVolume

◆ setupSystem()

MoFEMErrorCode Poisson::setupSystem ( )
private

Definition at line 64 of file poisson_prob.cpp.

64  {
66 
71 
73 }
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:672
Simple * simpleInterface
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::string domainField
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
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:269
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
continuous field
Definition: definitions.h:177
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:482

◆ solveSystem()

MoFEMErrorCode Poisson::solveSystem ( )
private

Definition at line 105 of file poisson_prob.cpp.

105  {
107 
108  Vec global_rhs, global_solution;
109  CHKERR DMCreateGlobalVector(dm, &global_rhs);
110  CHKERR VecDuplicate(global_rhs, &global_solution);
111 
112 
113  // ksp = createKSP(mField.get_comm());
114 
115  CHKERR KSPCreate(PETSC_COMM_WORLD, &ksp);
116 
117  // CHKERR KSPSetOperators(ksp, volPipelineLhs->getFEMethod()->ksp_B,
118  // volPipelineLhs->getFEMethod()->ksp_B);
119 
120  CHKERR KSPSetFromOptions(ksp);
121 
122  CHKERR KSPSetDM(ksp, dm);
123 
124  CHKERR KSPSetUp(ksp);
125 
126  CHKERR KSPSolve(ksp, global_rhs, global_solution);
127 
128  CHKERR DMoFEMMeshToGlobalVector(dm, global_solution, INSERT_VALUES,
129  SCATTER_REVERSE);
130 
131 
132 
133  CHKERR KSPDestroy(&ksp);
134  CHKERR VecDestroy(&global_solution);
135  CHKERR VecDestroy(&global_rhs);
136 
139 
140  CHKERR
141  boost::static_pointer_cast<PostProcFaceOnRefinedMesh>(postFaceVolume)
142  ->writeFile("out_vol.h5m");
143 
144  CHKERR DMDestroy(&dm);
145 
147 }
Simple * simpleInterface
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:507
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:457
#define CHKERR
Inline error check.
Definition: definitions.h:602
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition: Simple.hpp:281
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
boost::shared_ptr< VolEle > postFaceVolume

Member Data Documentation

◆ dm

DM Poisson::dm
private

Definition at line 33 of file poisson_prob.cpp.

◆ domainField

std::string Poisson::domainField
private

Definition at line 43 of file poisson_prob.cpp.

◆ invJac

MatrixDouble Poisson::invJac
private

Definition at line 36 of file poisson_prob.cpp.

◆ ksp

KSP Poisson::ksp
private

Definition at line 34 of file poisson_prob.cpp.

◆ mField

MoFEM::Interface& Poisson::mField
private

Definition at line 29 of file poisson_prob.cpp.

◆ moab

moab::Interface& Poisson::moab
private

Definition at line 30 of file poisson_prob.cpp.

◆ mpiComm

MPI_Comm Poisson::mpiComm
private

Definition at line 38 of file poisson_prob.cpp.

◆ mpiRank

const int Poisson::mpiRank
private

Definition at line 39 of file poisson_prob.cpp.

◆ order

int Poisson::order
private

Definition at line 41 of file poisson_prob.cpp.

◆ postFaceVolume

boost::shared_ptr<VolEle> Poisson::postFaceVolume
private

Definition at line 27 of file poisson_prob.cpp.

◆ simpleInterface

Simple* Poisson::simpleInterface
private

Definition at line 31 of file poisson_prob.cpp.

◆ volPipelineLhs

boost::shared_ptr<VolEle> Poisson::volPipelineLhs
private

Definition at line 25 of file poisson_prob.cpp.

◆ volPipelineRhs

boost::shared_ptr<VolEle> Poisson::volPipelineRhs
private

Definition at line 26 of file poisson_prob.cpp.


The documentation for this struct was generated from the following file: