v0.14.0
Functions | Variables
minimal_surface_area.cpp File Reference
#include <MoFEM.hpp>
#include <MethodForForceScaling.hpp>
#include <DirichletBC.hpp>
#include <PostProcOnRefMesh.hpp>
#include <MinimalSurfaceElement.hpp>

Go to the source code of this file.

Functions

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

Variables

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

Function Documentation

◆ main()

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

Definition at line 41 of file minimal_surface_area.cpp.

41  {
42 
43  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
44 
45  try {
46 
47  moab::Core mb_instance;
48  moab::Interface &moab = mb_instance;
49 
50  char mesh_file_name[255];
51  PetscInt order = 2;
52  PetscInt nb_sub_steps = 10;
53  PetscBool flg_file = PETSC_FALSE;
54  PetscBool is_partitioned = PETSC_FALSE;
55  PetscBool is_atom_test = PETSC_FALSE;
56 
57  {
58 
59  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "",
60  "Minimal surface area config", "none");
61  CHKERR PetscOptionsString("-my_file", "mesh file name", "", "mesh.h5m",
62  mesh_file_name, 255, &flg_file);
63  CHKERR PetscOptionsInt("-my_order", "default approximation order", "", 2,
64  &order, PETSC_NULL);
65  CHKERR PetscOptionsInt("-my_nb_sub_steps", "number of substeps", "", 10,
66  &nb_sub_steps, PETSC_NULL);
67  CHKERR PetscOptionsBool("-my_is_partitioned",
68  "set if mesh is partitioned (this result that "
69  "each process keeps only part of the mesh",
70  "", PETSC_FALSE, &is_partitioned, PETSC_NULL);
71  CHKERR PetscOptionsBool(
72  "-my_is_atom_test",
73  "is used with testing, exit with error when diverged", "",
74  PETSC_FALSE, &is_atom_test, PETSC_NULL);
75 
76  ierr = PetscOptionsEnd();
77  CHKERRG(ierr);
78 
79  // Reade parameters from line command
80  if (flg_file != PETSC_TRUE) {
81  SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
82  }
83 
84  if (is_partitioned == PETSC_TRUE) {
85  // Read mesh to MOAB
86  const char *option;
87  option = "PARALLEL=BCAST_DELETE;"
88  "PARALLEL_RESOLVE_SHARED_ENTS;"
89  "PARTITION=PARALLEL_PARTITION;";
90  CHKERR moab.load_file(mesh_file_name, 0, option);
91  } else {
92  const char *option;
93  option = "";
94  CHKERR moab.load_file(mesh_file_name, 0, option);
95  }
96  ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
97  if (pcomm == NULL)
98  pcomm = new ParallelComm(&moab, PETSC_COMM_WORLD);
99  }
100 
101  MoFEM::Core core(moab);
102  MoFEM::Interface &m_field = core;
103  // meshset consisting all entities in mesh
104  EntityHandle root_set = moab.get_root_set();
105 
106  // Seed all mesh entities to MoFEM database, those entities can be
107  // potentially used as finite elements or as entities which carry some
108  // approximation field.
109  BitRefLevel bit_level0;
110  bit_level0.set(0);
111  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
112  0, 2, bit_level0);
113 
114  // Add field
115  CHKERR m_field.add_field("U", H1, AINSWORTH_LEGENDRE_BASE, 1);
116  // Add entities to field (root_mesh, i.e. on all mesh entities fields are
117  // approx.) On those entities and lower dimension entities approximation is
118  // spaned,
119  CHKERR m_field.add_ents_to_field_by_type(root_set, MBTRI, "U");
120  // Set approximation oder
121  CHKERR m_field.set_field_order(root_set, MBTRI, "U", order);
122  CHKERR m_field.set_field_order(root_set, MBEDGE, "U", order);
123  CHKERR m_field.set_field_order(root_set, MBVERTEX, "U", 1);
124  CHKERR m_field.build_fields();
125 
126  // Add finite element (this defines element declaration comes later)
127  CHKERR m_field.add_finite_element("MINIMAL_SURFACE_ELEMENT");
128  // Set on which fields that element operates
130  "MINIMAL_SURFACE_ELEMENT", "U");
132  "MINIMAL_SURFACE_ELEMENT", "U");
134  "MINIMAL_SURFACE_ELEMENT", "U");
135  // Add entities to that element
137  root_set, MBTRI, "MINIMAL_SURFACE_ELEMENT");
138 
139  CHKERR m_field.add_finite_element("BC_ELEMENT");
140  CHKERR m_field.modify_finite_element_add_field_row("BC_ELEMENT", "U");
141  CHKERR m_field.modify_finite_element_add_field_col("BC_ELEMENT", "U");
142  CHKERR m_field.modify_finite_element_add_field_data("BC_ELEMENT", "U");
143  // Add entities to that element
144  Range bc_edges, bc_ents;
145  {
146  // Take a skin form the mesh
147  Skinner skin(&m_field.get_moab());
148  Range tris;
149  CHKERR moab.get_entities_by_type(root_set, MBTRI, tris, false);
150  Range bc_edges;
151  CHKERR skin.find_skin(root_set, tris, false, bc_edges);
152  Range bc_nodes;
153  CHKERR moab.get_connectivity(bc_edges, bc_nodes);
154  bc_ents.merge(bc_edges);
155  bc_ents.merge(bc_nodes);
156  // Add entities to that element
157  CHKERR m_field.add_ents_to_finite_element_by_type(bc_edges, MBEDGE,
158  "BC_ELEMENT");
159  }
160 
161  // build finite elements
162  CHKERR m_field.build_finite_elements();
163  // build adjacencies between elements and degrees of freedom
164  CHKERR m_field.build_adjacencies(bit_level0);
165 
166  // Create minimal surface area and elements
167  MinimalSurfaceElement minimal_surface_element(m_field);
168  // Data structure for common data
169  MinimalSurfaceElement::CommonData mse_common_data;
170 
171  // This class is used to fix dofs on entities
172  DirichletFixFieldAtEntitiesBc fix_edges_ents(m_field, "U", NULL, NULL, NULL,
173  bc_ents);
174 
175  CHKERR DMRegister_MoFEM("MoFEM");
176  // Solve BC problem
177  {
178  // define problems
179  DM bc_dm;
180  // Register DM problem
181  CHKERR DMCreate(PETSC_COMM_WORLD, &bc_dm);
182  CHKERR DMSetType(bc_dm, "MoFEM");
183  // Create DM instance
184  CHKERR DMMoFEMCreateMoFEM(bc_dm, &m_field, "DM_BC", bit_level0);
185  // Configure DM form line command options (DM itself, solvers,
186  // pre-conditioners, ... )
187  CHKERR DMSetFromOptions(bc_dm);
188  // Add elements to dm (only one here)
189  CHKERR DMMoFEMAddElement(bc_dm, "BC_ELEMENT");
190  // Set up problem (number dofs, partition mesh, etc.)
191  CHKERR DMSetUp(bc_dm);
192  Mat A;
193  Vec T, F;
194  CHKERR DMCreateGlobalVector(bc_dm, &T);
195  CHKERR VecDuplicate(T, &F);
196  CHKERR DMCreateMatrix(bc_dm, &A);
197 
198  // Boundary problem
199  minimal_surface_element.feBcEdge.getOpPtrVector().push_back(
201  minimal_surface_element.feBcEdge.getOpPtrVector().push_back(
203  CHKERR DMoFEMLoopFiniteElements(bc_dm, "BC_ELEMENT",
204  &minimal_surface_element.feBcEdge);
205  CHKERR VecAssemblyBegin(F);
206  CHKERR VecAssemblyEnd(F);
207  CHKERR MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
208  CHKERR MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
209 
210  {
211  KSP solver;
212  CHKERR KSPCreate(PETSC_COMM_WORLD, &solver);
213  CHKERR KSPSetOperators(solver, A, A);
214  CHKERR KSPSetFromOptions(solver);
215  CHKERR KSPSetUp(solver);
216  CHKERR KSPSolve(solver, F, T);
217  CHKERR KSPDestroy(&solver);
218  }
219 
220  // Scatter data on mesh
221  CHKERR VecGhostUpdateBegin(T, INSERT_VALUES, SCATTER_FORWARD);
222  CHKERR VecGhostUpdateEnd(T, INSERT_VALUES, SCATTER_FORWARD);
223  CHKERR DMoFEMMeshToGlobalVector(bc_dm, T, INSERT_VALUES, SCATTER_REVERSE);
224 
225  CHKERR VecDestroy(&T);
226  CHKERR VecDestroy(&F);
227  CHKERR MatDestroy(&A);
228  CHKERR DMDestroy(&bc_dm);
229  }
230 
231  // //define problems
232  // CHKERR m_field.add_problem("DM_MINIMAL_AREA");
233  // //set refinement level for problem
234  // CHKERR
235  // m_field.modify_problem_ref_level_add_bit("DM_MINIMAL_AREA",bit_level0);
236  // Create discrete manager instance
237  DM dm;
238  {
239  // Register problem
240  CHKERR DMCreate(PETSC_COMM_WORLD, &dm);
241  CHKERR DMSetType(dm, "MoFEM");
242  // Create DM instance
243  CHKERR DMMoFEMCreateMoFEM(dm, &m_field, "DM_MINIMAL_AREA", bit_level0);
244  // Get options fron line command
245  CHKERR DMSetFromOptions(dm);
246  // Add elements to dm (only one here)
247  CHKERR DMMoFEMAddElement(dm, "MINIMAL_SURFACE_ELEMENT");
248  // Set up data structures
249  CHKERR DMSetUp(dm);
250  }
251 
252  // Create matrices and vectors used for analysis
253  Vec T, F;
254  Mat A;
255  {
256  CHKERR DMCreateGlobalVector(dm, &T);
257  CHKERR VecDuplicate(T, &F);
258  CHKERR DMCreateMatrix(dm, &A);
259  }
260 
261  // Adding elements to DMSnes
262  {
263 
264  // Right hand side operators (run for each entity in element)
265  // Calculate inverse of jacobian
266  auto det_ptr = boost::make_shared<VectorDouble>();
267  auto jac_ptr = boost::make_shared<MatrixDouble>();
268  auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
269  minimal_surface_element.feRhs.getOpPtrVector().push_back(
270  new OpCalculateHOJacForFace(jac_ptr));
271  minimal_surface_element.feRhs.getOpPtrVector().push_back(
272  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
273  // Apply inv_jac_ptr to shape functions
274  minimal_surface_element.feRhs.getOpPtrVector().push_back(
275  new OpSetInvJacH1ForFace(inv_jac_ptr));
276  // Below operators are specific to minimal area problem
277  minimal_surface_element.feRhs.getOpPtrVector().push_back(
278  new MinimalSurfaceElement::OpGetDataAtGaussPts("U", mse_common_data));
279  minimal_surface_element.feRhs.getOpPtrVector().push_back(
281  "U", mse_common_data, true));
282  minimal_surface_element.feRhs.getOpPtrVector().push_back(
283  new MinimalSurfaceElement::OpAssembleResidaul("U", mse_common_data));
284 
285  // Left hand side operators (run for each entity in element)
286  // Calculate inverse of jacobian
287  minimal_surface_element.feLhs.getOpPtrVector().push_back(
288  new OpCalculateHOJacForFace(jac_ptr));
289  minimal_surface_element.feLhs.getOpPtrVector().push_back(
290  new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
291  // Apply inv_jac_ptr to shape functions
292  minimal_surface_element.feLhs.getOpPtrVector().push_back(
293  new OpSetInvJacH1ForFace(inv_jac_ptr));
294  // Below operators are specific to minimal area problem
295  minimal_surface_element.feLhs.getOpPtrVector().push_back(
296  new MinimalSurfaceElement::OpGetDataAtGaussPts("U", mse_common_data));
297  minimal_surface_element.feLhs.getOpPtrVector().push_back(
299  "U", mse_common_data, false));
300  minimal_surface_element.feLhs.getOpPtrVector().push_back(
301  new MinimalSurfaceElement::OpAssembleTangent("U", mse_common_data));
302 
303  // Rhs
304  CHKERR DMMoFEMSNESSetFunction(dm, DM_NO_ELEMENT, NULL, &fix_edges_ents,
305  NULL);
306  CHKERR DMMoFEMSNESSetFunction(dm, "MINIMAL_SURFACE_ELEMENT",
307  &minimal_surface_element.feRhs, PETSC_NULL,
308  PETSC_NULL);
310  &fix_edges_ents);
311 
312  // Lhs
313  CHKERR DMMoFEMSNESSetJacobian(dm, DM_NO_ELEMENT, NULL, &fix_edges_ents,
314  NULL);
315  CHKERR DMMoFEMSNESSetJacobian(dm, "MINIMAL_SURFACE_ELEMENT",
316  &minimal_surface_element.feLhs, PETSC_NULL,
317  PETSC_NULL);
319  &fix_edges_ents);
320  }
321 
322  SNESConvergedReason snes_reason;
323  // Solve problem
324  {
325  SNES snes;
326  SnesCtx *snes_ctx;
327  CHKERR SNESCreate(PETSC_COMM_WORLD, &snes);
328  // CHKERR SNESSetDM(snes,dm);
329  CHKERR DMMoFEMGetSnesCtx(dm, &snes_ctx);
330  CHKERR SNESSetFunction(snes, F, SnesRhs, snes_ctx);
331  CHKERR SNESSetJacobian(snes, A, A, SnesMat, snes_ctx);
332  CHKERR SNESSetFromOptions(snes);
333 
334  Vec T0;
335  CHKERR VecDuplicate(T, &T0);
336  CHKERR DMoFEMMeshToLocalVector(dm, T0, INSERT_VALUES, SCATTER_FORWARD);
337 
338  double step_size = 1. / nb_sub_steps;
339  for (int ss = 1; ss <= nb_sub_steps; ss++) {
340  CHKERR VecAXPY(T, step_size, T0);
341  CHKERR DMoFEMMeshToLocalVector(dm, T, INSERT_VALUES, SCATTER_REVERSE);
342  CHKERR SNESSolve(snes, PETSC_NULL, T);
343  CHKERR SNESGetConvergedReason(snes, &snes_reason);
344  int its;
345  CHKERR SNESGetIterationNumber(snes, &its);
346  CHKERR PetscPrintf(PETSC_COMM_WORLD,
347  "number of Newton iterations = %D\n\n", its);
348  if (snes_reason < 0) {
349  if (is_atom_test) {
350  SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
351  "atom test diverged");
352  } else {
353  break;
354  }
355  }
356  }
357  CHKERR VecGhostUpdateBegin(T, INSERT_VALUES, SCATTER_FORWARD);
358  CHKERR VecGhostUpdateEnd(T, INSERT_VALUES, SCATTER_FORWARD);
359  CHKERR DMoFEMMeshToGlobalVector(dm, T, INSERT_VALUES, SCATTER_REVERSE);
360 
361  // CHKERR MatView(A,PETSC_VIEWER_DRAW_SELF);
362  // std::string wait;
363  // std::cin >> wait;
364 
365  CHKERR VecDestroy(&T0);
366  CHKERR SNESDestroy(&snes);
367  }
368 
369  {
370  PostProcFaceOnRefinedMesh post_proc(m_field);
371  CHKERR post_proc.generateReferenceElementMesh();
372  CHKERR post_proc.addFieldValuesPostProc("U");
373  CHKERR DMoFEMLoopFiniteElements(dm, "MINIMAL_SURFACE_ELEMENT",
374  &post_proc);
375  // Save data on mesh
376  CHKERR post_proc.postProcMesh.write_file("out.h5m", "MOAB",
377  "PARALLEL=WRITE_PART");
378  }
379 
380  // Clean and destroy
381  {
382  CHKERR DMDestroy(&dm);
383  CHKERR VecDestroy(&T);
384  CHKERR VecDestroy(&F);
385  CHKERR MatDestroy(&A);
386  }
387  }
388  CATCH_ERRORS;
389 
391  return 0;
392 }

