v0.14.0
COR-8: Implementation of spring element

Introduction

This tutorial presents the implementation of spring elements in MoFEM. The spring elements are used to model the interaction between the structure being considered and other structures in surrounding environment. In contrast to other type of boundary conditions where the displacements/rotations at the surface of contact are fixed, the spring boundary condition allows the surface to move during the analysis which in turn creates reaction forces acting on the structure at the boundary surface. The magnitude of this force is changed propotionally to the spring stiffness and the extension or contraction in its length. Figure 1 illustrates a plate bearing load \( q \) sits on elastic foundation which is modelled by springs with stiffness \( k_s \).

Figure 1: Plate on elastic foundation modelled by springs

The remaining of this tutorial is organised as follows. The next section describes the finite element formulation of the spring element where its contributions to the right-hand side (RHS) and the left-hand side (LHS) are presented. It is followed by the description of the implementation of spring element in the platform of MoFEM. The following section closes this tutorial by showing sample implementations of spring elements to the elasticity and nonlinear elasticity modules.

Formulation

The finite element formulation of the spring element contributing to RHS and LHS can be respectively derived as follows

\[ \begin{align} f_s &= \int\limits_{\Gamma }^{} {{\psi ^T}{F^s}\left( u \right)d\Gamma } = \int\limits_{\Gamma }^{} {{\psi ^T}{k_s} u d \Gamma }, \\ {K^s} &= \int\limits_\Gamma ^{} {{\psi ^T}{k_s}\psi d\Gamma }. \end{align} \]

where \( \Gamma \) is the domain occupied by the surface with springs, \( \psi \) is approximation function, \( k_s \) denotes spring stiffness, and \( u \) is the extension/contraction of the spring which can be calculated by the displacement of the point on the surface where spring is applied. Physically, \( k_s \) can be considered as spring density on the surface.

Implementation

In the platform of MoFEM, the implementation of user data operators (UDO) for spring boundary conditions can be found in SpringElement.cpp and SpringElement.hpp. Particularly, the operators OpSpringFs and OpSpringKs are developed to calculate the contribution of springs to the RHS and LHS of the system of equations, respectively.

