v0.14.0
Classes | Functions | Variables
simple_interface.cpp File Reference
#include <MoFEM.hpp>

Go to the source code of this file.

Classes

struct  OpVolume
 
struct  OpFace
 
struct  VolRule
 Set integration rule to volume elements. More...
 
struct  FaceRule
 Set integration rule to boundary elements. More...
 

Functions

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

Variables

static char help [] = "...\n\n"
 

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)
Examples
simple_interface.cpp.

Definition at line 95 of file simple_interface.cpp.

95  {
96 
97  // initialize petsc
98  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
99 
100  try {
101 
102  // Create MoAB database
103  moab::Core moab_core;
104  moab::Interface &moab = moab_core;
105 
106  // Create MoFEM database and link it to MoAB
107  MoFEM::Core mofem_core(moab);
108  MoFEM::Interface &m_field = mofem_core;
109 
110  // Register DM Manager
111  DMType dm_name = "DMMOFEM";
112  CHKERR DMRegister_MoFEM(dm_name);
113 
114  // Simple interface
115  Simple *simple_interface;
116  CHKERR m_field.getInterface(simple_interface);
117  {
118  // get options from command line
119  CHKERR simple_interface->getOptions();
120  // load mesh file
121  CHKERR simple_interface->loadFile();
122  // add fields
123  CHKERR simple_interface->addDomainField("MESH_NODE_POSITIONS", H1,
125  CHKERR simple_interface->addBoundaryField("MESH_NODE_POSITIONS", H1,
127  CHKERR simple_interface->addSkeletonField("MESH_NODE_POSITIONS", H1,
129 
130  // set fields order
131  CHKERR simple_interface->setFieldOrder("MESH_NODE_POSITIONS", 1);
132  // setup problem
133  CHKERR simple_interface->setUp();
134  // Project mesh coordinate on mesh
135  Projection10NodeCoordsOnField ent_method(m_field, "MESH_NODE_POSITIONS");
136  CHKERR m_field.loop_dofs("MESH_NODE_POSITIONS", ent_method);
137  // get dm
138  auto dm = simple_interface->getDM();
139  // create elements
140  auto domain_fe =
141  boost::make_shared<VolumeElementForcesAndSourcesCore>(m_field);
142  auto boundary_fe =
143  boost::make_shared<FaceElementForcesAndSourcesCore>(m_field);
144 
145  // set integration rule
146  domain_fe->getRuleHook = VolRule();
147  boundary_fe->getRuleHook = FaceRule();
148  // create distributed vector to accumulate values from processors.
149  int ghosts[] = {0};
150 
151  auto vol = createGhostVector(
152  PETSC_COMM_WORLD, m_field.get_comm_rank() == 0 ? 1 : 0, 1,
153  m_field.get_comm_rank() == 0 ? 0 : 1, ghosts);
154  auto surf_vol = vectorDuplicate(vol);
155 
156  // set operator to the volume element
157  auto material_grad_mat = boost::make_shared<MatrixDouble>();
158  auto material_det_vec = boost::make_shared<VectorDouble>();
159 
160  domain_fe->getOpPtrVector().push_back(
161  new OpCalculateVectorFieldGradient<3, 3>("MESH_NODE_POSITIONS",
162  material_grad_mat));
163  domain_fe->getOpPtrVector().push_back(new OpInvertMatrix<3>(
164  material_grad_mat, material_det_vec, nullptr));
165  domain_fe->getOpPtrVector().push_back(
166  new OpSetHOWeights(material_det_vec));
167  domain_fe->getOpPtrVector().push_back(
168  new OpCalculateHOCoords("MESH_NODE_POSITIONS"));
169  domain_fe->getOpPtrVector().push_back(
170  new OpVolume("MESH_NODE_POSITIONS", vol));
171  // set operator to the face element
172  boundary_fe->getOpPtrVector().push_back(
173  new OpCalculateHOCoords("MESH_NODE_POSITIONS"));
174  boundary_fe->getOpPtrVector().push_back(
175  new OpFace("MESH_NODE_POSITIONS", surf_vol));
176  // make integration in volume (here real calculations starts)
177  CHKERR DMoFEMLoopFiniteElements(dm, simple_interface->getDomainFEName(),
178  domain_fe);
179  // make integration on boundary
180  CHKERR DMoFEMLoopFiniteElements(dm, simple_interface->getBoundaryFEName(),
181  boundary_fe);
182 
183  auto skeleton_fe = boost::make_shared<FEMethod>();
184  auto A = createDMMatrix(dm);
185  auto B = createDMMatrix(dm);
186  auto f = createDMVector(dm);
187  auto x = createDMVector(dm);
188  auto x_t = createDMVector(dm);
189  auto x_tt = createDMVector(dm);
190  skeleton_fe->f = f;
191  skeleton_fe->A = A;
192  skeleton_fe->B = B;
193  skeleton_fe->x = x;
194  skeleton_fe->x_t = x_t;
195  skeleton_fe->x_tt = x_tt;
196 
197  skeleton_fe->preProcessHook = [&]() {
199  if (f != skeleton_fe->ts_F)
200  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong ptr");
201  if (A != skeleton_fe->ts_A)
202  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong ptr");
203  if (B != skeleton_fe->ts_B)
204  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong ptr");
205  if (x != skeleton_fe->ts_u)
206  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong ptr");
207  if (x_t != skeleton_fe->ts_u_t)
208  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong ptr");
209  if (x_tt != skeleton_fe->ts_u_tt)
210  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong ptr");
212  };
213 
214  skeleton_fe->postProcessHook = []() { return 0; };
215  skeleton_fe->operatorHook = []() { return 0; };
216 
217  CHKERR DMoFEMLoopFiniteElements(dm, simple_interface->getSkeletonFEName(),
218  skeleton_fe);
219 
220  // assemble volumes from processors and accumulate on processor of rank 0
221  CHKERR VecAssemblyBegin(vol);
222  CHKERR VecAssemblyEnd(vol);
223  CHKERR VecGhostUpdateBegin(vol, ADD_VALUES, SCATTER_REVERSE);
224  CHKERR VecGhostUpdateEnd(vol, ADD_VALUES, SCATTER_REVERSE);
225  CHKERR VecAssemblyBegin(surf_vol);
226  CHKERR VecAssemblyEnd(surf_vol);
227  CHKERR VecGhostUpdateBegin(surf_vol, ADD_VALUES, SCATTER_REVERSE);
228  CHKERR VecGhostUpdateEnd(surf_vol, ADD_VALUES, SCATTER_REVERSE);
229  if (m_field.get_comm_rank() == 0) {
230  double *a_vol;
231  CHKERR VecGetArray(vol, &a_vol);
232  double *a_surf_vol;
233  CHKERR VecGetArray(surf_vol, &a_surf_vol);
234  cout << "Volume = " << a_vol[0] << endl;
235  cout << "Surf Volume = " << a_surf_vol[0] << endl;
236  if (fabs(a_vol[0] - a_surf_vol[0]) > 1e-12) {
237  SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "Should be zero");
238  }
239  CHKERR VecRestoreArray(vol, &a_vol);
240  CHKERR VecRestoreArray(vol, &a_surf_vol);
241  }
242 
243  }
244  }
245  CATCH_ERRORS;
246 
247  // finish work cleaning memory, getting statistics, etc.
249 
250  return 0;
251 }

Variable Documentation

◆ help

char help[] = "...\n\n"
static
Examples
simple_interface.cpp.

Definition at line 16 of file simple_interface.cpp.

MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
MoFEM::CoreInterface::loop_dofs
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
H1
@ H1
continuous field
Definition: definitions.h:85
MoFEM::CoreInterface::get_comm_rank
virtual int get_comm_rank() const =0
MoFEM::Simple::getSkeletonFEName
const std::string getSkeletonFEName() const
Get the Skeleton FE Name.
Definition: Simple.hpp:382
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
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
VolRule
Set integration rule to volume elements.
Definition: simple_interface.cpp:88
MoFEM::Projection10NodeCoordsOnField
Projection of edge entities with one mid-node on hierarchical basis.
Definition: Projection10NodeCoordsOnField.hpp:24
MoFEM::CoreTmp< 0 >::Finalize
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
MoFEM::Simple
Simple interface for fast problem set-up.
Definition: Simple.hpp:27
MoFEM::OpCalculateHOCoords
Calculate HO coordinates at gauss points.
Definition: HODataOperators.hpp:41
MoFEM::createDMMatrix
auto createDMMatrix(DM dm)
Get smart matrix from DM.
Definition: DMMoFEM.hpp:1056
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::Simple::addSkeletonField
MoFEMErrorCode addSkeletonField(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 skeleton.
Definition: Simple.cpp:372
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
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::createGhostVector
auto createGhostVector(MPI_Comm comm, PetscInt n, PetscInt N, PetscInt nghost, const PetscInt ghosts[])
Create smart ghost vector.
Definition: PetscSmartObj.hpp:179
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::DMRegister_MoFEM
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:43
OpFace
Definition: boundary_marker.cpp:19
MoFEM::Simple::getDomainFEName
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition: Simple.hpp:368
FaceRule
Set integration rule to boundary elements.
Definition: simple_interface.cpp:91
help
static char help[]
Definition: simple_interface.cpp:16
MoFEM::Simple::setFieldOrder
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:545
MoFEM::OpCalculateVectorFieldGradient
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Definition: UserDataOperators.hpp:1535
MoFEM::CoreTmp< 0 >::Initialize
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
MoFEM::vectorDuplicate
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
Definition: PetscSmartObj.hpp:221
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1148
HenckyOps::f
auto f
Definition: HenckyOps.hpp:15
MoFEM::OpSetHOWeights
Set inverse jacobian to base functions.
Definition: HODataOperators.hpp:159
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
MoFEM::Simple::getBoundaryFEName
const std::string getBoundaryFEName() const
Get the Boundary FE Name.
Definition: Simple.hpp:375
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
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
MOFEM_ATOM_TEST_INVALID
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40
MoFEM::Simple::setUp
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:683
MoFEM::DMoFEMLoopFiniteElements
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:586
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
OpVolume
Definition: simple_interface.cpp:18
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359