v0.10.0
mesh_smoothing.cpp
/** \file mesh_smoothing.cpp
* \brief Improve mesh quality using Volume-length quality measure with barrier
* \example mesh_smoothing.cpp
*
*/
/* This file is part of MoFEM.
* MoFEM is free software: you can redistribute it and/or modify it under
* the terms of the GNU Lesser General Public License as published by the
* Free Software Foundation, either version 3 of the License, or (at your
* option) any later version.
*
* MoFEM is distributed in the hope that it will be useful, but WITHOUT
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
* License for more details.
*
* You should have received a copy of the GNU Lesser General Public
* License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
#include <MoFEM.hpp>
using namespace MoFEM;
#include <BasicFiniteElements.hpp>
#include <Smoother.hpp>
using namespace MoFEM;
static char help[] = "mesh smoothing\n\n";
PetscBool flg_myfile = PETSC_TRUE;
char mesh_file_name[255];
PetscBool output_vtk = PETSC_TRUE;
double tol = 0.1;
double gamma_factor = 0.8;
int main(int argc, char *argv[]) {
MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
try {
CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "Mesh cut options", "none");
CHKERR PetscOptionsString("-my_file", "mesh file name", "", "mesh.h5m",
CHKERR PetscOptionsInt("-edges_block_set", "edges side set", "",
CHKERR PetscOptionsInt("-vertex_block_set", "vertex side set", "",
CHKERR PetscOptionsBool("-output_vtk", "if true outout vtk file", "",
output_vtk, &output_vtk, PETSC_NULL);
CHKERR PetscOptionsScalar("-quality_reduction_tol", "",
"Tolerance of quality reduction", tol, &tol,
PETSC_NULL);
CHKERR PetscOptionsScalar("-gamma_factor", "",
"Gamma factor", gamma_factor, &gamma_factor,
PETSC_NULL);
ierr = PetscOptionsEnd();
// Create MoAB database
moab::Core moab_core; // create database
moab::Interface &moab = moab_core; // create interface to database
// Create MoFEM database and link it to MoAB
MoFEM::Core mofem_core(moab); // create database
MoFEM::Interface &m_field = mofem_core; // create interface to database
// Register DM Manager
CHKERR DMRegister_MoFEM("DMMOFEM"); // register MoFEM DM in PETSc
// Get simple interface is simplified version enabling quick and
// easy construction of problem.
Simple *simple_interface;
// Query interface and get pointer to Simple interface
CHKERR m_field.getInterface(simple_interface);
// Build problem with simple interface
{
// Get options for simple interface from command line
CHKERR simple_interface->getOptions();
// Load mesh file to database
CHKERR simple_interface->loadFile();
// Add domain filed "U" in space H1 and Legendre base, Ainsworth recipe is
// used to construct base functions.
CHKERR simple_interface->addDomainField("MESH_NODE_POSITIONS", H1,
// Add Lagrange multiplier field on body boundary
CHKERR simple_interface->addBoundaryField("LAMBDA_SURFACE", H1,
// Set fields order domain and boundary fields.
CHKERR simple_interface->setFieldOrder("MESH_NODE_POSITIONS",
1); // to approximate function
CHKERR simple_interface->setFieldOrder("LAMBDA_SURFACE",
1); // to Lagrange multipliers
simple_interface->getDomainFEName() = "SMOOTHING";
simple_interface->getBoundaryFEName() = "SURFACE_SLIDING";
// Other fields and finite elements added to database directly
{
if (m_field.getInterface<MeshsetsManager>()->checkMeshset(
// Declare approximation fields
CHKERR m_field.add_field("LAMBDA_EDGE", H1, AINSWORTH_LOBATTO_BASE, 2,
MB_TAG_SPARSE, MF_ZERO);
Range edges;
CHKERR m_field.getInterface<MeshsetsManager>()
->getEntitiesByDimension(edges_block_set, BLOCKSET, 1, edges,
true);
CHKERR m_field.add_ents_to_field_by_type(edges, MBEDGE,
"LAMBDA_EDGE");
CHKERR m_field.getInterface<CommInterface>()
->synchroniseFieldEntities("LAMBDA_EDGE");
CHKERR m_field.set_field_order(0, MBVERTEX, "LAMBDA_EDGE", 1);
CHKERR m_field.add_finite_element("EDGE_SLIDING");
"EDGE_SLIDING", "MESH_NODE_POSITIONS");
"EDGE_SLIDING", "MESH_NODE_POSITIONS");
"EDGE_SLIDING", "MESH_NODE_POSITIONS");
"LAMBDA_EDGE");
"LAMBDA_EDGE");
"LAMBDA_EDGE");
"EDGE_SLIDING");
simple_interface->getOtherFiniteElements().push_back("EDGE_SLIDING");
}
}
CHKERR simple_interface->defineFiniteElements();
CHKERR simple_interface->defineProblem();
CHKERR simple_interface->buildFields();
// Remove vertices form LAMBDA_SURFACE which are on the edges
if (m_field.getInterface<MeshsetsManager>()->checkMeshset(
Range edges;
CHKERR m_field.getInterface<MeshsetsManager>()->getEntitiesByDimension(
edges_block_set, BLOCKSET, 1, edges, true);
Range verts;
CHKERR m_field.get_moab().get_connectivity(edges, verts, true);
CHKERR m_field.remove_ents_from_field("LAMBDA_SURFACE", verts);
}
CHKERR simple_interface->buildFiniteElements();
CHKERR simple_interface->buildProblem();
}
struct ElementsAndOperators {
Vec minQualityVec;
double *minQualityPtr;
ElementsAndOperators(MoFEM::Interface &m_field) : mField(m_field) {
ierr = VecCreateMPI(PETSC_COMM_WORLD, 1, m_field.get_comm_size(),
&minQualityVec);
CHKERRABORT(PETSC_COMM_WORLD, ierr);
ierr = VecGetArray(minQualityVec, &minQualityPtr);
CHKERRABORT(PETSC_COMM_WORLD, ierr);
}
virtual ~ElementsAndOperators() {
ierr = VecRestoreArray(minQualityVec, &minQualityPtr);
CHKERRABORT(PETSC_COMM_WORLD, ierr);
ierr = VecDestroy(&minQualityVec);
CHKERRABORT(PETSC_COMM_WORLD, ierr);
}
double getMinQuality() const { return *minQualityPtr; }
enum Tags {
SURFACE_CONSTRAINS_TAG,
EDGE_CONSTRAINS_TAG
};
boost::shared_ptr<Smoother> smootherFe;
boost::shared_ptr<Smoother::MyVolumeFE>
feSmootherRhs; ///< Integrate smoothing operators
boost::shared_ptr<Smoother::MyVolumeFE>
feSmootherLhs; ///< Integrate smoothing operators
boost::shared_ptr<VolumeLengthQuality<double> > volumeLengthDouble;
boost::shared_ptr<VolumeLengthQuality<adouble> > volumeLengthAdouble;
boost::shared_ptr<SurfaceSlidingConstrains> surfaceConstrain;
boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
skinOrientation;
boost::shared_ptr<EdgeSlidingConstrains> edgeConstrain;
boost::shared_ptr<DirichletFixFieldAtEntitiesBc> fixMaterialEnts;
boost::shared_ptr<MoFEM::VolumeElementForcesAndSourcesCore> minQualityFe;
double *minQualityPtr;
MinQualityOp(double *min_quality_ptr)
"MESH_NODE_POSITIONS", UserDataOperator::OPROW),
minQualityPtr(min_quality_ptr) {}
MoFEMErrorCode doWork(int side, EntityType type,
if (type != MBVERTEX)
double q = Tools::volumeLengthQuality(&*data.getFieldData().begin());
*minQualityPtr = fmin(*minQualityPtr, q);
}
};
MoFEMErrorCode createSmoothingFE() {
smootherFe = boost::shared_ptr<Smoother>(new Smoother(mField));
volumeLengthAdouble = boost::shared_ptr<VolumeLengthQuality<adouble> >(
volumeLengthDouble = boost::shared_ptr<VolumeLengthQuality<double> >(
Range tets;
CHKERR mField.get_moab().get_entities_by_type(0, MBTET, tets);
smootherFe->setOfBlocks[0].tEts.merge(tets);
smootherFe->setOfBlocks[0].materialDoublePtr = volumeLengthDouble;
smootherFe->setOfBlocks[0].materialAdoublePtr = volumeLengthAdouble;
// set element data
smootherFe->commonData.spatialPositions = "MESH_NODE_POSITIONS";
smootherFe->commonData.meshPositions = "NONE";
smootherFe->feRhs.meshPositionsFieldName = "NONE";
smootherFe->feLhs.meshPositionsFieldName = "NONE";
smootherFe->feRhs.addToRule = 0;
smootherFe->feLhs.addToRule = 0;
feSmootherRhs = smootherFe->feRhsPtr;
feSmootherLhs = smootherFe->feLhsPtr;
// Smoother right hand side
feSmootherRhs->getOpPtrVector().push_back(
"MESH_NODE_POSITIONS", smootherFe->commonData));
feSmootherRhs->getOpPtrVector().push_back(
"MESH_NODE_POSITIONS", smootherFe->setOfBlocks.at(0),
smootherFe->commonData, SMOOTHING_TAG, false));
feSmootherRhs->getOpPtrVector().push_back(new Smoother::OpRhsSmoother(
"MESH_NODE_POSITIONS", smootherFe->setOfBlocks[0],
smootherFe->commonData, smootherFe->smootherData));
// Smoother left hand side
feSmootherLhs->getOpPtrVector().push_back(
"MESH_NODE_POSITIONS", smootherFe->commonData));
feSmootherLhs->getOpPtrVector().push_back(
"MESH_NODE_POSITIONS", smootherFe->setOfBlocks.at(0),
smootherFe->commonData, SMOOTHING_TAG, true));
feSmootherLhs->getOpPtrVector().push_back(new Smoother::OpLhsSmoother(
"MESH_NODE_POSITIONS", "MESH_NODE_POSITIONS",
smootherFe->setOfBlocks.at(0), smootherFe->commonData,
smootherFe->smootherData, "LAMBDA_CRACKFRONT_AREA_TANGENT"));
minQualityFe =
boost::shared_ptr<MoFEM::VolumeElementForcesAndSourcesCore>(
minQualityFe->getOpPtrVector().push_back(
new MinQualityOp(minQualityPtr));
Range fixed_vertex;
if (mField.getInterface<MeshsetsManager>()->checkMeshset(
CHKERR mField.getInterface<MeshsetsManager>()->getEntitiesByDimension(
vertex_block_set, BLOCKSET, 0, fixed_vertex, true);
}
fixMaterialEnts = boost::shared_ptr<DirichletFixFieldAtEntitiesBc>(
new DirichletFixFieldAtEntitiesBc(mField, "MESH_NODE_POSITIONS",
fixed_vertex));
fixMaterialEnts->fieldNames.push_back("LAMBDA_SURFACE");
fixMaterialEnts->fieldNames.push_back("LAMBDA_EDGE");
}
MoFEMErrorCode createConstrians() {
skinOrientation = boost::shared_ptr<
surfaceConstrain = boost::shared_ptr<SurfaceSlidingConstrains>(
skinOrientation,
new SurfaceSlidingConstrains(mField, *skinOrientation));
surfaceConstrain->setOperators(SURFACE_CONSTRAINS_TAG, "LAMBDA_SURFACE",
"MESH_NODE_POSITIONS");
if (mField.getInterface<MeshsetsManager>()->checkMeshset(
Range edges;
CHKERR mField.getInterface<MeshsetsManager>()
->getEntitiesByDimension(edges_block_set, BLOCKSET, 1, edges,
true);
Range tets;
CHKERR mField.get_moab().get_entities_by_type(0, MBTET, tets);
Skinner skin(&mField.get_moab());
Range skin_faces; // skin faces from 3d ents
CHKERR skin.find_skin(0, tets, false, skin_faces);
edgeConstrain = boost::shared_ptr<EdgeSlidingConstrains>(
new EdgeSlidingConstrains(mField));
CHKERR edgeConstrain->setOperators(EDGE_CONSTRAINS_TAG, edges,
skin_faces, "LAMBDA_EDGE",
"MESH_NODE_POSITIONS");
// CHKERR EdgeSlidingConstrains::CalculateEdgeBase::saveEdges(
// mField.get_moab(), "out_edges.vtk", edges);
}
}
MoFEMErrorCode addFEtoDM(DM dm) {
boost::shared_ptr<ForcesAndSourcesCore> null;
CHKERR DMMoFEMSNESSetFunction(dm, DM_NO_ELEMENT, null, fixMaterialEnts,
null);
CHKERR DMMoFEMSNESSetFunction(dm, "SMOOTHING", feSmootherRhs, null,
null);
CHKERR DMMoFEMSNESSetFunction(dm, "SURFACE_SLIDING",
surfaceConstrain->feRhsPtr, null, null);
CHKERR DMMoFEMSNESSetFunction(dm, "EDGE_SLIDING",
edgeConstrain->feRhsPtr, null, null);
fixMaterialEnts);
CHKERR DMMoFEMSNESSetJacobian(dm, DM_NO_ELEMENT, null, fixMaterialEnts,
null);
CHKERR DMMoFEMSNESSetJacobian(dm, "SMOOTHING", feSmootherLhs, null,
null);
CHKERR DMMoFEMSNESSetJacobian(dm, "SURFACE_SLIDING",
surfaceConstrain->feLhsPtr, null, null);
CHKERR DMMoFEMSNESSetJacobian(dm, "EDGE_SLIDING",
edgeConstrain->feLhsPtr, null, null);
fixMaterialEnts);
// MoFEM::SnesCtx *snes_ctx;
// DMMoFEMGetSnesCtx(dm,&snes_ctx);
// snes_ctx->vErify = true;
}
MoFEMErrorCode calcuteMinQuality(DM dm) {
*minQualityPtr = 1;
CHKERR DMoFEMLoopFiniteElements(dm, "SMOOTHING", minQualityFe.get());
CHKERR VecMin(minQualityVec, PETSC_NULL, minQualityPtr);
}
};
ElementsAndOperators elements_and_operators(m_field);
CHKERR elements_and_operators.createSmoothingFE();
CHKERR elements_and_operators.createConstrians();
DM dm;
CHKERR simple_interface->getDM(&dm);
CHKERR elements_and_operators.addFEtoDM(dm);
struct Solve {
MoFEMErrorCode operator()(DM dm) const {
// Create the right hand side vector and vector of unknowns
Vec F, D;
CHKERR DMCreateGlobalVector(dm, &F);
// Create unknown vector by creating duplicate copy of F vector. only
// structure is duplicated no values.
CHKERR VecDuplicate(F, &D);
CHKERR zeroLambdaFields(dm);
CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_FORWARD);
// Create solver and link it to DM
SNES solver;
CHKERR SNESCreate(PETSC_COMM_WORLD, &solver);
CHKERR SNESSetFromOptions(solver);
CHKERR SNESSetDM(solver, dm);
// Set-up solver, is type of solver and pre-conditioners
CHKERR SNESSetUp(solver);
// At solution process, KSP solver using DM creates matrices, Calculate
// values of the left hand side and the right hand side vector. then
// solves system of equations. Results are stored in vector D.
CHKERR SNESSolve(solver, F, D);
// Scatter solution on the mesh. Stores unknown vector on field on the
// mesh.
CHKERR DMoFEMMeshToGlobalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
// Clean data. Solver and vector are not needed any more.
CHKERR SNESDestroy(&solver);
CHKERR VecDestroy(&D);
CHKERR VecDestroy(&F);
}
MoFEMErrorCode setCoordsFromField(DM dm) const {
MoFEM::Interface *m_field_ptr;
CHKERR DMoFEMGetInterfacePtr(dm, &m_field_ptr);
"MESH_NODE_POSITIONS", it)) {
if (it->get()->getEntType() != MBVERTEX)
continue;
VectorDouble3 coords(3);
for(int dd = 0;dd!=3;++dd)
coords[dd] = it->get()->getEntFieldData()[dd];
EntityHandle ent = it->get()->getEnt();
CHKERR m_field_ptr->get_moab().set_coords(&ent, 1, &*coords.begin());
}
}
MoFEMErrorCode setFieldFromCoords(DM dm) const {
MoFEM::Interface *m_field_ptr;
CHKERR DMoFEMGetInterfacePtr(dm, &m_field_ptr);
"MESH_NODE_POSITIONS", it)) {
if (it->get()->getEntType() != MBVERTEX)
continue;
EntityHandle ent = it->get()->getEnt();
VectorDouble3 coords(3);
CHKERR m_field_ptr->get_moab().get_coords(&ent, 1, &*coords.begin());
for(int dd = 0;dd!=3;++dd)
it->get()->getEntFieldData()[dd] = coords[dd];
}
}
private:
MoFEMErrorCode zeroLambdaFields(DM dm) const {
MoFEM::Interface *m_field_ptr;
CHKERR DMoFEMGetInterfacePtr(dm, &m_field_ptr);
CHKERR m_field_ptr->getInterface<FieldBlas>()->setField(
0, MBVERTEX, "LAMBDA_SURFACE");
}
};
Solve solve;
CHKERR solve.setFieldFromCoords(dm);
CHKERR elements_and_operators.calcuteMinQuality(dm);
double min_quality = elements_and_operators.getMinQuality();
PetscPrintf(PETSC_COMM_WORLD, "Min quality = %4.3f\n", min_quality);
double gamma = min_quality > 0 ? gamma_factor * min_quality
: min_quality / gamma_factor;
elements_and_operators.volumeLengthDouble->gAmma = gamma;
elements_and_operators.volumeLengthAdouble->gAmma = gamma;
double min_quality_p, eps;
do {
min_quality_p = min_quality;
CHKERR solve(dm);
CHKERR solve.setCoordsFromField(dm);
CHKERR elements_and_operators.calcuteMinQuality(dm);
min_quality = elements_and_operators.getMinQuality();
eps = (min_quality - min_quality_p) / min_quality;
PetscPrintf(PETSC_COMM_WORLD, "Min quality = %4.3f eps = %4.3f\n",
min_quality, eps);
double gamma = min_quality > 0 ? gamma_factor * min_quality
: min_quality / gamma_factor;
elements_and_operators.volumeLengthDouble->gAmma = gamma;
elements_and_operators.volumeLengthAdouble->gAmma = gamma;
} while (eps > tol);
// if (m_field.getInterface<MeshsetsManager>()->checkMeshset(edges_block_set,
// BLOCKSET)) {
// Range edges;
// CHKERR m_field.getInterface<MeshsetsManager>()->getEntitiesByDimension(
// edges_block_set, BLOCKSET, 1, edges, true);
// Range tets;
// CHKERR moab.get_entities_by_type(0,MBTET,tets);
// Skinner skin(&moab);
// Range skin_faces; // skin faces from 3d ents
// CHKERR skin.find_skin(0, tets, false, skin_faces);
// CHKERR EdgeSlidingConstrains::CalculateEdgeBase::setTags(moab, edges,
// skin_faces);
// CHKERR EdgeSlidingConstrains::CalculateEdgeBase::saveEdges(
// moab, "out_edges.vtk", edges);
// }
CHKERR m_field.getInterface<BitRefManager>()->writeBitLevelByType(
BitRefLevel().set(0), BitRefLevel().set(), MBTET, "out.vtk", "VTK",
"");
}
}
#define DM_NO_ELEMENT
Definition: DMMoFEM.hpp:22
Implementing surface sliding constrains.
Implementation of Volume-Length-Quality measure with barrier.
@ BARRIER_AND_QUALITY
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:441
@ MF_ZERO
Definition: definitions.h:189
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:154
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
@ H1
continuous field
Definition: definitions.h:177
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:552
@ BLOCKSET
Definition: definitions.h:217
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
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:637
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:507
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:445
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:678
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:48
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:336
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 _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:1775
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
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
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:92
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
CoreTmp< 0 > Core
Definition: Core.hpp:1129
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1943
const double D
diffusivity
DataForcesAndSourcesCore::EntData EntData
Fix dofs on entities.
virtual int get_comm_size() const =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.
virtual moab::Interface & get_moab()=0
Core (interface) class.
Definition: Core.hpp:77
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:60
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:100
Deprecated interface functions.
Data operator to do calculations at integration points.
static double volumeLengthQuality(const double *coords)
Calculate tetrahedron volume length quality.
Definition: Tools.cpp:34
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference 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.