Operator OpSpringFs for the RHS

  • The operator OpSpringFs is a class that inherits from another class of MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator. It is worth noting that, in the current implementation, springs are applied on face elements only.
  • Public members include (1) nF is a sub-vector that stores force at each degree of freedom, (2) commonDataPtr pointing to spring data associated with Gauss points, (3) dAta contains block data, and (4) is_spatial_position being a flag to check the field name, e.g. SPATIAL_POSITION or DISPLACEMENT, which later be used to determine the displacements at Gauss points.
    // vector used to store force vector for each degree of freedom
    boost::shared_ptr<DataAtIntegrationPtsSprings> commonDataPtr;
    BlockOptionDataSprings &dAta;
    bool is_spatial_position = true;
  • Then, the OpSpringFs constructor is implemented
    OpSpringFs(boost::shared_ptr<DataAtIntegrationPtsSprings> &common_data_ptr,
    BlockOptionDataSprings &data, const std::string field_name)
    field_name.c_str(), OPROW),
    commonDataPtr(common_data_ptr), dAta(data) {
    if (field_name.compare(0, 16, "SPATIAL_POSITION") != 0)
    is_spatial_position = false;
    }
  • Finally, the OpSpringFs::doWork method of the OpSpringFs operator is implemented as follows
    MoFEMErrorCode doWork(int side, EntityType type,
    // check that the faces have associated degrees of freedom
    const int nb_dofs = data.getIndices().size();
    if (nb_dofs == 0)
    if (dAta.tRis.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
    dAta.tRis.end()) {
    }
    CHKERR commonDataPtr->getBlockData(dAta);
    // size of force vector associated to the entity
    // set equal to the number of degrees of freedom of associated with the
    // entity
    nF.resize(nb_dofs, false);
    nF.clear();
    // get number of Gauss points
    const int nb_gauss_pts = data.getN().size1();
    // get intergration weights
    auto t_w = getFTensor0IntegrationWeight();
    // FTensor indices
    commonDataPtr->springStiffnessNormal, 0., 0., 0.,
    commonDataPtr->springStiffnessTangent, 0., 0., 0.,
    commonDataPtr->springStiffnessTangent);
    // create a 3d vector to be used as the normal to the face with length equal
    // to the face area
    auto t_normal_ptr = getFTensor1Normal();
    t_normal(i) = t_normal_ptr(i);
    // First tangent vector
    auto t_tangent1_ptr = getFTensor1Tangent1AtGaussPts();
    t_tangent1(i) = t_tangent1_ptr(i);
    // Second tangent vector, such that t_n = t_t1 x t_t2 | t_t2 = t_n x t_t1
    t_tangent2(i) = FTensor::levi_civita(i, j, k) * t_normal(j) * t_tangent1(k);
    // Spring stiffness in global coordinate
    t_spring_global = MetaSpringBC::transformLocalToGlobal(
    t_normal, t_tangent1, t_tangent2, t_spring_local);
    // Extract solution at Gauss points
    auto t_solution_at_gauss_point =
    getFTensor1FromMat<3>(*commonDataPtr->xAtPts);
    auto t_init_solution_at_gauss_point =
    getFTensor1FromMat<3>(*commonDataPtr->xInitAtPts);
    FTensor::Tensor1<double, 3> t_displacement_at_gauss_point;
    // loop over all Gauss points of the face
    for (int gg = 0; gg != nb_gauss_pts; ++gg) {
    // Calculate the displacement at the Gauss point
    if (is_spatial_position) { // "SPATIAL_POSITION"
    t_displacement_at_gauss_point(i) =
    t_solution_at_gauss_point(i) - t_init_solution_at_gauss_point(i);
    } else { // e.g. "DISPLACEMENT" or "U"
    t_displacement_at_gauss_point(i) = t_solution_at_gauss_point(i);
    }
    double w = t_w * getArea();
    auto t_base_func = data.getFTensor0N(gg, 0);
    // create a vector t_nf whose pointer points an array of 3 pointers
    // pointing to nF memory location of components
    &nF[2]);
    for (int rr = 0; rr != nb_dofs / 3; ++rr) { // loop over the nodes
    t_nf(i) += w * t_base_func * t_spring_global(i, j) *
    t_displacement_at_gauss_point(j);
    // move to next base function
    ++t_base_func;
    // move the pointer to next element of t_nf
    ++t_nf;
    }
    // move to next integration weight
    ++t_w;
    // move to the solutions at the next Gauss point
    ++t_solution_at_gauss_point;
    ++t_init_solution_at_gauss_point;
    }
    // add computed values of spring in the global right hand side vector
    Vec f = getFEMethod()->ksp_f != PETSC_NULL ? getFEMethod()->ksp_f
    : getFEMethod()->snes_f;
    CHKERR VecSetValues(f, nb_dofs, &data.getIndices()[0], &nF[0], ADD_VALUES);
    }

Operator OpSpringKs for the LHS

Similarly to the previous part, the operator OpSpringKs which calculates the contribution of the spring elements to the LHS is implemented as follows

