v0.8.14
Functions | Variables
mesh_smoothing.cpp File Reference

Improve mesh quality using Volume-length quality measure with barrier. More...

#include <MoFEM.hpp>
#include <BasicFiniteElements.hpp>
#include <Smoother.hpp>
#include <SurfaceSlidingConstrains.hpp>
#include <VolumeLengthQuality.hpp>

Go to the source code of this file.

Functions

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

Variables

static char help [] = "mesh smoothing\n\n"
 
PetscBool flg_myfile = PETSC_TRUE
 
char mesh_file_name [255]
 
int edges_block_set = 1
 
int vertex_block_set = 2
 
PetscBool output_vtk = PETSC_TRUE
 
double tol = 0.1
 
double gamma_factor = 0.8
 

Detailed Description

Improve mesh quality using Volume-length quality measure with barrier.

Definition in file mesh_smoothing.cpp.

Function Documentation

◆ main()

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

Definition at line 42 of file mesh_smoothing.cpp.

42  {
43 
44  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
45 
46  try {
47 
48  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "Mesh cut options", "none");
49  CHKERR PetscOptionsString("-my_file", "mesh file name", "", "mesh.h5m",
50  mesh_file_name, 255, &flg_myfile);
51  CHKERR PetscOptionsInt("-edges_block_set", "edges side set", "",
52  edges_block_set, &edges_block_set, PETSC_NULL);
53  CHKERR PetscOptionsInt("-vertex_block_set", "vertex side set", "",
54  vertex_block_set, &vertex_block_set, PETSC_NULL);
55  CHKERR PetscOptionsBool("-output_vtk", "if true outout vtk file", "",
56  output_vtk, &output_vtk, PETSC_NULL);
57  CHKERR PetscOptionsScalar("-quality_reduction_tol", "",
58  "Tolerance of quality reduction", tol, &tol,
59  PETSC_NULL);
60  CHKERR PetscOptionsScalar("-gamma_factor", "",
61  "Gamma factor", gamma_factor, &gamma_factor,
62  PETSC_NULL);
63  ierr = PetscOptionsEnd();
64  CHKERRG(ierr);
65 
66  // Create MoAB database
67  moab::Core moab_core; // create database
68  moab::Interface &moab = moab_core; // create interface to database
69  // Create MoFEM database and link it to MoAB
70  MoFEM::Core mofem_core(moab); // create database
71  MoFEM::Interface &m_field = mofem_core; // create interface to database
72  // Register DM Manager
73  CHKERR DMRegister_MoFEM("DMMOFEM"); // register MoFEM DM in PETSc
74 
75  // Get simple interface is simplified version enabling quick and
76  // easy construction of problem.
77  Simple *simple_interface;
78  // Query interface and get pointer to Simple interface
79  CHKERR m_field.getInterface(simple_interface);
80 
81  // Build problem with simple interface
82  {
83  // Get options for simple interface from command line
84  CHKERR simple_interface->getOptions();
85  // Load mesh file to database
86  CHKERR simple_interface->loadFile();
87 
88  // Add domain filed "U" in space H1 and Legendre base, Ainsworth recipe is
89  // used to construct base functions.
90  CHKERR simple_interface->addDomainField("MESH_NODE_POSITIONS", H1,
92  // Add Lagrange multiplier field on body boundary
93  CHKERR simple_interface->addBoundaryField("LAMBDA_SURFACE", H1,
95 
96  // Set fields order domain and boundary fields.
97  CHKERR simple_interface->setFieldOrder("MESH_NODE_POSITIONS",
98  1); // to approximate function
99  CHKERR simple_interface->setFieldOrder("LAMBDA_SURFACE",
100  1); // to Lagrange multipliers
101 
102  simple_interface->getDomainFEName() = "SMOOTHING";
103  simple_interface->getBoundaryFEName() = "SURFACE_SLIDING";
104 
105  // Other fields and finite elements added to database directly
106  {
107  if (m_field.getInterface<MeshsetsManager>()->checkMeshset(
109  // Declare approximation fields
110  CHKERR m_field.add_field("LAMBDA_EDGE", H1, AINSWORTH_LOBATTO_BASE, 2,
111  MB_TAG_SPARSE, MF_ZERO);
112  Range edges;
114  ->getEntitiesByDimension(edges_block_set, BLOCKSET, 1, edges,
115  true);
116  CHKERR m_field.add_ents_to_field_by_type(edges, MBEDGE,
117  "LAMBDA_EDGE");
118  CHKERR m_field.synchronise_field_entities("LAMBDA_EDGE");
119  CHKERR m_field.set_field_order(0, MBVERTEX, "LAMBDA_EDGE", 1);
120 
121  CHKERR m_field.add_finite_element("EDGE_SLIDING");
123  "EDGE_SLIDING", "MESH_NODE_POSITIONS");
125  "EDGE_SLIDING", "MESH_NODE_POSITIONS");
127  "EDGE_SLIDING", "MESH_NODE_POSITIONS");
128  CHKERR m_field.modify_finite_element_add_field_row("EDGE_SLIDING",
129  "LAMBDA_EDGE");
130  CHKERR m_field.modify_finite_element_add_field_col("EDGE_SLIDING",
131  "LAMBDA_EDGE");
132  CHKERR m_field.modify_finite_element_add_field_data("EDGE_SLIDING",
133  "LAMBDA_EDGE");
134  CHKERR m_field.add_ents_to_finite_element_by_type(edges, MBEDGE,
135  "EDGE_SLIDING");
136  simple_interface->getOtherFiniteElements().push_back("EDGE_SLIDING");
137  }
138  }
139 
140  CHKERR simple_interface->defineFiniteElements();
141  CHKERR simple_interface->defineProblem();
142  CHKERR simple_interface->buildFields();
143  // Remove vertices form LAMBDA_SURFACE which are on the edges
144  if (m_field.getInterface<MeshsetsManager>()->checkMeshset(
146  Range edges;
147  CHKERR m_field.getInterface<MeshsetsManager>()->getEntitiesByDimension(
148  edges_block_set, BLOCKSET, 1, edges, true);
149  Range verts;
150  CHKERR m_field.get_moab().get_connectivity(edges, verts, true);
151  CHKERR m_field.remove_ents_from_field("LAMBDA_SURFACE", verts);
152  }
153  CHKERR simple_interface->buildFiniteElements();
154  CHKERR simple_interface->buildProblem();
155  }
156 
157  struct ElementsAndOperators {
158 
159  MoFEM::Interface &mField;
160  Vec minQualityVec;
161  double *minQualityPtr;
162 
163  ElementsAndOperators(MoFEM::Interface &m_field) : mField(m_field) {
164  ierr = VecCreateMPI(PETSC_COMM_WORLD, 1, m_field.get_comm_size(),
165  &minQualityVec);
166  CHKERRABORT(PETSC_COMM_WORLD, ierr);
167  ierr = VecGetArray(minQualityVec, &minQualityPtr);
168  CHKERRABORT(PETSC_COMM_WORLD, ierr);
169  }
170 
171  virtual ~ElementsAndOperators() {
172  ierr = VecRestoreArray(minQualityVec, &minQualityPtr);
173  CHKERRABORT(PETSC_COMM_WORLD, ierr);
174  ierr = VecDestroy(&minQualityVec);
175  CHKERRABORT(PETSC_COMM_WORLD, ierr);
176  }
177 
178  double getMinQuality() const { return *minQualityPtr; }
179 
180  enum Tags {
181  SMOOTHING_TAG = 1,
182  SURFACE_CONSTRAINS_TAG,
183  EDGE_CONSTRAINS_TAG
184  };
185 
186  boost::shared_ptr<Smoother> smootherFe;
187  boost::shared_ptr<Smoother::MyVolumeFE>
188  feSmootherRhs; ///< Integrate smoothing operators
189  boost::shared_ptr<Smoother::MyVolumeFE>
190  feSmootherLhs; ///< Integrate smoothing operators
191  boost::shared_ptr<VolumeLengthQuality<double> > volumeLengthDouble;
192  boost::shared_ptr<VolumeLengthQuality<adouble> > volumeLengthAdouble;
193 
194  boost::shared_ptr<SurfaceSlidingConstrains> surfaceConstrain;
195  boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
196  skinOrientation;
197  boost::shared_ptr<EdgeSlidingConstrains> edgeConstrain;
198 
199  boost::shared_ptr<DirichletFixFieldAtEntitiesBc> fixMaterialEnts;
200 
201  boost::shared_ptr<MoFEM::VolumeElementForcesAndSourcesCore> minQualityFe;
202 
203  struct MinQualityOp : MoFEM::ForcesAndSourcesCore::UserDataOperator {
204  double *minQualityPtr;
205  MinQualityOp(double *min_quality_ptr)
207  "MESH_NODE_POSITIONS", UserDataOperator::OPROW),
208  minQualityPtr(min_quality_ptr) {}
209  MoFEMErrorCode doWork(int side, EntityType type,
212  if (type != MBVERTEX)
214  double q = Tools::volumeLengthQuality(&*data.getFieldData().begin());
215  *minQualityPtr = fmin(*minQualityPtr, q);
217  }
218  };
219 
220  MoFEMErrorCode createSmoothingFE() {
222  smootherFe = boost::shared_ptr<Smoother>(new Smoother(mField));
223  volumeLengthAdouble = boost::shared_ptr<VolumeLengthQuality<adouble> >(
225  volumeLengthDouble = boost::shared_ptr<VolumeLengthQuality<double> >(
227 
228  Range tets;
229  CHKERR mField.get_moab().get_entities_by_type(0, MBTET, tets);
230  smootherFe->setOfBlocks[0].tEts.merge(tets);
231 
232  smootherFe->setOfBlocks[0].materialDoublePtr = volumeLengthDouble;
233  smootherFe->setOfBlocks[0].materialAdoublePtr = volumeLengthAdouble;
234 
235  // set element data
236  smootherFe->commonData.spatialPositions = "MESH_NODE_POSITIONS";
237  smootherFe->commonData.meshPositions = "NONE";
238 
239  smootherFe->feRhs.meshPositionsFieldName = "NONE";
240  smootherFe->feLhs.meshPositionsFieldName = "NONE";
241  smootherFe->feRhs.addToRule = 0;
242  smootherFe->feLhs.addToRule = 0;
243 
244  feSmootherRhs = smootherFe->feRhsPtr;
245  feSmootherLhs = smootherFe->feLhsPtr;
246 
247  // Smoother right hand side
248  feSmootherRhs->getOpPtrVector().push_back(
250  "MESH_NODE_POSITIONS", smootherFe->commonData));
251  feSmootherRhs->getOpPtrVector().push_back(
253  "MESH_NODE_POSITIONS", smootherFe->setOfBlocks.at(0),
254  smootherFe->commonData, SMOOTHING_TAG, false));
255  feSmootherRhs->getOpPtrVector().push_back(new Smoother::OpRhsSmoother(
256  "MESH_NODE_POSITIONS", smootherFe->setOfBlocks[0],
257  smootherFe->commonData, smootherFe->smootherData));
258 
259  // Smoother left hand side
260  feSmootherLhs->getOpPtrVector().push_back(
262  "MESH_NODE_POSITIONS", smootherFe->commonData));
263  feSmootherLhs->getOpPtrVector().push_back(
265  "MESH_NODE_POSITIONS", smootherFe->setOfBlocks.at(0),
266  smootherFe->commonData, SMOOTHING_TAG, true));
267  feSmootherLhs->getOpPtrVector().push_back(new Smoother::OpLhsSmoother(
268  "MESH_NODE_POSITIONS", "MESH_NODE_POSITIONS",
269  smootherFe->setOfBlocks.at(0), smootherFe->commonData,
270  smootherFe->smootherData, "LAMBDA_CRACKFRONT_AREA_TANGENT"));
271 
272  minQualityFe =
273  boost::shared_ptr<MoFEM::VolumeElementForcesAndSourcesCore>(
275  minQualityFe->getOpPtrVector().push_back(
276  new MinQualityOp(minQualityPtr));
277 
278  Range fixed_vertex;
281  CHKERR mField.getInterface<MeshsetsManager>()->getEntitiesByDimension(
282  vertex_block_set, BLOCKSET, 0, fixed_vertex, true);
283  }
284  fixMaterialEnts = boost::shared_ptr<DirichletFixFieldAtEntitiesBc>(
285  new DirichletFixFieldAtEntitiesBc(mField, "MESH_NODE_POSITIONS",
286  fixed_vertex));
287  fixMaterialEnts->fieldNames.push_back("LAMBDA_SURFACE");
288  fixMaterialEnts->fieldNames.push_back("LAMBDA_EDGE");
289 
291  }
292 
293  MoFEMErrorCode createConstrians() {
295  skinOrientation = boost::shared_ptr<
296  SurfaceSlidingConstrains::DriverElementOrientation>(
297  new SurfaceSlidingConstrains::DriverElementOrientation());
298  surfaceConstrain = boost::shared_ptr<SurfaceSlidingConstrains>(
299  skinOrientation,
300  new SurfaceSlidingConstrains(mField, *skinOrientation));
301  surfaceConstrain->setOperators(SURFACE_CONSTRAINS_TAG, "LAMBDA_SURFACE",
302  "MESH_NODE_POSITIONS");
303 
306  Range edges;
308  ->getEntitiesByDimension(edges_block_set, BLOCKSET, 1, edges,
309  true);
310 
311  Range tets;
312  CHKERR mField.get_moab().get_entities_by_type(0, MBTET, tets);
313  Skinner skin(&mField.get_moab());
314  Range skin_faces; // skin faces from 3d ents
315  CHKERR skin.find_skin(0, tets, false, skin_faces);
316 
317  edgeConstrain = boost::shared_ptr<EdgeSlidingConstrains>(
318  new EdgeSlidingConstrains(mField));
319  CHKERR edgeConstrain->setOperators(EDGE_CONSTRAINS_TAG, edges,
320  skin_faces, "LAMBDA_EDGE",
321  "MESH_NODE_POSITIONS");
322 
323  // CHKERR EdgeSlidingConstrains::CalculateEdgeBase::saveEdges(
324  // mField.get_moab(), "out_edges.vtk", edges);
325 
326  }
328  }
329 
330  MoFEMErrorCode addFEtoDM(DM dm) {
332  boost::shared_ptr<ForcesAndSourcesCore> null;
333 
334  CHKERR DMMoFEMSNESSetFunction(dm, DM_NO_ELEMENT, null, fixMaterialEnts,
335  null);
336  CHKERR DMMoFEMSNESSetFunction(dm, "SMOOTHING", feSmootherRhs, null,
337  null);
338  CHKERR DMMoFEMSNESSetFunction(dm, "SURFACE_SLIDING",
339  surfaceConstrain->feRhsPtr, null, null);
340  CHKERR DMMoFEMSNESSetFunction(dm, "EDGE_SLIDING",
341  edgeConstrain->feRhsPtr, null, null);
343  fixMaterialEnts);
344 
345  CHKERR DMMoFEMSNESSetJacobian(dm, DM_NO_ELEMENT, null, fixMaterialEnts,
346  null);
347  CHKERR DMMoFEMSNESSetJacobian(dm, "SMOOTHING", feSmootherLhs, null,
348  null);
349  CHKERR DMMoFEMSNESSetJacobian(dm, "SURFACE_SLIDING",
350  surfaceConstrain->feLhsPtr, null, null);
351  CHKERR DMMoFEMSNESSetJacobian(dm, "EDGE_SLIDING",
352  edgeConstrain->feLhsPtr, null, null);
354  fixMaterialEnts);
355 
356  // MoFEM::SnesCtx *snes_ctx;
357  // DMMoFEMGetSnesCtx(dm,&snes_ctx);
358  // snes_ctx->vErify = true;
359 
361  }
362 
363  MoFEMErrorCode calcuteMinQuality(DM dm) {
365  *minQualityPtr = 1;
366  CHKERR DMoFEMLoopFiniteElements(dm, "SMOOTHING", minQualityFe.get());
367  CHKERR VecMin(minQualityVec, PETSC_NULL, minQualityPtr);
369  }
370  };
371 
372  ElementsAndOperators elements_and_operators(m_field);
373  CHKERR elements_and_operators.createSmoothingFE();
374  CHKERR elements_and_operators.createConstrians();
375 
376  DM dm;
377  CHKERR simple_interface->getDM(&dm);
378  CHKERR elements_and_operators.addFEtoDM(dm);
379 
380  struct Solve {
381 
382  MoFEMErrorCode operator()(DM dm) const {
384 
385  // Create the right hand side vector and vector of unknowns
386  Vec F, D;
387  CHKERR DMCreateGlobalVector(dm, &F);
388  // Create unknown vector by creating duplicate copy of F vector. only
389  // structure is duplicated no values.
390  CHKERR VecDuplicate(F, &D);
391 
392  CHKERR zeroLambdaFields(dm);
393  CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_FORWARD);
394 
395  // Create solver and link it to DM
396  SNES solver;
397  CHKERR SNESCreate(PETSC_COMM_WORLD, &solver);
398  CHKERR SNESSetFromOptions(solver);
399  CHKERR SNESSetDM(solver, dm);
400  // Set-up solver, is type of solver and pre-conditioners
401  CHKERR SNESSetUp(solver);
402  // At solution process, KSP solver using DM creates matrices, Calculate
403  // values of the left hand side and the right hand side vector. then
404  // solves system of equations. Results are stored in vector D.
405  CHKERR SNESSolve(solver, F, D);
406 
407  // Scatter solution on the mesh. Stores unknown vector on field on the
408  // mesh.
409  CHKERR DMoFEMMeshToGlobalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
410  // Clean data. Solver and vector are not needed any more.
411  CHKERR SNESDestroy(&solver);
412  CHKERR VecDestroy(&D);
413  CHKERR VecDestroy(&F);
414 
416  }
417 
418  MoFEMErrorCode setCoordsFromField(DM dm) const {
420  MoFEM::Interface *m_field_ptr;
421  CHKERR DMoFEMGetInterfacePtr(dm, &m_field_ptr);
422  for (_IT_GET_ENT_FIELD_BY_NAME_FOR_LOOP_(*m_field_ptr,
423  "MESH_NODE_POSITIONS", it)) {
424  if (it->get()->getEntType() != MBVERTEX)
425  continue;
426  VectorDouble3 coords(3);
427  for(int dd = 0;dd!=3;++dd)
428  coords[dd] = it->get()->getEntFieldData()[dd];
429  EntityHandle ent = it->get()->getEnt();
430  CHKERR m_field_ptr->get_moab().set_coords(&ent, 1, &*coords.begin());
431  }
433  }
434 
435  MoFEMErrorCode setFieldFromCoords(DM dm) const {
437  MoFEM::Interface *m_field_ptr;
438  CHKERR DMoFEMGetInterfacePtr(dm, &m_field_ptr);
439  for (_IT_GET_ENT_FIELD_BY_NAME_FOR_LOOP_(*m_field_ptr,
440  "MESH_NODE_POSITIONS", it)) {
441  if (it->get()->getEntType() != MBVERTEX)
442  continue;
443  EntityHandle ent = it->get()->getEnt();
444  VectorDouble3 coords(3);
445  CHKERR m_field_ptr->get_moab().get_coords(&ent, 1, &*coords.begin());
446  for(int dd = 0;dd!=3;++dd)
447  it->get()->getEntFieldData()[dd] = coords[dd];
448  }
450  }
451 
452  private:
453  MoFEMErrorCode zeroLambdaFields(DM dm) const {
455  MoFEM::Interface *m_field_ptr;
456  CHKERR DMoFEMGetInterfacePtr(dm, &m_field_ptr);
457  CHKERR m_field_ptr->getInterface<FieldBlas>()->setField(
458  0, MBVERTEX, "LAMBDA_SURFACE");
460  }
461 
462  };
463 
464  Solve solve;
465  CHKERR solve.setFieldFromCoords(dm);
466 
467  CHKERR elements_and_operators.calcuteMinQuality(dm);
468  double min_quality = elements_and_operators.getMinQuality();
469  PetscPrintf(PETSC_COMM_WORLD, "Min quality = %4.3f\n", min_quality);
470 
471  double gamma = min_quality > 0 ? gamma_factor * min_quality
472  : min_quality / gamma_factor;
473  elements_and_operators.volumeLengthDouble->gAmma = gamma;
474  elements_and_operators.volumeLengthAdouble->gAmma = gamma;
475 
476  double min_quality_p, eps;
477  do {
478 
479  min_quality_p = min_quality;
480 
481  CHKERR solve(dm);
482 
483  CHKERR solve.setCoordsFromField(dm);
484  CHKERR elements_and_operators.calcuteMinQuality(dm);
485  min_quality = elements_and_operators.getMinQuality();
486 
487  eps = (min_quality - min_quality_p) / min_quality;
488  PetscPrintf(PETSC_COMM_WORLD, "Min quality = %4.3f eps = %4.3f\n",
489  min_quality, eps);
490 
491  double gamma = min_quality > 0 ? gamma_factor * min_quality
492  : min_quality / gamma_factor;
493  elements_and_operators.volumeLengthDouble->gAmma = gamma;
494  elements_and_operators.volumeLengthAdouble->gAmma = gamma;
495 
496  } while (eps > tol);
497 
498  // if (m_field.getInterface<MeshsetsManager>()->checkMeshset(edges_block_set,
499  // BLOCKSET)) {
500  // Range edges;
501  // CHKERR m_field.getInterface<MeshsetsManager>()->getEntitiesByDimension(
502  // edges_block_set, BLOCKSET, 1, edges, true);
503 
504  // Range tets;
505  // CHKERR moab.get_entities_by_type(0,MBTET,tets);
506  // Skinner skin(&moab);
507  // Range skin_faces; // skin faces from 3d ents
508  // CHKERR skin.find_skin(0, tets, false, skin_faces);
509 
510  // CHKERR EdgeSlidingConstrains::CalculateEdgeBase::setTags(moab, edges,
511  // skin_faces);
512  // CHKERR EdgeSlidingConstrains::CalculateEdgeBase::saveEdges(
513  // moab, "out_edges.vtk", edges);
514  // }
515 
516  if (output_vtk)
517  CHKERR m_field.getInterface<BitRefManager>()->writeBitLevelByType(
518  BitRefLevel().set(0), BitRefLevel().set(), MBTET, "out.vtk", "VTK",
519  "");
520  }
521  CATCH_ERRORS;
522 
524 }
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
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:553
#define _IT_GET_ENT_FIELD_BY_NAME_FOR_LOOP_(MFIELD, NAME, IT)
loop over all dofs from a moFEM field and particular field
Definition: Interface.hpp:1680
#define DM_NO_ELEMENT
Definition: DMMoFEM.hpp:22
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: DMMMoFEM.cpp:644
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
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition: Simple.cpp:605
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:186
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Common.hpp:60
virtual moab::Interface & get_moab()=0
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
static char help[]
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:452
virtual int get_comm_size() const =0
PetscBool output_vtk
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:526
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.
Interface for managing meshsets containing materials and boundary conditions.
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.
Core (interface) class.
Definition: Core.hpp:50
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
Volume finite elementUser is implementing own operator at Gauss point level, by class derived from Vo...
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:514
const std::string & getDomainFEName() const
Definition: Simple.hpp:197
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
char mesh_file_name[255]
Simple interface for fast problem set-up.
Definition: Simple.hpp:33
Volume Length Quality.
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:51
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Common.hpp:147
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
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:139
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:579
double tol
Basic algebra on fields.
Definition: FieldBlas.hpp:34
int edges_block_set
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:320
double gamma_factor
const VectorDouble & getFieldData() const
get dofs values
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.
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:464
bool checkMeshset(const int ms_id, const CubitBCType cubit_bc_type) const
check for CUBIT Id and CUBIT type
Managing BitRefLevels.
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:91
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:142
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:123
virtual MoFEMErrorCode synchronise_field_entities(const BitFieldId id, int verb=DEFAULT_VERBOSITY)=0
virtual MoFEMErrorCode remove_ents_from_field(const std::string &name, const EntityHandle meshset, const EntityType type, int verb=DEFAULT_VERBOSITY)=0
remove entities from field
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: DMMMoFEM.cpp:685
MoFEMErrorCode loadFile()
Load mesh file.
Definition: Simple.cpp:75
int vertex_block_set
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define CHKERR
Inline error check.
Definition: definitions.h:578
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
ForcesAndSourcesCore::UserDataOperator UserDataOperator
Fix dofs on entities.
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:347
Data on single entity (This is passed as argument to DataOperator::doWork)
structure to get information form mofem into DataForcesAndSourcesCore
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:60
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
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:429
std::vector< std::string > & getOtherFiniteElements()
Definition: Simple.hpp:205
continuous field
Definition: definitions.h:168
static const double eps
Data operator to do calculations at integration points.Is inherited and implemented by user to do cal...
const std::string & getBoundaryFEName() const
Definition: Simple.hpp:198
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:312
PetscBool flg_myfile
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:285
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Common.hpp:220
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:61

Variable Documentation

◆ edges_block_set

int edges_block_set = 1
Examples:
mesh_cut.cpp, and mesh_smoothing.cpp.

Definition at line 36 of file mesh_smoothing.cpp.

◆ flg_myfile

PetscBool flg_myfile = PETSC_TRUE
Examples:
mesh_cut.cpp, and mesh_smoothing.cpp.

Definition at line 34 of file mesh_smoothing.cpp.

◆ gamma_factor

double gamma_factor = 0.8
Examples:
mesh_smoothing.cpp.

Definition at line 40 of file mesh_smoothing.cpp.

◆ help

char help[] = "mesh smoothing\n\n"
static
Examples:
mesh_smoothing.cpp.

Definition at line 32 of file mesh_smoothing.cpp.

◆ mesh_file_name

char mesh_file_name[255]

◆ output_vtk

PetscBool output_vtk = PETSC_TRUE
Examples:
mesh_cut.cpp, mesh_smoothing.cpp, and split_sideset.cpp.

Definition at line 38 of file mesh_smoothing.cpp.

◆ tol

double tol = 0.1
Examples:
MagneticElement.hpp, and mesh_smoothing.cpp.

Definition at line 39 of file mesh_smoothing.cpp.

◆ vertex_block_set

int vertex_block_set = 2
Examples:
mesh_cut.cpp, and mesh_smoothing.cpp.

Definition at line 37 of file mesh_smoothing.cpp.