13#ifndef __BASICBOUNDARYCONDIONSINTERFACE_HPP__
14#define __BASICBOUNDARYCONDIONSINTERFACE_HPP__
22 using DomainEle = VolumeElementForcesAndSourcesCore;
57 :
sCale(
scale), TimeScaleVector3(file_name, false) {}
117 string mesh_pos_field_name =
"MESH_NODE_POSITIONS",
118 string problem_name =
"ELASTIC",
119 string domain_element_name =
"ELASTIC_FE",
121 double *snes_load_factor =
nullptr,
bool is_partitioned =
true)
135 PetscBool quasi_static = PETSC_FALSE;
136 PetscBool is_linear = PETSC_FALSE;
137 CHKERR PetscOptionsGetBool(PETSC_NULLPTR,
"-is_quasi_static", &quasi_static,
139 CHKERR PetscOptionsGetBool(PETSC_NULLPTR,
"-is_linear", &is_linear,
160 dirichletBcPtr = boost::make_shared<DirichletSpatialRemoveDofsBc>(
164 dirichletBcPtr = boost::make_shared<DirichletDisplacementRemoveDofsBc>(
220 boost::make_shared<VolumeElementForcesAndSourcesCore>(
mField);
222 boost::make_shared<VolumeElementForcesAndSourcesCore>(
mField);
238 auto get_id_block_param = [&](
string base_name,
int id) {
239 char load_hist_file[255] =
"hist.in";
240 PetscBool ctg_flag = PETSC_FALSE;
241 string param_name_with_id =
"-" + base_name +
"_" + to_string(
id);
242 CHKERR PetscOptionsGetString(PETSC_NULLPTR, PETSC_NULLPTR,
243 param_name_with_id.c_str(), load_hist_file,
246 return param_name_with_id;
248 param_name_with_id =
"-" + base_name;
249 CHKERR PetscOptionsGetString(PETSC_NULLPTR, PETSC_NULLPTR,
250 param_name_with_id.c_str(), load_hist_file,
254 <<
"Setting one accelerogram for all blocks!";
255 return param_name_with_id;
261 auto get_adj_ents = [&](
const Range &ents) {
264 for (
size_t d = 1; d < 3; ++d)
266 moab::Interface::UNION);
273 const std::string block_name =
"BODY_FORCE";
274 if (it->getName().compare(0, block_name.size(), block_name) == 0) {
275 std::vector<double> attr;
276 CHKERR it->getAttributes(attr);
277 if (attr.size() > 3) {
279 const int id = it->getMeshsetId();
284 auto bc_ents_ptr = boost::make_shared<Range>(get_adj_ents(bc_ents));
288 VectorDouble acc({attr[1], attr[2], attr[3]});
289 double density = attr[0];
291 attr.size() > 4 ?
bool(std::floor(attr[4])) :
true;
294 std::vector<boost::shared_ptr<TimeScaleVector3>> methods_for_scaling;
295 std::string param_name_for_scaling =
296 get_id_block_param(
"accelerogram",
id);
298 if (!param_name_for_scaling.empty())
299 methods_for_scaling.push_back(
300 boost::make_shared<BasicBCVectorScale>(density,
301 param_name_for_scaling));
303 methods_for_scaling.push_back(
304 boost::make_shared<BasicBCVectorConst>(
312 return density * fe_domain_lhs->ts_aa;
324 pipeline_rhs.push_back(
330 pipeline_rhs,
mField,
"U", {}, methods_for_scaling,
331 it->getName(), Sev::inform);
334 pipeline_lhs.push_back(
336 pipeline_lhs.push_back(
338 auto mat_acceleration = boost::make_shared<MatrixDouble>();
339 pipeline_rhs.push_back(
new OpCalculateVectorFieldValuesDotDot<3>(
343 [&](
double,
double,
double) {
return density; }, bc_ents_ptr));
351 "There should be (1 density + 3 accelerations ) attributes in "
352 "BODY_FORCE blockset, but is %ld. Optionally, you can set 5th "
353 "parameter to inertia flag.",
359 auto integration_rule_vol = [](int, int,
int approx_order) {
362 auto integration_rule_boundary = [](int, int,
int approx_order) {
400 bc_mng->getMergedBlocksMarker(vector<string>{
"FIX_",
"ROTATE"});
409 vector<const char *> element_list{
"FORCE_FE",
"PRESSURE_FE",
413 for (
auto &el : element_list) {
415 simple->getOtherFiniteElements().push_back(el);
430 char load_hist_file[255] =
"hist.in";
431 PetscBool ctg_flag = PETSC_FALSE;
432 string new_param_file = string(
"-") + prefix + string(
"_history");
433 CHKERR PetscOptionsGetString(PETSC_NULLPTR, PETSC_NULLPTR, new_param_file.c_str(),
434 load_hist_file, 255, &ctg_flag);
436 return new_param_file;
437 return string(
"-load_history");
440 template <
typename T>
446 template <
typename T>
452 template <
typename T,
bool RHS>
458 auto set_solver_pipelines =
459 [&](PetscErrorCode (*function)(DM,
const char fe_name[],
MoFEM::FEMethod *,
462 PetscErrorCode (*jacobian)(DM,
const char fe_name[],
MoFEM::FEMethod *,
468 if (std::is_same_v<T, TS>)
478 auto set_neumann_methods = [&](
auto &neumann_el,
string hist_name,
481 for (
auto &&mit : neumann_el) {
482 if constexpr (std::is_same_v<T, SNES>)
483 mit->second->methodsOp.push_back(
485 if constexpr (std::is_same_v<T, TS>)
486 mit->second->methodsOp.push_back(
488 string element_name = mit->first;
492 mit->second->getLoopFe().getOpPtrVector(), {},
497 mit->second->getLoopFe().getOpPtrVector(), {},
507 CHKERR function(
dM, element_name.c_str(),
508 &mit->second->getLoopFe(), NULL, NULL);
524 CHKERR function(
dM,
"FLUID_PRESSURE_FE",
541 if constexpr (std::is_same_v<T, SNES>) {
545 "SNES lambda factor pointer not set in the module constructor");
547 CHKERR set_solver_pipelines(&DMMoFEMSNESSetFunction,
548 &DMMoFEMSNESSetJacobian);
550 }
else if constexpr (std::is_same_v<T, TS>) {
554 CHKERR set_solver_pipelines(&DMMoFEMTSSetIFunction,
555 &DMMoFEMTSSetIJacobian);
558 CHKERR set_solver_pipelines(&DMMoFEMTSSetI2Function,
559 &DMMoFEMTSSetI2Jacobian);
562 CHKERR set_solver_pipelines(&DMMoFEMTSSetRHSFunction,
563 &DMMoFEMTSSetRHSJacobian);
567 "This TS is not yet implemented for basic BCs");
572 static_assert(!std::is_same_v<T, KSP>,
573 "this solver has not been implemented for basic BCs yet");
void simple(double P1[], double P2[], double P3[], double c[], const int N)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
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 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 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_row(const std::string &fe_name, const std::string name_row)=0
set field row which finite element use
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_field)=0
set finite element field data
#define MOFEM_LOG(channel, severity)
Log.
#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
PetscBool is_quasi_static
static constexpr int approx_order
FTensor::Tensor1< double, 3 > getVector(const double time)
FTensor::Tensor1< double, 3 > tForce
BasicBCVectorConst(double scale, FTensor::Tensor1< double, 3 > t_vec)
BasicBCVectorScale(double scale, std::string file_name)
FTensor::Tensor1< double, 3 > getVector(const double time)
MoFEMErrorCode scaleNf(const FEMethod *fe, VectorDouble &nf)
LoadScale(double *my_lambda)
Set of functions declaring elements and setting operators for basic boundary conditions interface.
MoFEMErrorCode setupSolverFunctionSNES() override
BasicBoundaryConditionsInterface(MoFEM::Interface &m_field, string postion_field, string mesh_pos_field_name="MESH_NODE_POSITIONS", string problem_name="ELASTIC", string domain_element_name="ELASTIC_FE", bool is_displacement_field=true, bool is_quasi_static=true, double *snes_load_factor=nullptr, bool is_partitioned=true)
double * snesLambdaLoadFactorPtr
boost::shared_ptr< FaceElementForcesAndSourcesCore > springLhsPtr
MoFEMErrorCode getCommandLineParameters() override
MoFEMErrorCode addElementFields() override
MoFEMErrorCode setupSolverFunction(const TSType type=IM)
~BasicBoundaryConditionsInterface()
const string domainElementName
MoFEMErrorCode setupSolverJacobianTS(const TSType type) override
MoFEM::Interface & mField
boost::shared_ptr< VolumeElementForcesAndSourcesCore > bodyForceLhsPtr
boost::shared_ptr< FluidPressure > fluidPressureElementPtr
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, 3 > OpMass
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, 3 > OpBodyForce
boost::shared_ptr< FaceElementForcesAndSourcesCore > springRhsPtr
boost::ptr_map< std::string, EdgeForce > edge_forces
MoFEMErrorCode setupSolverImpl(const TSType type=IM)
MoFEMErrorCode setOperators() override
boost::ptr_map< std::string, NodalForce > nodal_forces
MoFEMErrorCode createElements() override
string getHistoryParam(string prefix)
boost::ptr_map< std::string, NeumannForcesSurface > neumann_forces
MoFEMErrorCode setupSolverFunctionTS(const TSType type) override
boost::shared_ptr< VolumeElementForcesAndSourcesCore > bodyForceRhsPtr
MoFEMErrorCode addElementsToDM(SmartPetscObj< DM > dm) override
MoFEMErrorCode setupSolverJacobian(const TSType type=IM)
boost::shared_ptr< KelvinVoigtDamper > damperElementPtr
boost::shared_ptr< DirichletDisplacementBc > dirichletBcPtr
NaturalBC< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS > DomainNaturalBC
const string domainProblemName
MoFEMErrorCode updateElementVariables() override
VolumeElementForcesAndSourcesCore DomainEle
MoFEMErrorCode postProcessElement(int step) override
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesVector< 1, 3, 1 > OpInertiaForce
MoFEMErrorCode setupSolverJacobianSNES() override
Set of functions declaring elements and setting operators for generic element interface.
BcMarkerPtr mBoundaryMarker
Constitutive model functions.
Class used to scale loads, f.e. in arc-length control.
Data structure to exchange data between mofem and User Loop Methods.
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
structure for User Loop Methods on finite elements
Force scale operator for reading four columns (time and vector)
virtual FTensor::Tensor1< double, SPACE_DIM > getVector(const double time)
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Force scale operator for reading two columns.