boost::shared_ptr<DataAtIntegrationPtsSprings> commonDataPtr;
BlockOptionDataSprings &dAta;
MatrixDouble locKs;
MatrixDouble transLocKs;
OpSpringKs(boost::shared_ptr<DataAtIntegrationPtsSprings> &common_data_ptr,
BlockOptionDataSprings &data, const std::string field_name)
field_name.c_str(), field_name.c_str(), OPROWCOL),
commonDataPtr(common_data_ptr), dAta(data) {}
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
EntityType col_type,
// check if the volumes have associated degrees of freedom
const int row_nb_dofs = row_data.getIndices().size();
if (!row_nb_dofs)
const int col_nb_dofs = col_data.getIndices().size();
if (!col_nb_dofs)
if (dAta.tRis.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
dAta.tRis.end()) {
}
CHKERR commonDataPtr->getBlockData(dAta);
// size associated to the entity
locKs.resize(row_nb_dofs, col_nb_dofs, false);
locKs.clear();
// get number of Gauss points
const int row_nb_gauss_pts = row_data.getN().size1();
if (!row_nb_gauss_pts) // check if number of Gauss point <> 0
const int row_nb_base_functions = row_data.getN().size2();
// get intergration weights
auto get_tensor2 = [](MatrixDouble &m, const int r, const int c) {
&m(3 * r + 0, 3 * c + 0), &m(3 * r + 0, 3 * c + 1),
&m(3 * r + 0, 3 * c + 2), &m(3 * r + 1, 3 * c + 0),
&m(3 * r + 1, 3 * c + 1), &m(3 * r + 1, 3 * c + 2),
&m(3 * r + 2, 3 * c + 0), &m(3 * r + 2, 3 * c + 1),
&m(3 * r + 2, 3 * c + 2));
};
commonDataPtr->springStiffnessNormal, 0., 0., 0.,
commonDataPtr->springStiffnessTangent, 0., 0., 0.,
commonDataPtr->springStiffnessTangent);
// create a 3d vector to be used as the normal to the face with length equal
// to the face area
auto t_normal_ptr = getFTensor1Normal();
t_normal(i) = t_normal_ptr(i);
// First tangent vector
auto t_tangent1_ptr = getFTensor1Tangent1AtGaussPts();
t_tangent1(i) = t_tangent1_ptr(i);
// Second tangent vector, such that t_n = t_t1 x t_t2 | t_t2 = t_n x t_t1
t_tangent2(i) = FTensor::levi_civita(i, j, k) * t_normal(j) * t_tangent1(k);
// Spring stiffness in global coordinate
t_spring_global = MetaSpringBC::transformLocalToGlobal(
t_normal, t_tangent1, t_tangent2, t_spring_local);
// loop over the Gauss points
for (int gg = 0; gg != row_nb_gauss_pts; gg++) {
// get area and integration weight
double w = t_w * getArea();
auto t_row_base_func = row_data.getFTensor0N(gg, 0);
for (int rr = 0; rr != row_nb_dofs / 3; rr++) {
auto t_col_base_func = col_data.getFTensor0N(gg, 0);
for (int cc = 0; cc != col_nb_dofs / 3; cc++) {
auto assemble_m = get_tensor2(locKs, rr, cc);
assemble_m(i, j) +=
w * t_row_base_func * t_col_base_func * t_spring_global(i, j);
++t_col_base_func;
}
++t_row_base_func;
}
// move to next integration weight
++t_w;
}
// Add computed values of spring stiffness to the global LHS matrix
Mat B = getFEMethod()->ksp_B != PETSC_NULL ? getFEMethod()->ksp_B
CHKERR MatSetValues(B, row_nb_dofs, &*row_data.getIndices().begin(),
col_nb_dofs, &*col_data.getIndices().begin(),
&locKs(0, 0), ADD_VALUES);
// is symmetric
if (row_side != col_side || row_type != col_type) {
transLocKs.resize(col_nb_dofs, row_nb_dofs, false);
noalias(transLocKs) = trans(locKs);
CHKERR MatSetValues(B, col_nb_dofs, &*col_data.getIndices().begin(),
row_nb_dofs, &*row_data.getIndices().begin(),
&transLocKs(0, 0), ADD_VALUES);
}
}
};

Other functions

Declaration of spring elements

This function declares the spring element allowing the implementation of the spring element in the main program of a module become more straightforward.

const std::string field_name,
const std::string mesh_nodals_positions) {
// Define boundary element that operates on rows, columns and data of a
// given field
CHKERR m_field.add_finite_element("SPRING");
if (m_field.check_field(mesh_nodals_positions)) {
mesh_nodals_positions);
}
// Add entities to that element, here we add all triangles with SPRING_BC
// from cubit
if (bit->getName().compare(0, 9, "SPRING_BC") == 0) {
MBTRI, "SPRING");
}
}
}

Pushing operators to finite element instances

With same motivation as the previous function, this function which pushes implemented operators (OpSpringKs, OpSpringFs) to the finite element instances makes the implementation in different module easier.

MoFEM::Interface &m_field,
boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_lhs_ptr,
boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_rhs_ptr,
const std::string field_name, const std::string mesh_nodals_positions) {
// Push operators to instances for springs
// loop over blocks
boost::shared_ptr<DataAtIntegrationPtsSprings> commonDataPtr =
boost::make_shared<DataAtIntegrationPtsSprings>(m_field);
CHKERR commonDataPtr->getParameters();
for (auto &sitSpring : commonDataPtr->mapSpring) {
fe_spring_lhs_ptr->getOpPtrVector().push_back(
new OpSpringKs(commonDataPtr, sitSpring.second, field_name));
fe_spring_rhs_ptr->getOpPtrVector().push_back(
new OpCalculateVectorFieldValues<3>(field_name, commonDataPtr->xAtPts));
fe_spring_rhs_ptr->getOpPtrVector().push_back(
new OpCalculateVectorFieldValues<3>(mesh_nodals_positions,
commonDataPtr->xInitAtPts));
fe_spring_rhs_ptr->getOpPtrVector().push_back(
new OpSpringFs(commonDataPtr, sitSpring.second, field_name));
}
}

Transformation of the spring

stiffnesses between local and global coordinates

This functions allows the transformation of the spring stiffnesses from local coordinate associated with the element face, where normal ( \(k_n\)) and in-plane (tangent, \(k_t\)) stiffnesses are initially defined, to the global coordinate.

FTensor::Tensor2<double, 3, 3> MetaSpringBC::transformLocalToGlobal(
FTensor::Tensor1<double, 3> t_tangent1_local,
FTensor::Tensor1<double, 3> t_tangent2_local,
// FTensor indices
// Global base vectors
FTensor::Tensor1<double, 3> t_e1(1., 0., 0.);
FTensor::Tensor1<double, 3> t_e2(0., 1., 0.);
FTensor::Tensor1<double, 3> t_e3(0., 0., 1.);
// Direction cosines
auto get_cosine = [&](auto x, auto xp) {
return (x(i) * xp(i)) / (sqrt(x(i) * x(i)) * sqrt(xp(i) * xp(i)));
};
// Transformation matrix (tensor)
FTensor::Tensor2<double, 3, 3> t_transformation_matrix(
get_cosine(t_e1, t_normal_local), get_cosine(t_e1, t_tangent1_local),
get_cosine(t_e1, t_tangent2_local), get_cosine(t_e2, t_normal_local),
get_cosine(t_e2, t_tangent1_local), get_cosine(t_e2, t_tangent2_local),
get_cosine(t_e3, t_normal_local), get_cosine(t_e3, t_tangent1_local),
get_cosine(t_e3, t_tangent2_local));
// Spring stiffness in global coordinate, Q*ls*Q^T
t_spring_global(i, j) = t_transformation_matrix(i, k) * t_spring_local(k, l) *
t_transformation_matrix(j, l);
return t_spring_global;
};

Sample applications

With the spring element developed, the boundary condition employing elastic spring to model the interaction between the structure of interest and surrounding environment can be applied to any existing solid problems. For demonstration purposes, the following part presents how spring elements are added to the linear and nonlinear elasticity problems, i.e. elasticity.cpp and arc_length_nonlinear_elasticity.cpp, and how the inclusion of spring boundary condition affects the responses of the structure bearing external loads.

For numerical examples, a rod with rectangular cross section will be used. The input for the geometry and mesh of a rod with spring boundary condition as follows

  • Simple journal script which can be used in Cubit to generate the rod geometry and mesh
    reset
    brick x 1.2 y 1.4 z 4
    #rotate Volume 1 angle 70 about X include_merged
    block 1 volume all
    block 1 name 'MAT_ELASTIC'
    block 1 attribute count 2
    block 1 attribute index 1 {young_modulus}
    block 1 attribute index 2 {poisson_ratio}
    # Boundary conditions
    #create displacement on vertex 6 dof 1 dof 2 fix 0
    #create displacement on vertex 8 dof 1 fix 0
    create pressure on surface 1 magnitude 12
    # Spring
    {spring_stiffness_normal = 5}
    {spring_stiffness_tangent = 100}
    block 4 surface 2
    block 4 name 'SPRING_BC'
    block 4 attribute count 2
    block 4 attribute index 1 {spring_stiffness_normal}
    block 4 attribute index 2 {spring_stiffness_tangent}
    volume all scheme Tetmesh
    volume all size auto factor 10
    mesh volume all
    save as "/Users/username/mofem_install/mofem-cephas/mofem/users_modules/basic_finite_elements/nonlinear_elasticity/spring_rod.cub" overwrite
  • The mesh of the rod can be found in Figure 2 where the translation in x and y directions of two vertices at the bottom surface are fixed to avoid rigid body rotation and the spring is applied at the bottom surface.

