v0.13.1
mesh_smoothing.cpp
Go to the documentation of this file.
1/** \file mesh_smoothing.cpp
2 * \brief Improve mesh quality using Volume-length quality measure with barrier
3 * \example mesh_smoothing.cpp
4 *
5 */
6
7/* This file is part of MoFEM.
8 * MoFEM is free software: you can redistribute it and/or modify it under
9 * the terms of the GNU Lesser General Public License as published by the
10 * Free Software Foundation, either version 3 of the License, or (at your
11 * option) any later version.
12 *
13 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
14 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
15 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16 * License for more details.
17 *
18 * You should have received a copy of the GNU Lesser General Public
19 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
20
21#include <MoFEM.hpp>
22using namespace MoFEM;
23
24#include <BasicFiniteElements.hpp>
25
26#include <Smoother.hpp>
29
30using namespace MoFEM;
31
32static char help[] = "mesh smoothing\n\n";
33
34PetscBool flg_myfile = PETSC_TRUE;
38PetscBool output_vtk = PETSC_TRUE;
39double tol = 0.1;
40double gamma_factor = 0.8;
41
42int main(int argc, char *argv[]) {
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",
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", "",
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();
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 {
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");
119 ->synchroniseFieldEntities("LAMBDA_EDGE");
120 CHKERR m_field.set_field_order(0, MBVERTEX, "LAMBDA_EDGE", 1);
121
122 CHKERR m_field.add_finite_element("EDGE_SLIDING");
124 "EDGE_SLIDING", "MESH_NODE_POSITIONS");
126 "EDGE_SLIDING", "MESH_NODE_POSITIONS");
128 "EDGE_SLIDING", "MESH_NODE_POSITIONS");
129 CHKERR m_field.modify_finite_element_add_field_row("EDGE_SLIDING",
130 "LAMBDA_EDGE");
131 CHKERR m_field.modify_finite_element_add_field_col("EDGE_SLIDING",
132 "LAMBDA_EDGE");
133 CHKERR m_field.modify_finite_element_add_field_data("EDGE_SLIDING",
134 "LAMBDA_EDGE");
135 CHKERR m_field.add_ents_to_finite_element_by_type(edges, MBEDGE,
136 "EDGE_SLIDING");
137 simple_interface->getOtherFiniteElements().push_back("EDGE_SLIDING");
138 }
139 }
140
141 CHKERR simple_interface->defineFiniteElements();
142 CHKERR simple_interface->defineProblem();
143 CHKERR simple_interface->buildFields();
144 // Remove vertices form LAMBDA_SURFACE which are on the edges
147 Range edges;
148 CHKERR m_field.getInterface<MeshsetsManager>()->getEntitiesByDimension(
149 edges_block_set, BLOCKSET, 1, edges, true);
150 Range verts;
151 CHKERR m_field.get_moab().get_connectivity(edges, verts, true);
152 CHKERR m_field.remove_ents_from_field("LAMBDA_SURFACE", verts);
153 }
154 CHKERR simple_interface->buildFiniteElements();
155 CHKERR simple_interface->buildProblem();
156 }
157
158 struct ElementsAndOperators {
159
160 MoFEM::Interface &mField;
161 Vec minQualityVec;
162 double *minQualityPtr;
163
164 ElementsAndOperators(MoFEM::Interface &m_field) : mField(m_field) {
165 ierr = VecCreateMPI(PETSC_COMM_WORLD, 1, m_field.get_comm_size(),
166 &minQualityVec);
167 CHKERRABORT(PETSC_COMM_WORLD, ierr);
168 ierr = VecGetArray(minQualityVec, &minQualityPtr);
169 CHKERRABORT(PETSC_COMM_WORLD, ierr);
170 }
171
172 virtual ~ElementsAndOperators() {
173 ierr = VecRestoreArray(minQualityVec, &minQualityPtr);
174 CHKERRABORT(PETSC_COMM_WORLD, ierr);
175 ierr = VecDestroy(&minQualityVec);
176 CHKERRABORT(PETSC_COMM_WORLD, ierr);
177 }
178
179 double getMinQuality() const { return *minQualityPtr; }
180
181 enum Tags {
182 SMOOTHING_TAG = 1,
183 SURFACE_CONSTRAINS_TAG,
184 EDGE_CONSTRAINS_TAG
185 };
186
187 boost::shared_ptr<Smoother> smootherFe;
188 boost::shared_ptr<Smoother::MyVolumeFE>
189 feSmootherRhs; ///< Integrate smoothing operators
190 boost::shared_ptr<Smoother::MyVolumeFE>
191 feSmootherLhs; ///< Integrate smoothing operators
192 boost::shared_ptr<VolumeLengthQuality<double> > volumeLengthDouble;
193 boost::shared_ptr<VolumeLengthQuality<adouble> > volumeLengthAdouble;
194
195 boost::shared_ptr<SurfaceSlidingConstrains> surfaceConstrain;
196 boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
197 skinOrientation;
198 boost::shared_ptr<EdgeSlidingConstrains> edgeConstrain;
199
200 boost::shared_ptr<DirichletFixFieldAtEntitiesBc> fixMaterialEnts;
201
202 boost::shared_ptr<MoFEM::VolumeElementForcesAndSourcesCore> minQualityFe;
203
205 double *minQualityPtr;
206 MinQualityOp(double *min_quality_ptr)
208 "MESH_NODE_POSITIONS", UserDataOperator::OPROW),
209 minQualityPtr(min_quality_ptr) {}
210 MoFEMErrorCode doWork(int side, EntityType type,
213 if (type != MBVERTEX)
215 double q = Tools::volumeLengthQuality(&*data.getFieldData().begin());
216 *minQualityPtr = fmin(*minQualityPtr, q);
218 }
219 };
220
221 MoFEMErrorCode createSmoothingFE() {
223 smootherFe = boost::shared_ptr<Smoother>(new Smoother(mField));
224 volumeLengthAdouble = boost::shared_ptr<VolumeLengthQuality<adouble> >(
226 volumeLengthDouble = boost::shared_ptr<VolumeLengthQuality<double> >(
228
229 Range tets;
230 CHKERR mField.get_moab().get_entities_by_type(0, MBTET, tets);
231 smootherFe->setOfBlocks[0].tEts.merge(tets);
232
233 smootherFe->setOfBlocks[0].materialDoublePtr = volumeLengthDouble;
234 smootherFe->setOfBlocks[0].materialAdoublePtr = volumeLengthAdouble;
235
236 // set element data
237 smootherFe->commonData.spatialPositions = "MESH_NODE_POSITIONS";
238 smootherFe->commonData.meshPositions = "NONE";
239
240 smootherFe->feRhs.meshPositionsFieldName = "NONE";
241 smootherFe->feLhs.meshPositionsFieldName = "NONE";
242 smootherFe->feRhs.addToRule = 0;
243 smootherFe->feLhs.addToRule = 0;
244
245 feSmootherRhs = smootherFe->feRhsPtr;
246 feSmootherLhs = smootherFe->feLhsPtr;
247
248 // Smoother right hand side
249 feSmootherRhs->getOpPtrVector().push_back(
251 "MESH_NODE_POSITIONS", smootherFe->commonData));
252 feSmootherRhs->getOpPtrVector().push_back(
254 "MESH_NODE_POSITIONS", smootherFe->setOfBlocks.at(0),
255 smootherFe->commonData, SMOOTHING_TAG, false));
256 feSmootherRhs->getOpPtrVector().push_back(new Smoother::OpRhsSmoother(
257 "MESH_NODE_POSITIONS", smootherFe->setOfBlocks[0],
258 smootherFe->commonData, smootherFe->smootherData));
259
260 // Smoother left hand side
261 feSmootherLhs->getOpPtrVector().push_back(
263 "MESH_NODE_POSITIONS", smootherFe->commonData));
264 feSmootherLhs->getOpPtrVector().push_back(
266 "MESH_NODE_POSITIONS", smootherFe->setOfBlocks.at(0),
267 smootherFe->commonData, SMOOTHING_TAG, true));
268 feSmootherLhs->getOpPtrVector().push_back(new Smoother::OpLhsSmoother(
269 "MESH_NODE_POSITIONS", "MESH_NODE_POSITIONS",
270 smootherFe->setOfBlocks.at(0), smootherFe->commonData,
271 smootherFe->smootherData, "LAMBDA_CRACKFRONT_AREA_TANGENT"));
272
273 minQualityFe =
274 boost::shared_ptr<MoFEM::VolumeElementForcesAndSourcesCore>(
276 minQualityFe->getOpPtrVector().push_back(
277 new MinQualityOp(minQualityPtr));
278
279 Range fixed_vertex;
282 CHKERR mField.getInterface<MeshsetsManager>()->getEntitiesByDimension(
283 vertex_block_set, BLOCKSET, 0, fixed_vertex, true);
284 }
285 fixMaterialEnts = boost::shared_ptr<DirichletFixFieldAtEntitiesBc>(
286 new DirichletFixFieldAtEntitiesBc(mField, "MESH_NODE_POSITIONS",
287 fixed_vertex));
288 fixMaterialEnts->fieldNames.push_back("LAMBDA_SURFACE");
289 fixMaterialEnts->fieldNames.push_back("LAMBDA_EDGE");
290
292 }
293
294 MoFEMErrorCode createConstrians() {
296 skinOrientation = boost::shared_ptr<
299 surfaceConstrain = boost::shared_ptr<SurfaceSlidingConstrains>(
300 skinOrientation,
301 new SurfaceSlidingConstrains(mField, *skinOrientation));
302 surfaceConstrain->setOperators(SURFACE_CONSTRAINS_TAG, "LAMBDA_SURFACE",
303 "MESH_NODE_POSITIONS");
304
307 Range edges;
309 ->getEntitiesByDimension(edges_block_set, BLOCKSET, 1, edges,
310 true);
311
312 Range tets;
313 CHKERR mField.get_moab().get_entities_by_type(0, MBTET, tets);
314 Skinner skin(&mField.get_moab());
315 Range skin_faces; // skin faces from 3d ents
316 CHKERR skin.find_skin(0, tets, false, skin_faces);
317
318 edgeConstrain = boost::shared_ptr<EdgeSlidingConstrains>(
319 new EdgeSlidingConstrains(mField));
320 CHKERR edgeConstrain->setOperators(EDGE_CONSTRAINS_TAG, edges,
321 skin_faces, "LAMBDA_EDGE",
322 "MESH_NODE_POSITIONS");
323
324 // CHKERR EdgeSlidingConstrains::CalculateEdgeBase::saveEdges(
325 // mField.get_moab(), "out_edges.vtk", edges);
326
327 }
329 }
330
331 MoFEMErrorCode addFEtoDM(DM dm) {
333 boost::shared_ptr<ForcesAndSourcesCore> null;
334
335 CHKERR DMMoFEMSNESSetFunction(dm, DM_NO_ELEMENT, null, fixMaterialEnts,
336 null);
337 CHKERR DMMoFEMSNESSetFunction(dm, "SMOOTHING", feSmootherRhs, null,
338 null);
339 CHKERR DMMoFEMSNESSetFunction(dm, "SURFACE_SLIDING",
340 surfaceConstrain->feRhsPtr, null, null);
341 CHKERR DMMoFEMSNESSetFunction(dm, "EDGE_SLIDING",
342 edgeConstrain->feRhsPtr, null, null);
344 fixMaterialEnts);
345
346 CHKERR DMMoFEMSNESSetJacobian(dm, DM_NO_ELEMENT, null, fixMaterialEnts,
347 null);
348 CHKERR DMMoFEMSNESSetJacobian(dm, "SMOOTHING", feSmootherLhs, null,
349 null);
350 CHKERR DMMoFEMSNESSetJacobian(dm, "SURFACE_SLIDING",
351 surfaceConstrain->feLhsPtr, null, null);
352 CHKERR DMMoFEMSNESSetJacobian(dm, "EDGE_SLIDING",
353 edgeConstrain->feLhsPtr, null, null);
355 fixMaterialEnts);
356
357 // MoFEM::SnesCtx *snes_ctx;
358 // DMMoFEMGetSnesCtx(dm,&snes_ctx);
359 // snes_ctx->vErify = true;
360
362 }
363
364 MoFEMErrorCode calcuteMinQuality(DM dm) {
366 *minQualityPtr = 1;
367 CHKERR DMoFEMLoopFiniteElements(dm, "SMOOTHING", minQualityFe.get());
368 CHKERR VecMin(minQualityVec, PETSC_NULL, minQualityPtr);
370 }
371 };
372
373 ElementsAndOperators elements_and_operators(m_field);
374 CHKERR elements_and_operators.createSmoothingFE();
375 CHKERR elements_and_operators.createConstrians();
376
377 DM dm;
378 CHKERR simple_interface->getDM(&dm);
379 CHKERR elements_and_operators.addFEtoDM(dm);
380
381 struct Solve {
382
383 MoFEMErrorCode operator()(DM dm) const {
385
386 // Create the right hand side vector and vector of unknowns
387 Vec F, D;
388 CHKERR DMCreateGlobalVector(dm, &F);
389 // Create unknown vector by creating duplicate copy of F vector. only
390 // structure is duplicated no values.
391 CHKERR VecDuplicate(F, &D);
392
393 CHKERR zeroLambdaFields(dm);
394 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_FORWARD);
395
396 // Create solver and link it to DM
397 SNES solver;
398 CHKERR SNESCreate(PETSC_COMM_WORLD, &solver);
399 CHKERR SNESSetFromOptions(solver);
400 CHKERR SNESSetDM(solver, dm);
401 // Set-up solver, is type of solver and pre-conditioners
402 CHKERR SNESSetUp(solver);
403 // At solution process, KSP solver using DM creates matrices, Calculate
404 // values of the left hand side and the right hand side vector. then
405 // solves system of equations. Results are stored in vector D.
406 CHKERR SNESSolve(solver, F, D);
407
408 // Scatter solution on the mesh. Stores unknown vector on field on the
409 // mesh.
410 CHKERR DMoFEMMeshToGlobalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
411 // Clean data. Solver and vector are not needed any more.
412 CHKERR SNESDestroy(&solver);
413 CHKERR VecDestroy(&D);
414 CHKERR VecDestroy(&F);
415
417 }
418
419 MoFEMErrorCode setCoordsFromField(DM dm) const {
421 MoFEM::Interface *m_field_ptr;
422 CHKERR DMoFEMGetInterfacePtr(dm, &m_field_ptr);
423 for (_IT_GET_ENT_FIELD_BY_NAME_FOR_LOOP_(*m_field_ptr,
424 "MESH_NODE_POSITIONS", it)) {
425 if (it->get()->getEntType() != MBVERTEX)
426 continue;
427 VectorDouble3 coords(3);
428 for(int dd = 0;dd!=3;++dd)
429 coords[dd] = it->get()->getEntFieldData()[dd];
430 EntityHandle ent = it->get()->getEnt();
431 CHKERR m_field_ptr->get_moab().set_coords(&ent, 1, &*coords.begin());
432 }
434 }
435
436 MoFEMErrorCode setFieldFromCoords(DM dm) const {
438 MoFEM::Interface *m_field_ptr;
439 CHKERR DMoFEMGetInterfacePtr(dm, &m_field_ptr);
440 for (_IT_GET_ENT_FIELD_BY_NAME_FOR_LOOP_(*m_field_ptr,
441 "MESH_NODE_POSITIONS", it)) {
442 if (it->get()->getEntType() != MBVERTEX)
443 continue;
444 EntityHandle ent = it->get()->getEnt();
445 VectorDouble3 coords(3);
446 CHKERR m_field_ptr->get_moab().get_coords(&ent, 1, &*coords.begin());
447 for(int dd = 0;dd!=3;++dd)
448 it->get()->getEntFieldData()[dd] = coords[dd];
449 }
451 }
452
453 private:
454 MoFEMErrorCode zeroLambdaFields(DM dm) const {
456 MoFEM::Interface *m_field_ptr;
457 CHKERR DMoFEMGetInterfacePtr(dm, &m_field_ptr);
458 CHKERR m_field_ptr->getInterface<FieldBlas>()->setField(
459 0, MBVERTEX, "LAMBDA_SURFACE");
461 }
462
463 };
464
465 Solve solve;
466 CHKERR solve.setFieldFromCoords(dm);
467
468 CHKERR elements_and_operators.calcuteMinQuality(dm);
469 double min_quality = elements_and_operators.getMinQuality();
470 PetscPrintf(PETSC_COMM_WORLD, "Min quality = %4.3f\n", min_quality);
471
472 double gamma = min_quality > 0 ? gamma_factor * min_quality
473 : min_quality / gamma_factor;
474 elements_and_operators.volumeLengthDouble->gAmma = gamma;
475 elements_and_operators.volumeLengthAdouble->gAmma = gamma;
476
477 double min_quality_p, eps;
478 do {
479
480 min_quality_p = min_quality;
481
482 CHKERR solve(dm);
483
484 CHKERR solve.setCoordsFromField(dm);
485 CHKERR elements_and_operators.calcuteMinQuality(dm);
486 min_quality = elements_and_operators.getMinQuality();
487
488 eps = (min_quality - min_quality_p) / min_quality;
489 PetscPrintf(PETSC_COMM_WORLD, "Min quality = %4.3f eps = %4.3f\n",
490 min_quality, eps);
491
492 double gamma = min_quality > 0 ? gamma_factor * min_quality
493 : min_quality / gamma_factor;
494 elements_and_operators.volumeLengthDouble->gAmma = gamma;
495 elements_and_operators.volumeLengthAdouble->gAmma = gamma;
496
497 } while (eps > tol);
498
499 // if (m_field.getInterface<MeshsetsManager>()->checkMeshset(edges_block_set,
500 // BLOCKSET)) {
501 // Range edges;
502 // CHKERR m_field.getInterface<MeshsetsManager>()->getEntitiesByDimension(
503 // edges_block_set, BLOCKSET, 1, edges, true);
504
505 // Range tets;
506 // CHKERR moab.get_entities_by_type(0,MBTET,tets);
507 // Skinner skin(&moab);
508 // Range skin_faces; // skin faces from 3d ents
509 // CHKERR skin.find_skin(0, tets, false, skin_faces);
510
511 // CHKERR EdgeSlidingConstrains::CalculateEdgeBase::setTags(moab, edges,
512 // skin_faces);
513 // CHKERR EdgeSlidingConstrains::CalculateEdgeBase::saveEdges(
514 // moab, "out_edges.vtk", edges);
515 // }
516
517 if (output_vtk)
518 CHKERR m_field.getInterface<BitRefManager>()->writeBitLevelByType(
519 BitRefLevel().set(0), BitRefLevel().set(), MBTET, "out.vtk", "VTK",
520 "");
521 }
523
525}
#define DM_NO_ELEMENT
Definition: DMMoFEM.hpp:22
Implementing surface sliding constrains.
Implementation of Volume-Length-Quality measure with barrier.
@ BARRIER_AND_QUALITY
static const double eps
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
@ MF_ZERO
Definition: definitions.h:111
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:75
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
@ H1
continuous field
Definition: definitions.h:98
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
@ BLOCKSET
Definition: definitions.h:161
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
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:677
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:482
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:718
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:545
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:373
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:494
#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:1797
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
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
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
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
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
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.
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.
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
bool checkMeshset(const int ms_id, const CubitBCType cubit_bc_type) const
check for CUBIT Id and CUBIT type
int main(int argc, char *argv[])
double tol
PetscBool output_vtk
static char help[]
PetscBool flg_myfile
double gamma_factor
int vertex_block_set
int edges_block_set
char mesh_file_name[255]
const FTensor::Tensor2< T, Dim, Dim > Vec
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 MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:103
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
CoreTmp< 0 > Core
Definition: Core.hpp:1096
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
const double D
diffusivity
Fix dofs on entities.
Managing BitRefLevels.
Managing BitRefLevels.
virtual int get_comm_size() const =0
virtual moab::Interface & get_moab()=0
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.
Core (interface) class.
Definition: Core.hpp:92
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:85
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:125
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
const VectorDouble & getFieldData() const
get dofs values
Basic algebra on fields.
Definition: FieldBlas.hpp:31
Interface for managing meshsets containing materials and boundary conditions.
Simple interface for fast problem set-up.
Definition: Simple.hpp:35
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:739
const std::string getBoundaryFEName() const
Get the Boundary FE Name.
Definition: Simple.hpp:328
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:497
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:713
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:378
std::vector< std::string > & getOtherFiniteElements()
Get the Other Finite Elements.
Definition: Simple.hpp:380
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:294
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition: Simple.cpp:810
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:601
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:308
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:593
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:396
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:569
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition: Simple.hpp:321
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Class implemented by user to detect face orientation.
Shape preserving constrains, i.e. nodes sliding on body surface.
Volume Length Quality.