Variable Documentation

◆ help

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

Definition at line 39 of file minimal_surface_area.cpp.

MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
help
static char help[]
Definition: minimal_surface_area.cpp:39
MYPCOMM_INDEX
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:228
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
H1
@ H1
continuous field
Definition: definitions.h:85
EntityHandle
MoFEM::CoreInterface::modify_finite_element_add_field_row
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string name_row)=0
set field row which finite element use
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
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
MinimalSurfaceEquation::MinimalSurfaceElement
Implementation of minimal area element.
Definition: MinimalSurfaceElement.hpp:37
MoFEM::CoreTmp< 0 >::Finalize
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
MoFEM::CoreInterface::add_ents_to_field_by_type
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
MoFEM::DMMoFEMAddElement
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition: DMMoFEM.cpp:497
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
DirichletFixFieldAtEntitiesBc
Fix dofs on entities.
Definition: DirichletBC.hpp:258
MinimalSurfaceEquation::MinimalSurfaceElement::OpAssmebleBcRhs
Integrate vector on rhs,.
Definition: MinimalSurfaceElement.hpp:91
MoFEM::CoreInterface::add_ents_to_finite_element_by_type
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::CoreInterface::add_finite_element
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
MoFEM::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
MoFEM::CoreInterface::modify_finite_element_add_field_col
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string name_row)=0
set field col which finite element use
MoFEM::CoreInterface::build_finite_elements
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elements.
DM_NO_ELEMENT
#define DM_NO_ELEMENT
Definition: DMMoFEM.hpp:10
MoFEM::DMoFEMMeshToGlobalVector
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMoFEM.cpp:535
MoFEM::DMMoFEMGetSnesCtx
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMoFEM.cpp:1094
MoFEM::DMRegister_MoFEM
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:43
MoFEM::DMMoFEMSNESSetJacobian
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Definition: DMMoFEM.cpp:759
MoFEM::DMMoFEMCreateMoFEM
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMoFEM.cpp:114
MoFEM::CoreInterface::modify_finite_element_add_field_data
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_field)=0
set finite element field data
PostProcFaceOnRefinedMesh
Postprocess on face.
Definition: PostProcOnRefMesh.hpp:1032
mesh_file_name
char mesh_file_name[255]
Definition: mesh_smoothing.cpp:22
MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleTangent
Assemble tangent matrix.
Definition: MinimalSurfaceElement.hpp:328
MoFEM::OpSetInvJacH1ForFace
Definition: UserDataOperators.hpp:2918
MoFEM::SnesRhs
PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx)
This is MoFEM implementation for the right hand side (residual vector) evaluation in SNES solver.
Definition: SnesCtx.cpp:27
MinimalSurfaceEquation::MinimalSurfaceElement::OpCalculateCoefficientsAtGaussPts
Evaluate function values and gradients at Gauss Pts.
Definition: MinimalSurfaceElement.hpp:228
Range
MinimalSurfaceEquation::MinimalSurfaceElement::OpAssmebleBcLhs
Integrate vector on lhs,.
Definition: MinimalSurfaceElement.hpp:139
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
MinimalSurfaceEquation::MinimalSurfaceElement::CommonData
Keep date between operators.
Definition: MinimalSurfaceElement.hpp:75
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1148
MinimalSurfaceEquation::MinimalSurfaceElement::OpGetDataAtGaussPts
Evaluate function values and gradients at Gauss Pts.
Definition: MinimalSurfaceElement.hpp:190
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
EigenMatrix::Vec
const FTensor::Tensor2< T, Dim, Dim > Vec
Definition: MatrixFunction.hpp:66
MoFEM::OpInvertMatrix
Definition: UserDataOperators.hpp:3249
MoFEM::CoreInterface::build_fields
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
MoFEM::OpCalculateHOJacForFace
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
Definition: HODataOperators.hpp:264
MoFEM::BitRefManager
Managing BitRefLevels.
Definition: BitRefManager.hpp:21
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
MOFEM_ATOM_TEST_INVALID
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40
MoFEM::CoreInterface::build_adjacencies
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
MoFEM::SnesCtx
Interface for nonlinear (SNES) solver.
Definition: SnesCtx.hpp:13
MoFEM::CoreInterface::set_field_order
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)=0
Set order approximation of the entities in the field.
MoFEM::SnesMat
PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx)
This is MoFEM implementation for the left hand side (tangent matrix) evaluation in SNES solver.
Definition: SnesCtx.cpp:139
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
MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleResidaul
Assemble residual.
Definition: MinimalSurfaceElement.hpp:295
CHKERRG
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
MoFEM::CoreInterface::add_field
virtual MoFEMErrorCode add_field(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_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add field.
F
@ F
Definition: free_surface.cpp:394
MoFEM::DMMoFEMSNESSetFunction
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
Definition: DMMoFEM.cpp:718