};
};
GAUSS>::OpGradTensorGrad<1, SPACE_DIM, SPACE_DIM, 1>;
GAUSS>::OpSource<1, SPACE_DIM>;
private:
boost::shared_ptr<MatrixDouble> matGradPtr;
boost::shared_ptr<MatrixDouble> matStrainPtr;
boost::shared_ptr<MatrixDouble> matStressPtr;
boost::shared_ptr<MatrixDouble> matDPtr;
boost::shared_ptr<HenckyOps::CommonData> commonHenckyDataPtr;
};
auto set_material_stiffness = [&]() {
: 1;
auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(*matDPtr);
};
matGradPtr = boost::make_shared<MatrixDouble>();
matStrainPtr = boost::make_shared<MatrixDouble>();
matStressPtr = boost::make_shared<MatrixDouble>();
matDPtr = boost::make_shared<MatrixDouble>();
commonHenckyDataPtr = boost::make_shared<HenckyOps::CommonData>();
commonHenckyDataPtr->matGradPtr = matGradPtr;
commonHenckyDataPtr->matDPtr = matDPtr;
matDPtr->resize(size_symm * size_symm, 1);
CHKERR set_material_stiffness();
}
}
auto simple = mField.getInterface<Simple>();
}
Simple *
simple = mField.getInterface<Simple>();
}
auto *pipeline_mng = mField.getInterface<PipelineManager>();
auto simple = mField.getInterface<Simple>();
auto bc_mng = mField.getInterface<BcManager>();
CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_X",
"U", 0, 0);
CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_Y",
"U", 1, 1);
CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_Z",
"U", 2, 2);
CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_ALL",
"U", 0, 3);
auto get_time = [&](double, double, double) {
auto *pipeline_mng = mField.getInterface<PipelineManager>();
auto &fe_domain_rhs = pipeline_mng->getDomainRhsFE();
return fe_domain_rhs->ts_t;
};
if (it->getName().compare(0, 5, "FORCE") == 0) {
Range force_edges;
std::vector<double> attr_vec;
CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 1,
force_edges, true);
it->getAttributes(attr_vec);
"Wrong size of boundary attributes vector. Set right block "
"size attributes.");
auto force_vec_ptr = boost::make_shared<MatrixDouble>(
SPACE_DIM, 1);
std::copy(&attr_vec[0], &attr_vec[
SPACE_DIM],
force_vec_ptr->data().begin());
pipeline_mng->getOpBoundaryRhsPipeline().push_back(
boost::make_shared<Range>(force_edges)));
}
}
auto get_body_force = [this](const double, const double, const double) {
auto *pipeline_mng = mField.getInterface<PipelineManager>();
auto fe_domain_rhs = pipeline_mng->getDomainRhsFE();
const auto time = fe_domain_rhs->ts_t;
t_source(1) = 1 * time;
return t_source;
};
pipeline_mng->getOpDomainRhsPipeline().push_back(
}
auto *
simple = mField.getInterface<Simple>();
auto *pipeline_mng = mField.getInterface<PipelineManager>();
auto add_domain_base_ops = [&](auto &pipeline) {
auto det_ptr = boost::make_shared<VectorDouble>();
auto jac_ptr = boost::make_shared<MatrixDouble>();
auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
pipeline.push_back(new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
pipeline.push_back(new OpSetInvJacH1ForFace(inv_jac_ptr));
}
"U", matGradPtr));
pipeline.push_back(
new OpCalculateEigenVals<SPACE_DIM>("U", commonHenckyDataPtr));
pipeline.push_back(
new OpCalculateLogC<SPACE_DIM>("U", commonHenckyDataPtr));
pipeline.push_back(
new OpCalculateLogC_dC<SPACE_DIM>("U", commonHenckyDataPtr));
pipeline.push_back(
new OpCalculateHenckyStress<SPACE_DIM>("U", commonHenckyDataPtr));
pipeline.push_back(
new OpCalculatePiolaStress<SPACE_DIM>("U", commonHenckyDataPtr));
};
auto add_domain_ops_lhs = [&](auto &pipeline) {
pipeline.push_back(
new OpHenckyTangent<SPACE_DIM>("U", commonHenckyDataPtr));
pipeline.push_back(
new OpK(
"U",
"U", commonHenckyDataPtr->getMatTangent()));
};
auto add_domain_ops_rhs = [&](auto &pipeline) {
"U", commonHenckyDataPtr->getMatFirstPiolaStress()));
};
CHKERR add_domain_base_ops(pipeline_mng->getOpDomainLhsPipeline());
CHKERR add_domain_ops_lhs(pipeline_mng->getOpDomainLhsPipeline());
CHKERR add_domain_base_ops(pipeline_mng->getOpDomainRhsPipeline());
CHKERR add_domain_ops_rhs(pipeline_mng->getOpDomainRhsPipeline());
};
}
Monitor(SmartPetscObj<DM> dm, boost::shared_ptr<PostProcEle> post_proc)
: dM(dm), postProc(post_proc){};
"out_step_" + boost::lexical_cast<std::string>(ts_step) + ".h5m");
}
}
private:
SmartPetscObj<DM> dM;
boost::shared_ptr<PostProcEle> postProc;
};
auto *
simple = mField.getInterface<Simple>();
auto *pipeline_mng = mField.getInterface<PipelineManager>();
auto ts = pipeline_mng->createTSIM();
auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
post_proc_fe->generateReferenceElementMesh();
auto det_ptr = boost::make_shared<VectorDouble>();
auto jac_ptr = boost::make_shared<MatrixDouble>();
auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
post_proc_fe->getOpPtrVector().push_back(
post_proc_fe->getOpPtrVector().push_back(
new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
post_proc_fe->getOpPtrVector().push_back(
new OpSetInvJacH1ForFace(inv_jac_ptr));
}
post_proc_fe->getOpPtrVector().push_back(
"U", commonHenckyDataPtr->matGradPtr));
post_proc_fe->getOpPtrVector().push_back(
new OpCalculateEigenVals<SPACE_DIM>("U", commonHenckyDataPtr));
post_proc_fe->getOpPtrVector().push_back(
new OpCalculateLogC<SPACE_DIM>("U", commonHenckyDataPtr));
post_proc_fe->getOpPtrVector().push_back(
new OpCalculateLogC_dC<SPACE_DIM>("U", commonHenckyDataPtr));
post_proc_fe->getOpPtrVector().push_back(
new OpCalculateHenckyStress<SPACE_DIM>("U", commonHenckyDataPtr));
post_proc_fe->getOpPtrVector().push_back(
new OpCalculatePiolaStress<SPACE_DIM>("U", commonHenckyDataPtr));
post_proc_fe->getOpPtrVector().push_back(new OpPostProcHencky<SPACE_DIM>(
"U", post_proc_fe->postProcMesh, post_proc_fe->mapGaussPts,
commonHenckyDataPtr));
post_proc_fe->addFieldValuesPostProc("U");
boost::shared_ptr<FEMethod> null_fe;
auto monitor_ptr = boost::make_shared<Monitor>(dm, post_proc_fe);
null_fe, monitor_ptr);
double ftime = 1;
CHKERR TSSetDuration(ts, PETSC_DEFAULT, ftime);
CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
PetscInt steps, snesfails, rejects, nonlinits, linits;
CHKERR TSGetTimeStepNumber(ts, &steps);
CHKERR TSGetSNESFailures(ts, &snesfails);
CHKERR TSGetStepRejections(ts, &rejects);
CHKERR TSGetSNESIterations(ts, &nonlinits);
CHKERR TSGetKSPIterations(ts, &linits);
"steps %D (%D rejected, %D SNES fails), ftime %g, nonlinits "
"%D, linits %D\n",
steps, rejects, snesfails, ftime, nonlinits, linits);
}
PetscBool test_flg = PETSC_FALSE;
if (test_flg) {
auto *
simple = mField.getInterface<Simple>();
SCATTER_FORWARD);
double nrm2;
MOFEM_LOG(
"EXAMPLE", Sev::inform) <<
"Regression norm " << nrm2;
constexpr double regression_value = 1.09572;
if (fabs(nrm2 - regression_value) > 1e-2)
"Regression test faileed; wrong norm value.");
}
}
}
static char help[] =
"...\n\n";
int main(
int argc,
char *argv[]) {
auto core_log = logging::core::get();
core_log->add_sink(
try {
DMType dm_name = "DMMOFEM";
}
}
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_LOG_C(channel, severity, format,...)
void simple(double P1[], double P2[], double P3[], double c[], const int N)
int main(int argc, char *argv[])
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
Kronecker Delta class symmetric.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForce
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, SPACE_DIM > OpBodyForce
#define CATCH_ERRORS
Catch errors.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_ATOM_TEST_INVALID
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
auto smartCreateDMVector
Get smart vector from DM.
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
FaceElementForcesAndSourcesCoreSwitch< 0 > FaceElementForcesAndSourcesCore
Face finite element default.
VolumeElementForcesAndSourcesCoreSwitch< 0 > VolumeElementForcesAndSourcesCore
Volume finite element default.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
#define MOFEM_LOG(channel, severity)
Log.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
OpCalculateHOJacForFaceImpl< 2 > OpCalculateHOJacForFace
DeprecatedCoreInterface Interface
const double D
diffusivity
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesVector< 1, SPACE_DIM, 0 > OpBoundaryVec
static constexpr int approx_order
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
MoFEMErrorCode boundaryCondition()
[Set up problem]
MoFEMErrorCode assembleSystem()
[Boundary condition]
MoFEMErrorCode readMesh()
[Run problem]
MoFEMErrorCode checkResults()
[Postprocess results]
MoFEMErrorCode solveSystem()
[Solve]
MoFEMErrorCode createCommonData()
[Create common data]
MoFEMErrorCode runProblem()
[Create common data]
MoFEMErrorCode setupProblem()
[Read mesh]
MoFEMErrorCode outputResults()
[Solve]
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
default operator for TRI element
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
MoFEM::FaceElementForcesAndSourcesCore FaceEle
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
VolumeElementForcesAndSourcesCoreBase::UserDataOperator UserDataOperator
[Push operators to pipeline]
EntitiesFieldData::EntData EntData
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpGradTimesTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForce
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpGradTensorGrad< 1, SPACE_DIM, SPACE_DIM, 1 > OpK
constexpr double poisson_ratio
constexpr double shear_modulus_G
constexpr double bulk_modulus_K
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesVector< 1, SPACE_DIM, 0 > OpBoundaryVec
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, SPACE_DIM > OpBodyForce
constexpr double young_modulus