11using namespace boost::numeric;
31 boost::shared_ptr<MatrixDouble>
xAtPts =
45 CHKERRABORT(PETSC_COMM_WORLD,
ierr);
50 PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"Problem",
"none");
70 if (
bit->getName().compare(0, 3,
"ROD") == 0) {
72 const int id =
bit->getMeshsetId();
77 std::vector<double> attributes;
78 bit->getAttributes(attributes);
79 if (attributes.size() < 3) {
82 "Input mesh for ROD should have 3 attributes but there is %ld",
91 CHKERR PetscPrintf(PETSC_COMM_WORLD,
"\nSimple rod block %d\n",
id);
92 CHKERR PetscPrintf(PETSC_COMM_WORLD,
"\tYoung's modulus %3.4g\n",
94 CHKERR PetscPrintf(PETSC_COMM_WORLD,
"\tCross-section area %3.4g\n",
96 CHKERR PetscPrintf(PETSC_COMM_WORLD,
"\tPrestress %3.4g\n",
120 boost::shared_ptr<DataAtIntegrationPtsSimpleRods> &common_data_ptr,
135 const int row_nb_dofs = row_data.
getIndices().size();
139 const int col_nb_dofs = col_data.
getIndices().size();
149 locK.resize(row_nb_dofs, col_nb_dofs,
false);
152 double tension_stiffness =
commonDataPtr->simpleRodYoungModulus *
158 double coeff = tension_stiffness / (
L *
L *
L);
160 double x21 = coords(3) - coords(0);
161 double y21 = coords(4) - coords(1);
162 double z21 = coords(5) - coords(2);
165 locK(0, 0) = coeff * x21 * x21;
166 locK(0, 1) = coeff * x21 * y21;
167 locK(0, 2) = coeff * x21 * z21;
168 locK(0, 3) = -coeff * x21 * x21;
169 locK(0, 4) = -coeff * x21 * y21;
170 locK(0, 5) = -coeff * x21 * z21;
173 locK(1, 1) = coeff * y21 * y21;
174 locK(1, 2) = coeff * y21 * z21;
175 locK(1, 3) = -coeff * y21 * x21;
176 locK(1, 4) = -coeff * y21 * y21;
177 locK(1, 5) = -coeff * y21 * z21;
181 locK(2, 2) = coeff * z21 * z21;
182 locK(2, 3) = -coeff * z21 * x21;
183 locK(2, 4) = -coeff * z21 * y21;
184 locK(2, 5) = -coeff * z21 * z21;
189 locK(3, 3) = coeff * x21 * x21;
190 locK(3, 4) = coeff * x21 * y21;
191 locK(3, 5) = coeff * x21 * z21;
197 locK(4, 4) = coeff * y21 * y21;
198 locK(4, 5) = coeff * y21 * z21;
205 locK(5, 5) = coeff * z21 * z21;
225 boost::shared_ptr<DataAtIntegrationPtsSimpleRods> &common_data_ptr,
250 nF.resize(nb_dofs,
false);
258 for (
auto d : {0, 1, 2}) {
259 nF(d) = -axial_force * dir[d];
260 nF(d + 3) = axial_force * dir[d];
271 const std::string mesh_nodals_positions) {
282 mesh_nodals_positions);
287 if (
bit->getName().compare(0, 3,
"ROD") == 0) {
289 MBEDGE,
"SIMPLE_ROD");
299 boost::shared_ptr<EdgeElementForcesAndSourcesCore>
300 fe_simple_rod_lhs_ptr,
301 boost::shared_ptr<EdgeElementForcesAndSourcesCore>
302 fe_simple_rod_rhs_ptr,
303 const std::string
field_name,
const std::string mesh_nodals_positions) {
308 boost::shared_ptr<DataAtIntegrationPtsSimpleRods> commonDataPtr =
309 boost::make_shared<DataAtIntegrationPtsSimpleRods>(m_field);
310 CHKERR commonDataPtr->getParameters();
312 for (
auto &sitSimpleRod : commonDataPtr->mapSimpleRod) {
313 fe_simple_rod_lhs_ptr->getOpPtrVector().push_back(
317 commonDataPtr, sitSimpleRod.second,
field_name));
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#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 ...
@ MOFEM_ATOM_TEST_INVALID
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
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 build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elements.
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
virtual bool check_field(const std::string &name) const =0
check if field is in database
#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.
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasMatrix< double > MatrixDouble
implementation of Data Operators for Forces and Sources
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
constexpr auto field_name
double simpleRodYoungModulus
BlockOptionDataSimpleRods()
double simpleRodSectionArea
double simpleRodPreStress
double simpleRodYoungModulus
DataAtIntegrationPtsSimpleRods(MoFEM::Interface &m_field)
MoFEMErrorCode getBlockData(BlockOptionDataSimpleRods &data)
boost::shared_ptr< MatrixDouble > xAtPts
double simpleRodPreStress
double simpleRodSectionArea
MoFEMErrorCode setBlocks()
MoFEM::Interface & mField
MoFEMErrorCode getParameters()
boost::shared_ptr< MatrixDouble > gradDispPtr
boost::shared_ptr< MatrixDouble > xInitAtPts
std::map< int, BlockOptionDataSimpleRods > mapSimpleRod
virtual moab::Interface & get_moab()=0
bool sYmm
If true assume that matrix is symmetric structure.
Deprecated interface functions.
default operator for EDGE element
VectorDouble & getCoords()
get edge node coordinates
double getLength()
get edge length
VectorDouble & getDirection()
get edge direction
Data on single entity (This is passed as argument to DataOperator::doWork)
const VectorInt & getIndices() const
Get global indices of dofs on entity.
EntityHandle getFEEntityHandle() const
Return finite element entity handle.
@ OPROW
operator doWork function is executed on FE rows
@ OPROWCOL
operator doWork is executed on FE rows &columns
Assemble contribution of SimpleRod element to LHS *.
BlockOptionDataSimpleRods & dAta
boost::shared_ptr< DataAtIntegrationPtsSimpleRods > commonDataPtr
OpSimpleRodK(boost::shared_ptr< DataAtIntegrationPtsSimpleRods > &common_data_ptr, BlockOptionDataSimpleRods &data, const std::string field_name)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
Add ROD pre-stress to the RHS *.
OpSimpleRodPreStress(boost::shared_ptr< DataAtIntegrationPtsSimpleRods > &common_data_ptr, BlockOptionDataSimpleRods &data, const std::string field_name)
boost::shared_ptr< DataAtIntegrationPtsSimpleRods > commonDataPtr
BlockOptionDataSimpleRods & dAta
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.