|
| v0.14.0
|
Definition at line 125 of file contact.cpp.
◆ Contact()
◆ bC()
[Create common data]
[Boundary condition]
Definition at line 435 of file contact.cpp.
440 for (
auto f : {
"U",
"SIGMA"}) {
441 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
442 "REMOVE_X",
f, 0, 0);
443 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
444 "REMOVE_Y",
f, 1, 1);
445 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
446 "REMOVE_Z",
f, 2, 2);
447 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
448 "REMOVE_ALL",
f, 0, 3);
451 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_X",
452 "SIGMA", 0, 0,
false,
true);
453 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_Y",
454 "SIGMA", 1, 1,
false,
true);
455 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_Z",
456 "SIGMA", 2, 2,
false,
true);
457 CHKERR bc_mng->removeBlockDOFsOnEntities(
simple->getProblemName(),
"FIX_ALL",
458 "SIGMA", 0, 3,
false,
true);
459 CHKERR bc_mng->removeBlockDOFsOnEntities(
460 simple->getProblemName(),
"NO_CONTACT",
"SIGMA", 0, 3,
false,
true);
465 simple->getProblemName(),
"U");
◆ checkResults()
[Solve]
[Check]
Definition at line 839 of file contact.cpp.
846 double analytical_active_area = 1.0;
848 double tol_force = 1e-3;
849 double tol_norm = 7.5;
850 double tol_area = 3e-2;
851 double fem_active_area = t_ptr[3];
856 fem_force = t_ptr[1];
861 fem_force = t_ptr[1];
868 fem_force = t_ptr[2];
869 analytical_active_area = M_PI / 4;
879 hertz_force = 15.873;
881 fem_force = t_ptr[1];
883 analytical_active_area = M_PI;
888 fem_force = t_ptr[1];
892 hertz_force = 0.5289;
893 fem_force = t_ptr[2];
898 "atom test %d does not exist",
atom_test);
900 if (fabs(fem_force - hertz_force) / hertz_force > tol_force) {
902 "atom test %d failed: Wrong FORCE output: %3.4e != %3.4e",
905 if (norm > tol_norm) {
907 "atom test %d failed: Wrong NORM output: %3.4e > %3.4e",
910 if (fabs(fem_active_area - analytical_active_area) > tol_area) {
912 "atom test %d failed: AREA computed %3.4e but should be %3.4e",
913 atom_test, fem_active_area, analytical_active_area);
◆ createCommonData()
[Set up problem]
[Create common data]
Definition at line 314 of file contact.cpp.
317 PetscBool use_mfront = PETSC_FALSE;
343 MOFEM_LOG(
"CONTACT", Sev::inform) <<
"Using MFront for material model";
354 PetscBool use_scale = PETSC_FALSE;
369 char sdf_file_name[255];
371 sdf_file_name, 255, PETSC_NULL);
373 sdfPythonPtr = boost::make_shared<SDFPython>();
374 CHKERR sdfPythonPtr->sdfInit(sdf_file_name);
375 sdfPythonWeakPtr = sdfPythonPtr;
381 "Use executable contact_2d with axisymmetric model");
385 "Axisymmetric model is only available with MFront (set "
388 MOFEM_LOG(
"CONTACT", Sev::inform) <<
"Using axisymmetric model";
393 MOFEM_LOG(
"CONTACT", Sev::inform) <<
"Using plane strain model";
398 #ifndef WITH_MODULE_MFRONT_INTERFACE
401 "MFrontInterface module was not found while use_mfront was set to 1");
405 boost::make_shared<MFrontMoFEMInterface<TRIDIMENSIONAL>>(
410 boost::make_shared<MFrontMoFEMInterface<AXISYMMETRICAL>>(
413 mfrontInterface = boost::make_shared<MFrontMoFEMInterface<PLANESTRAIN>>(
422 auto dm =
simple->getDM();
◆ OPs()
[Boundary condition]
[Push operators to pip]
[Only used for dynamics]
[Only used for dynamics]
[Only used for dynamics]
[Only used for dynamics]
[Operators used for contact]
[Operators used for contact]
[Operators used for contact]
[Operators used for contact]
Definition at line 472 of file contact.cpp.
477 auto time_scale = boost::make_shared<ScaledTimeScale>();
478 auto body_force_time_scale =
479 boost::make_shared<ScaledTimeScale>(
"body_force_hist.txt");
488 auto add_domain_base_ops = [&](
auto &pip) {
495 auto henky_common_data_ptr = boost::make_shared<HenckyOps::CommonData>();
496 henky_common_data_ptr->matDPtr = boost::make_shared<MatrixDouble>();
497 henky_common_data_ptr->matGradPtr = boost::make_shared<MatrixDouble>();
499 auto add_domain_ops_lhs = [&](
auto &pip) {
509 auto fe_domain_lhs = pip_mng->getDomainLhsFE();
511 auto get_inertia_and_mass_damping =
513 return (
rho *
scale) * fe_domain_lhs->ts_aa +
516 pip.push_back(
new OpMass(
"U",
"U", get_inertia_and_mass_damping));
520 auto fe_domain_lhs = pip_mng->getDomainLhsFE();
522 auto get_inertia_and_mass_damping =
526 pip.push_back(
new OpMass(
"U",
"U", get_inertia_and_mass_damping));
530 CHKERR HenckyOps::opFactoryDomainLhs<SPACE_DIM, AT, IT, DomainEleOp>(
531 mField, pip,
"U",
"MAT_ELASTIC", Sev::verbose,
scale);
539 auto add_domain_ops_rhs = [&](
auto &pip) {
543 pip,
mField,
"U", {body_force_time_scale}, Sev::inform);
547 AT>::LinearForm<IT>::OpBaseTimesVector<1,
SPACE_DIM, 1>;
552 auto mat_acceleration = boost::make_shared<MatrixDouble>();
554 "U", mat_acceleration));
556 new OpInertiaForce(
"U", mat_acceleration, [](
double,
double,
double) {
563 auto mat_velocity = boost::make_shared<MatrixDouble>();
567 new OpInertiaForce(
"U", mat_velocity, [](
double,
double,
double) {
573 CHKERR HenckyOps::opFactoryDomainRhs<SPACE_DIM, AT, IT, DomainEleOp>(
574 mField, pip,
"U",
"MAT_ELASTIC", Sev::inform,
scale);
579 CHKERR ContactOps::opFactoryDomainRhs<SPACE_DIM, AT, IT, DomainEleOp>(
585 auto add_boundary_base_ops = [&](
auto &pip) {
595 auto add_boundary_ops_lhs = [&](
auto &pip) {
605 pip,
mField,
"U", Sev::inform);
610 auto fe_boundary_lhs = pip_mng->getBoundaryLhsFE();
615 [
this, fe_boundary_lhs](
double,
double,
double) {
624 ContactOps::opFactoryBoundaryLhs<SPACE_DIM, AT, GAUSS, BoundaryEleOp>(
628 mField, pip,
simple->getDomainFEName(),
"SIGMA",
"U",
"GEOMETRY",
634 auto add_boundary_ops_rhs = [&](
auto &pip) {
639 AT>::LinearForm<IT>::OpBaseTimesVector<1,
SPACE_DIM, 1>;
644 pip,
mField,
"U", {time_scale}, Sev::inform);
647 auto u_disp = boost::make_shared<MatrixDouble>();
648 auto dot_u_disp = boost::make_shared<MatrixDouble>();
653 new OpSpringRhs(
"U", u_disp, [
this](
double,
double,
double) {
657 new OpSpringRhs(
"U", dot_u_disp, [
this](
double,
double,
double) {
663 ContactOps::opFactoryBoundaryRhs<SPACE_DIM, AT, GAUSS, BoundaryEleOp>(
669 CHKERR add_domain_base_ops(pip_mng->getOpDomainLhsPipeline());
670 CHKERR add_domain_base_ops(pip_mng->getOpDomainRhsPipeline());
671 CHKERR add_domain_ops_lhs(pip_mng->getOpDomainLhsPipeline());
672 CHKERR add_domain_ops_rhs(pip_mng->getOpDomainRhsPipeline());
674 CHKERR add_boundary_base_ops(pip_mng->getOpBoundaryLhsPipeline());
675 CHKERR add_boundary_base_ops(pip_mng->getOpBoundaryRhsPipeline());
676 CHKERR add_boundary_ops_lhs(pip_mng->getOpBoundaryLhsPipeline());
677 CHKERR add_boundary_ops_rhs(pip_mng->getOpBoundaryRhsPipeline());
683 CHKERR pip_mng->setDomainRhsIntegrationRule(integration_rule_vol);
684 CHKERR pip_mng->setDomainLhsIntegrationRule(integration_rule_vol);
685 CHKERR pip_mng->setBoundaryRhsIntegrationRule(integration_rule_boundary);
686 CHKERR pip_mng->setBoundaryLhsIntegrationRule(integration_rule_boundary);
◆ runProblem()
◆ setupProblem()
[Run problem]
[Set up problem]
< blocs interation is collective, so that is set irrespective if there are entities in given rank or not in the block
Definition at line 174 of file contact.cpp.
193 enum bases { AINSWORTH, DEMKOWICZ, LASBASETOPT };
194 const char *list_bases[LASBASETOPT] = {
"ainsworth",
"demkowicz"};
195 PetscInt choice_base_value = AINSWORTH;
197 LASBASETOPT, &choice_base_value, PETSC_NULL);
200 switch (choice_base_value) {
204 <<
"Set AINSWORTH_LEGENDRE_BASE for displacements";
209 <<
"Set DEMKOWICZ_JACOBI_BASE for displacements";
229 auto get_skin = [&]() {
234 CHKERR skin.find_skin(0, body_ents,
false, skin_ents);
238 auto filter_blocks = [&](
auto skin) {
239 bool is_contact_block =
false;
244 (boost::format(
"%s(.*)") %
"CONTACT").str()
253 <<
"Find contact block set: " <<
m->getName();
254 auto meshset =
m->getMeshset();
255 Range contact_meshset_range;
257 meshset,
SPACE_DIM - 1, contact_meshset_range,
true);
260 contact_meshset_range);
261 contact_range.merge(contact_meshset_range);
263 if (is_contact_block) {
265 <<
"Nb entities in contact surface: " << contact_range.size();
267 skin = intersect(skin, contact_range);
272 auto filter_true_skin = [&](
auto skin) {
274 ParallelComm *pcomm =
276 CHKERR pcomm->filter_pstatus(skin, PSTATUS_SHARED | PSTATUS_MULTISHARED,
277 PSTATUS_NOT, -1, &boundary_ents);
278 return boundary_ents;
281 auto boundary_ents = filter_true_skin(filter_blocks(get_skin()));
289 moab::Interface::UNION);
298 auto project_ho_geometry = [&]() {
303 PetscBool project_geometry = PETSC_TRUE;
305 &project_geometry, PETSC_NULL);
306 if (project_geometry) {
307 CHKERR project_ho_geometry();
◆ tsSolve()
Definition at line 702 of file contact.cpp.
709 auto set_section_monitor = [&](
auto solver) {
712 CHKERR TSGetSNES(solver, &snes);
713 PetscViewerAndFormat *vf;
714 CHKERR PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD,
715 PETSC_VIEWER_DEFAULT, &vf);
718 (
MoFEMErrorCode(*)(SNES, PetscInt, PetscReal,
void *))SNESMonitorFields,
723 auto scatter_create = [&](
auto D,
auto coeff) {
725 CHKERR is_manager->isCreateProblemFieldAndRank(
simple->getProblemName(),
726 ROW,
"U", coeff, coeff, is);
728 CHKERR ISGetLocalSize(is, &loc_size);
732 CHKERR VecScatterCreate(
D, is,
v, PETSC_NULL, &scatter);
737 auto set_time_monitor = [&](
auto dm,
auto solver) {
740 boost::shared_ptr<ForcesAndSourcesCore>
null;
746 auto set_essential_bc = [&]() {
750 auto pre_proc_ptr = boost::make_shared<FEMethod>();
751 auto post_proc_rhs_ptr = boost::make_shared<FEMethod>();
752 auto post_proc_lhs_ptr = boost::make_shared<FEMethod>();
755 auto time_scale = boost::make_shared<TimeScale>();
757 auto get_bc_hook_rhs = [&]() {
759 {time_scale},
false);
762 pre_proc_ptr->preProcessHook = get_bc_hook_rhs();
764 auto get_post_proc_hook_rhs = [&]() {
766 mField, post_proc_rhs_ptr, 1.);
768 auto get_post_proc_hook_lhs = [&]() {
770 mField, post_proc_lhs_ptr, 1.);
772 post_proc_rhs_ptr->postProcessHook = get_post_proc_hook_rhs();
774 ts_ctx_ptr->getPreProcessIFunction().push_front(pre_proc_ptr);
775 ts_ctx_ptr->getPreProcessIJacobian().push_front(pre_proc_ptr);
776 ts_ctx_ptr->getPostProcessIFunction().push_back(post_proc_rhs_ptr);
777 post_proc_lhs_ptr->postProcessHook = get_post_proc_hook_lhs();
778 ts_ctx_ptr->getPostProcessIJacobian().push_back(post_proc_lhs_ptr);
783 auto set_schur_pc = [&](
auto solver) {
784 boost::shared_ptr<SetUpSchur> schur_ptr;
793 auto dm =
simple->getDM();
805 CHKERR set_essential_bc();
808 auto solver = pip_mng->createTSIM();
809 CHKERR TSSetFromOptions(solver);
810 auto schur_pc_ptr = set_schur_pc(solver);
813 CHKERR set_section_monitor(solver);
814 CHKERR set_time_monitor(dm, solver);
815 CHKERR TSSetSolution(solver,
D);
817 CHKERR TSSolve(solver, NULL);
819 auto solver = pip_mng->createTSIM2();
820 CHKERR TSSetFromOptions(solver);
821 auto schur_pc_ptr = set_schur_pc(solver);
823 auto dm =
simple->getDM();
827 CHKERR set_section_monitor(solver);
828 CHKERR set_time_monitor(dm, solver);
829 CHKERR TS2SetSolution(solver,
D, DD);
831 CHKERR TSSolve(solver, NULL);
◆ mField
◆ mfrontInterface
◆ monitorPtr
boost::shared_ptr<Monitor> Contact::monitorPtr |
|
private |
◆ uXScatter
◆ uYScatter
◆ uZScatter
The documentation for this struct was generated from the following file:
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
#define MYPCOMM_INDEX
default communicator number PCOMM
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
Specialization for DisplacementCubitBcData.
virtual MPI_Comm & get_comm() const =0
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Get values at integration pts for tensor filed rank 1, i.e. vector field.
virtual int get_comm_rank() const =0
Specialization for DisplacementCubitBcData.
PipelineManager interface.
Definition of the displacement bc data structure.
Projection of edge entities with one mid-node on hierarchical basis.
Simple interface for fast problem set-up.
Approximate field values for given petsc vector.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesVector< 1, SPACE_DIM, 1 > OpInertiaForce
#define CHKERR
Inline error check.
auto createDMVector(DM dm)
Get smart vector from DM.
Specialization for DisplacementCubitBcData.
virtual moab::Interface & get_moab()=0
Section manager is used to create indexes and sections.
Simple interface for fast problem set-up.
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
[Only used with Hooke equation (linear material model)]
void simple(double P1[], double P2[], double P3[], double c[], const int N)
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Add operators pushing bases from local to physical configuration.
const double v
phase velocity of light in medium (cm/ns)
#define MOFEM_LOG(channel, severity)
Log.
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
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.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
const FTensor::Tensor2< T, Dim, Dim > Vec
@ MOFEM_DATA_INCONSISTENCY
PetscErrorCode PetscOptionsGetString(PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
Interface for managing meshsets containing materials and boundary conditions.
PetscErrorCode PetscOptionsGetEList(PetscOptions *, const char pre[], const char name[], const char *const *list, PetscInt next, PetscInt *value, PetscBool *set)
FieldApproximationBase
approximation base
const double D
diffusivity
FTensor::Index< 'm', 3 > m
@ MOFEM_ATOM_TEST_INVALID
PetscErrorCode PetscOptionsGetScalar(PetscOptions *, const char pre[], const char name[], PetscScalar *dval, PetscBool *set)
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
static boost::shared_ptr< SetUpSchur > createSetUpSchur(MoFEM::Interface &m_field)
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ HDIV
field with continuous normal traction
auto getDMTsCtx(DM dm)
Get TS context data structure used by DM.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)