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

Assemble contribution of springs to RHS *

\[ f_s = \int\limits_{\partial \Omega }^{} {{\psi ^T}{F^s}\left( u \right)d\partial \Omega } = \int\limits_{\partial \Omega }^{} {{\psi ^T}{k_s}ud\partial \Omega } \]

. More...

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

Public Member Functions

 OpSpringFs (boost::shared_ptr< DataAtIntegrationPtsSprings > &common_data_ptr, BlockOptionDataSprings &data, const std::string field_name)
 
MoFEMErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
 

Public Attributes

VectorDouble nF
 
boost::shared_ptr< DataAtIntegrationPtsSpringscommonDataPtr
 
BlockOptionDataSpringsdAta
 
bool is_spatial_position = true
 

Detailed Description

Assemble contribution of springs to RHS *

\[ f_s = \int\limits_{\partial \Omega }^{} {{\psi ^T}{F^s}\left( u \right)d\partial \Omega } = \int\limits_{\partial \Omega }^{} {{\psi ^T}{k_s}ud\partial \Omega } \]

.

Definition at line 113 of file SpringElement.cpp.

Constructor & Destructor Documentation

◆ OpSpringFs()

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

Definition at line 122 of file SpringElement.cpp.

125  field_name.c_str(), OPROW),
126  commonDataPtr(common_data_ptr), dAta(data) {
127  if (field_name.compare(0, 16, "SPATIAL_POSITION") != 0)
128  is_spatial_position = false;
129  }
boost::shared_ptr< DataAtIntegrationPtsSprings > commonDataPtr
bool is_spatial_position
BlockOptionDataSprings & dAta
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

MoFEMErrorCode OpSpringFs::doWork ( int  side,
EntityType  type,
DataForcesAndSourcesCore::EntData data 
)

Definition at line 131 of file SpringElement.cpp.

132  {
133 
135 
136  // check that the faces have associated degrees of freedom
137  const int nb_dofs = data.getIndices().size();
138  if (nb_dofs == 0)
140 
141  if (dAta.tRis.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
142  dAta.tRis.end()) {
144  }
145 
146  CHKERR commonDataPtr->getBlockData(dAta);
147 
148  // size of force vector associated to the entity
149  // set equal to the number of degrees of freedom of associated with the
150  // entity
151  nF.resize(nb_dofs, false);
152  nF.clear();
153 
154  // get number of Gauss points
155  const int nb_gauss_pts = data.getN().size1();
156 
157  // get intergration weights
158  auto t_w = getFTensor0IntegrationWeight();
159 
160  // FTensor indices
164 
165  FTensor::Tensor2<double, 3, 3> t_spring_local(
166  commonDataPtr->springStiffnessNormal, 0., 0., 0.,
167  commonDataPtr->springStiffnessTangent, 0., 0., 0.,
168  commonDataPtr->springStiffnessTangent);
169  // create a 3d vector to be used as the normal to the face with length equal
170  // to the face area
171  auto t_normal_ptr = getFTensor1Normal();
172 
174  t_normal(i) = t_normal_ptr(i);
175 
176  // First tangent vector
177  auto t_tangent1_ptr = getFTensor1Tangent1AtGaussPts();
178  FTensor::Tensor1<double, 3> t_tangent1;
179  t_tangent1(i) = t_tangent1_ptr(i);
180 
181  // Second tangent vector, such that t_n = t_t1 x t_t2 | t_t2 = t_n x t_t1
182  FTensor::Tensor1<double, 3> t_tangent2;
183  t_tangent2(i) = FTensor::levi_civita(i, j, k) * t_normal(j) * t_tangent1(k);
184 
185  // Spring stiffness in global coordinate
186  FTensor::Tensor2<double, 3, 3> t_spring_global;
187  t_spring_global = MetaSpringBC::transformLocalToGlobal(
188  t_normal, t_tangent1, t_tangent2, t_spring_local);
189 
190  // Extract solution at Gauss points
191  auto t_solution_at_gauss_point =
192  getFTensor1FromMat<3>(*commonDataPtr->xAtPts);
193  auto t_init_solution_at_gauss_point =
194  getFTensor1FromMat<3>(*commonDataPtr->xInitAtPts);
195  FTensor::Tensor1<double, 3> t_displacement_at_gauss_point;
196 
197  // loop over all Gauss points of the face
198  for (int gg = 0; gg != nb_gauss_pts; ++gg) {
199 
200  // Calculate the displacement at the Gauss point
201  if (is_spatial_position) { // "SPATIAL_POSITION"
202  t_displacement_at_gauss_point(i) =
203  t_solution_at_gauss_point(i) - t_init_solution_at_gauss_point(i);
204  } else { // e.g. "DISPLACEMENT" or "U"
205  t_displacement_at_gauss_point(i) = t_solution_at_gauss_point(i);
206  }
207 
208  double w = t_w * getArea();
209 
210  auto t_base_func = data.getFTensor0N(gg, 0);
211 
212  // create a vector t_nf whose pointer points an array of 3 pointers
213  // pointing to nF memory location of components
215  &nF[2]);
216 
217  for (int rr = 0; rr != nb_dofs / 3; ++rr) { // loop over the nodes
218  t_nf(i) += w * t_base_func * t_spring_global(i, j) *
219  t_displacement_at_gauss_point(j);
220 
221  // move to next base function
222  ++t_base_func;
223  // move the pointer to next element of t_nf
224  ++t_nf;
225  }
226  // move to next integration weight
227  ++t_w;
228  // move to the solutions at the next Gauss point
229  ++t_solution_at_gauss_point;
230  ++t_init_solution_at_gauss_point;
231  }
232  // add computed values of spring in the global right hand side vector
233  Vec f = getFEMethod()->ksp_f != PETSC_NULL ? getFEMethod()->ksp_f
234  : getFEMethod()->snes_f;
235  CHKERR VecSetValues(f, nb_dofs, &data.getIndices()[0], &nF[0], ADD_VALUES);
237  }
boost::shared_ptr< DataAtIntegrationPtsSprings > commonDataPtr
#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
const VectorInt & getIndices() const
Get global indices of dofs on entity.
bool is_spatial_position
MoFEMErrorCode VecSetValues(Vec V, const DataForcesAndSourcesCore::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
BlockOptionDataSprings & dAta
#define CHKERR
Inline error check.
Definition: definitions.h:596
VectorDouble nF
#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> OpSpringFs::commonDataPtr

Definition at line 118 of file SpringElement.cpp.

◆ dAta

BlockOptionDataSprings& OpSpringFs::dAta

Definition at line 119 of file SpringElement.cpp.

◆ is_spatial_position

bool OpSpringFs::is_spatial_position = true

Definition at line 120 of file SpringElement.cpp.

◆ nF

VectorDouble OpSpringFs::nF

Definition at line 116 of file SpringElement.cpp.


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