v0.9.0
Public Member Functions | Public Attributes | List of all members
OpSpringKs Struct Reference

Assemble contribution of spring to LHS *

\[ {K^s} = \int\limits_\Omega ^{} {{\psi ^T}{k_s}\psi d\Omega } \]

. More...

Inheritance diagram for OpSpringKs:
[legend]
Collaboration diagram for OpSpringKs:
[legend]

Public Member Functions

 OpSpringKs (boost::shared_ptr< DataAtIntegrationPtsSprings > &common_data_ptr, BlockOptionDataSprings &data, const std::string field_name)
 
MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 

Public Attributes

boost::shared_ptr< DataAtIntegrationPtsSpringscommonDataPtr
 
BlockOptionDataSpringsdAta
 
MatrixDouble locKs
 
MatrixDouble transLocKs
 

Detailed Description

Assemble contribution of spring to LHS *

\[ {K^s} = \int\limits_\Omega ^{} {{\psi ^T}{k_s}\psi d\Omega } \]

.

Definition at line 246 of file SpringElement.cpp.

Constructor & Destructor Documentation

◆ OpSpringKs()

OpSpringKs::OpSpringKs ( boost::shared_ptr< DataAtIntegrationPtsSprings > &  common_data_ptr,
BlockOptionDataSprings data,
const std::string  field_name 
)

Definition at line 254 of file SpringElement.cpp.

257  field_name.c_str(), field_name.c_str(), OPROWCOL),
258  commonDataPtr(common_data_ptr), dAta(data) {}
boost::shared_ptr< DataAtIntegrationPtsSprings > commonDataPtr
BlockOptionDataSprings & dAta
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

MoFEMErrorCode OpSpringKs::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
DataForcesAndSourcesCore::EntData row_data,
DataForcesAndSourcesCore::EntData col_data 
)

Definition at line 260 of file SpringElement.cpp.