Figure 2: Mesh of a rod with pressure on the top surface and the boundary condition imposed on the bottom surface.

Elasticity problem

In order to include spring element in the linear elesticity problem elasticity.cpp, the following lines of code are needed in the main program.

  • The header files for spring element implementation are embedded in BasicFiniteElements.hpp (directly included) and All.cpp (via CMake)
  • Declaration of the spring element
    // Add spring boundary condition applied on surfaces.
    // This is only declaration not implementation.
    CHKERR MetaSpringBC::addSpringElements(m_field, "DISPLACEMENT",
    "MESH_NODE_POSITIONS");
  • Implementation of spring element by creating new instances and calling function MetaSpringBC::setSpringOperators to push operators to them
    // Implementation of spring element
    // Create new instances of face elements for springs
    boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_lhs_ptr(
    boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_rhs_ptr(
    CHKERR MetaSpringBC::setSpringOperators(m_field, fe_spring_lhs_ptr,
    fe_spring_rhs_ptr, "DISPLACEMENT",
    "MESH_NODE_POSITIONS");
  • Assign global matrix and vector to which springs contribute
    fe_spring_lhs_ptr->ksp_B = Aij;
    fe_spring_rhs_ptr->ksp_f = F;
  • Loop over the elements to calculate the contribution of the spring boundary condition to the RHS and LHS and assemble them
    // Assemble springs
    CHKERR DMoFEMLoopFiniteElements(dm, "SPRING", fe_spring_lhs_ptr);
    CHKERR DMoFEMLoopFiniteElements(dm, "SPRING", fe_spring_rhs_ptr);

Once the additional lines of code are added to include spring element to the linear elasticity module elasticity.cpp, the test using above rod model with springs is conducted. The displacement results using finite element analysis is presented in Figure 3. The command line that has been used to run this test is as follows

./elasticity -my_file spring_rod.cub -my_order 2 -ksp_type gmres -pc_type asm -sub_pc_type lu

Figure 3: Linear displacement result of the rod.

The hand calculation for the displacement of a point on the bottom surface where the spring is applied can be done as follows

\[ u = \frac{{F\left( {{A_F}} \right)}}{{k_s\left( {{A_k}} \right)}} = \frac{{p \times {A_F}}}{{{k_sz} \times {A_k}}} = \frac{{12 \times \left( {1.2 \times 1.4} \right)}}{{5 \times \left( {1.2 \times 1.4} \right)}} = 2.4 \]

where \( p \) and \( k_z \) are the pressure and spring stiffness in z direction, respectively. Meanwhile, \( A_F \) and \( A_k \) represent the areas where pressure and spring are applied, respectively. As expected, the hand calculated result is consistent with that shown in Figure 3.

Nonlinear elasticity problem

Similar to the previous application, in order to include spring element in the nonlinear elesticity problem arc_length_nonlinear_elasticity.cpp, the following lines of code are needed in the main program.

  • The header files for spring element implementation are embedded in BasicFiniteElements.hpp (directly included) and All.cpp (via CMake)
  • Declaration of the spring element
    // Add spring boundary condition applied on surfaces.
    // This is only declaration not implementation.
    CHKERR MetaSpringBC::addSpringElements(m_field, "SPATIAL_POSITION",
    "MESH_NODE_POSITIONS");
  • Set finite element for spring boundary condition
    CHKERR m_field.modify_problem_add_finite_element("ELASTIC_MECHANICS",
    "SPRING");
  • Implementation of spring element by creating new instances and calling function MetaSpringBC::setSpringOperators to push operators to them
    // Implementation of spring element
    // Create new instances of face elements for springs
    boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_lhs_ptr(
    boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_rhs_ptr(
    m_field, fe_spring_lhs_ptr, fe_spring_rhs_ptr, "SPATIAL_POSITION",
    "MESH_NODE_POSITIONS");
  • Loop over the elements to calculate the contribution of the spring boundary condition to the RHS and LHS
    loops_to_do_Rhs.push_back(
    SnesCtx::PairNameFEMethodPtr("SPRING", fe_spring_rhs_ptr.get()));
    loops_to_do_Mat.push_back(
    SnesCtx::PairNameFEMethodPtr("SPRING", fe_spring_lhs_ptr.get()));

The nonlinear analysis of the rod is run using the command line with the output of the first two load steps as follows

./arc_length_nonlinear_elasticity -my_file spring_rod.cub -my_sr 1e-6 -my_ms 5 -snes_converged_reason -my_order 2
dlambda = 9.1847e-08 dx2 = 0.0000e+00
res_lambda = 0.0000e+00
lambda = 9.1847e-08
0 SNES Function norm 1.000000000004e-06
diag = 2.1775e-05
Flambda2 = 1.1854e+02
lambda = 9.1847e-08
fnorm = 2.3303e-13
dlambda = 9.1847e-08 dx2 = 4.1486e-12
res_lambda = 0.0000e+00
lambda = 9.1847e-08
1 SNES Function norm 2.330272528283e-13
Nonlinear solve converged due to CONVERGED_FNORM_ABS iterations 1
number of Newton iterations = 1
Set alpha = 1.0000e+00 beta = 0.0000e+00
dlambda = 9.1847e-08 dx2 = 4.1486e-12
Set s = 2.0368e-06
Load Step 2 step_size = 2.0368e-06 dlambda0 = 9.1847e-08 dx_nrm = 2.0368e-06 dx2 = 4.1486e-12
lambda = 9.1847e-08, 1.8369e-07 (9.1847e-08)
Flambda2 = 1.1854e+02
lambda = 1.8369e-07
fnorm = 9.3481e-13
dlambda = 9.1847e-08 dx2 = 4.1486e-12
res_lambda = 0.0000e+00
lambda = 1.8369e-07
0 SNES Function norm 9.348064473394e-13
Nonlinear solve converged due to CONVERGED_FNORM_ABS iterations 0
number of Newton iterations = 0

Figure 4 presents the nonlinear displacement illustration of the rod after the first load step.

Figure 4: Nonlinear displacement result of the rod.

The hand calculation for the nonlinear displacement after the first load step of a point on the bottom surface where the spring is applied can be done as follows

\[ u = \frac{{F\left( {{A_F}} \right)}} {{k_s\left( {{A_k}} \right)}} = \frac{{\lambda \times p \times {A_F}}}{{{k_sz} \times {A_k}}} = \frac{{9.1847 \times 10^{-08} \times 12 \times \left({1.2 \times 1.4} \right)}}{{5 \times \left( {1.2 \times 1.4} \right)}} = 2.204\times 10^{-07} \]

where \( \lambda \) is the load factor being used in each load step which can be found in the output as \( lambda \). As expected, the numerical and hand-calculated results are consistent.

As an additional demonstration of the effect of spring boundary conditions on the structural response, the analysis result for a L-shaped structure is conducted. Figure 5 shows the geometry and the mesh while Figure 6 presents the comparison of the results for different boundary conditions applied at a face of the L-shaped structure. There are two analyses where springs with different stiffness magnitude in normal and in-plane local direction are applied. Additionally, another case in which clamp boundary condition is used is also presented. It is worth commenting that the deformed shapes of the structure get bigger are due to the small strain assumption and the solution scaling.

Figure 5: Mesh of the L-shaped structure with the boundary condition imposed on the highlighted surface.

Figure 6: Effects of different boundary conditions on the displacement (E = 10, v = 0, pressure = 25).
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::MatSetValues
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
Definition: EntitiesFieldData.hpp:1644
MoFEM::KspMethod::ksp_B
Mat & ksp_B
Definition: LoopMethods.hpp:91
surface
Definition: surface.py:1
FTensor::Tensor1< double, 3 >
young_modulus
double young_modulus
Young modulus.
Definition: plastic.cpp:121
MoFEM::CoreInterface::modify_finite_element_add_field_row
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
paraview_deform_macro.Function
Function
Definition: paraview_deform_macro.py:15
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::Types::MatrixDouble
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
MoFEM::SnesMethod::snes_B
Mat & snes_B
preconditioner of jacobian matrix
Definition: LoopMethods.hpp:124
FTensor::levi_civita
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
Definition: Levi_Civita.hpp:617
MoFEM::DataOperator::doWork
virtual 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.
Definition: DataOperators.hpp:40
MoFEM::VecSetValues
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
Definition: EntitiesFieldData.hpp:1589
sdf.r
int r
Definition: sdf.py:8
FTensor::Tensor2< double, 3, 3 >
MoFEM::ForcesAndSourcesCore::UserDataOperator::getFTensor0IntegrationWeight
auto getFTensor0IntegrationWeight()
Get integration weights.
Definition: ForcesAndSourcesCore.hpp:1240
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
c
const double c
speed of light (cm/ns)
Definition: initial_diffusion.cpp:39
MoFEM::CoreInterface::add_ents_to_finite_element_by_type
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
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::CoreInterface::add_finite_element
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::CoreInterface::modify_finite_element_add_field_col
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
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Normal
auto getFTensor1Normal()
get normal as tensor
Definition: FaceElementForcesAndSourcesCore.hpp:255
MoFEM::ForcesAndSourcesCore::UserDataOperator
Definition: ForcesAndSourcesCore.hpp:549
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator
default operator for TRI element
Definition: FaceElementForcesAndSourcesCore.hpp:94
bit
auto bit
set bit
Definition: hanging_node_approx.cpp:75
convert.type
type
Definition: convert.py:64
poisson_ratio
double poisson_ratio
Poisson ratio.
Definition: plastic.cpp:122
MoFEM::solve
static PetscErrorCode solve(Mat mat, Vec x, Vec y)
Definition: Schur.cpp:1223
MoFEM::CoreInterface::check_field
virtual bool check_field(const std::string &name) const =0
check if field is in database
MoFEM::CoreInterface::modify_finite_element_add_field_data
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_field)=0
set finite element field data
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
FaceElementForcesAndSourcesCore
EntData
EntitiesFieldData::EntData EntData
Definition: child_and_parent.cpp:37
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
FTensor::Index< 'i', 3 >
MoFEM::ForcesAndSourcesCore::UserDataOperator::getFEMethod
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
Definition: ForcesAndSourcesCore.hpp:1042
MetaSpringBC::setSpringOperators
static MoFEMErrorCode setSpringOperators(MoFEM::Interface &m_field, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_lhs_ptr, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_rhs_ptr, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS", double stiffness_scale=1.)
Implementation of spring element. Set operators to calculate LHS and RHS.
Definition: SpringElement.cpp:1178
_IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_
#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.
Definition: MeshsetsManager.hpp:71
HenckyOps::f
auto f
Definition: HenckyOps.hpp:15
BLOCKSET
@ BLOCKSET
Definition: definitions.h:161
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1Tangent1AtGaussPts
auto getFTensor1Tangent1AtGaussPts()
get tangent 1 at integration points
Definition: FaceElementForcesAndSourcesCore.hpp:323
lambda
static double lambda
Definition: incompressible_elasticity.cpp:199
EigenMatrix::Vec
const FTensor::Tensor2< T, Dim, Dim > Vec
Definition: MatrixFunction.hpp:66
MoFEM::Types::VectorDouble
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
sdf_wavy_2d.w
int w
Definition: sdf_wavy_2d.py:6
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
MoFEM::CoreInterface::modify_problem_add_finite_element
virtual MoFEMErrorCode modify_problem_add_finite_element(const std::string name_problem, const std::string &fe_name)=0
add finite element to problem, this add entities assigned to finite element to a particular problem
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator::getArea
double getArea()
get area of face
Definition: FaceElementForcesAndSourcesCore.hpp:239
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
MoFEM::DMoFEMLoopFiniteElements
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:586
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MetaSpringBC::addSpringElements
static MoFEMErrorCode addSpringElements(MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
Declare spring element.
Definition: SpringElement.cpp:1127
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
F
@ F
Definition: free_surface.cpp:394