263  {
265 
266  // check if the volumes have associated degrees of freedom
267  const int row_nb_dofs = row_data.getIndices().size();
268  if (!row_nb_dofs)
270 
271  const int col_nb_dofs = col_data.getIndices().size();
272  if (!col_nb_dofs)
274 
275  if (dAta.tRis.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
276  dAta.tRis.end()) {
278  }
279 
280  CHKERR commonDataPtr->getBlockData(dAta);
281  // size associated to the entity
282  locKs.resize(row_nb_dofs, col_nb_dofs, false);
283  locKs.clear();
284 
285  // get number of Gauss points
286  const int row_nb_gauss_pts = row_data.getN().size1();
287  if (!row_nb_gauss_pts) // check if number of Gauss point <> 0
289 
290  // get intergration weights
291  auto t_w = getFTensor0IntegrationWeight();
292 
296  auto get_tensor2 = [](MatrixDouble &m, const int r, const int c) {
298  &m(3 * r + 0, 3 * c + 0), &m(3 * r + 0, 3 * c + 1),
299  &m(3 * r + 0, 3 * c + 2), &m(3 * r + 1, 3 * c + 0),
300  &m(3 * r + 1, 3 * c + 1), &m(3 * r + 1, 3 * c + 2),
301  &m(3 * r + 2, 3 * c + 0), &m(3 * r + 2, 3 * c + 1),
302  &m(3 * r + 2, 3 * c + 2));
303  };
304 
305  FTensor::Tensor2<double, 3, 3> t_spring_local(
306  commonDataPtr->springStiffnessNormal, 0., 0., 0.,
307  commonDataPtr->springStiffnessTangent, 0., 0., 0.,
308  commonDataPtr->springStiffnessTangent);
309  // create a 3d vector to be used as the normal to the face with length equal
310  // to the face area
311  auto t_normal_ptr = getFTensor1Normal();
312 
314  t_normal(i) = t_normal_ptr(i);
315 
316  // First tangent vector
317  auto t_tangent1_ptr = getFTensor1Tangent1AtGaussPts();
318  FTensor::Tensor1<double, 3> t_tangent1;
319  t_tangent1(i) = t_tangent1_ptr(i);
320 
321  // Second tangent vector, such that t_n = t_t1 x t_t2 | t_t2 = t_n x t_t1
322  FTensor::Tensor1<double, 3> t_tangent2;
323  t_tangent2(i) = FTensor::levi_civita(i, j, k) * t_normal(j) * t_tangent1(k);
324 
325  // Spring stiffness in global coordinate
326  FTensor::Tensor2<double, 3, 3> t_spring_global;
327  t_spring_global = MetaSpringBC::transformLocalToGlobal(
328  t_normal, t_tangent1, t_tangent2, t_spring_local);
329 
330  // loop over the Gauss points
331  for (int gg = 0; gg != row_nb_gauss_pts; gg++) {
332  // get area and integration weight
333  double w = t_w * getArea();
334 
335  auto t_row_base_func = row_data.getFTensor0N(gg, 0);
336 
337  for (int rr = 0; rr != row_nb_dofs / 3; rr++) {
338  auto t_col_base_func = col_data.getFTensor0N(gg, 0);
339  for (int cc = 0; cc != col_nb_dofs / 3; cc++) {
340  auto assemble_m = get_tensor2(locKs, rr, cc);
341  assemble_m(i, j) +=
342  w * t_row_base_func * t_col_base_func * t_spring_global(i, j);
343  ++t_col_base_func;
344  }
345  ++t_row_base_func;
346  }
347  // move to next integration weight
348  ++t_w;
349  }
350 
351  // Add computed values of spring stiffness to the global LHS matrix
352  Mat B = getFEMethod()->ksp_B != PETSC_NULL ? getFEMethod()->ksp_B
353  : getFEMethod()->snes_B;
354  CHKERR MatSetValues(B, row_nb_dofs, &*row_data.getIndices().begin(),
355  col_nb_dofs, &*col_data.getIndices().begin(),
356  &locKs(0, 0), ADD_VALUES);
357 
358  // is symmetric
359  if (row_side != col_side || row_type != col_type) {
360  transLocKs.resize(col_nb_dofs, row_nb_dofs, false);
361  noalias(transLocKs) = trans(locKs);
362 
363  CHKERR MatSetValues(B, col_nb_dofs, &*col_data.getIndices().begin(),
364  row_nb_dofs, &*row_data.getIndices().begin(),
365  &transLocKs(0, 0), ADD_VALUES);
366  }
367 
369  }
MatrixDouble transLocKs
MatrixDouble locKs
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
boost::shared_ptr< DataAtIntegrationPtsSprings > commonDataPtr
const VectorInt & getIndices() const
Get global indices of dofs on entity.
MoFEMErrorCode MatSetValues(Mat M, const DataForcesAndSourcesCore::EntData &row_data, const DataForcesAndSourcesCore::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
BlockOptionDataSprings & dAta
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
static FTensor::Tensor2< double, 3, 3 > transformLocalToGlobal(FTensor::Tensor1< double, 3 > t_normal, FTensor::Tensor1< double, 3 > t_tangent1, FTensor::Tensor1< double, 3 > t_tangent2, FTensor::Tensor2< double, 3, 3 > t_spring_local)
Implementation of spring element. Set operators to calculate LHS and RHS.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
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

Member Data Documentation

◆ commonDataPtr

boost::shared_ptr<DataAtIntegrationPtsSprings> OpSpringKs::commonDataPtr

Definition at line 248 of file SpringElement.cpp.

◆ dAta

BlockOptionDataSprings& OpSpringKs::dAta

Definition at line 249 of file SpringElement.cpp.

◆ locKs

MatrixDouble OpSpringKs::locKs

Definition at line 251 of file SpringElement.cpp.

◆ transLocKs

MatrixDouble OpSpringKs::transLocKs

Definition at line 252 of file SpringElement.cpp.


The documentation for this struct was generated from the following file: