v0.10.0
Public Member Functions | Public Attributes | List of all members
OpAssembleP Struct Reference

Assemble P *

\[ {\bf{P}} = - \int\limits_\Omega {{\bf{N}}_p^T\frac{1}{\lambda }{{\bf{N}}_p}d\Omega } \]

. More...

#include <users_modules/tutorials/low/elasticity_mixed_formulation/src/ElasticityMixedFormulation.hpp>

Inherits UserDataOperator.

Collaboration diagram for OpAssembleP:
[legend]

Public Member Functions

 OpAssembleP (DataAtIntegrationPts &common_data, BlockData &data)
 
PetscErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 

Public Attributes

DataAtIntegrationPtscommonData
 
MatrixDouble locP
 
MatrixDouble translocP
 
BlockDatadAta
 

Detailed Description

Assemble P *

\[ {\bf{P}} = - \int\limits_\Omega {{\bf{N}}_p^T\frac{1}{\lambda }{{\bf{N}}_p}d\Omega } \]

.

Examples
elasticity_mixed_formulation.cpp, and ElasticityMixedFormulation.hpp.

Definition at line 127 of file ElasticityMixedFormulation.hpp.

Constructor & Destructor Documentation

◆ OpAssembleP()

OpAssembleP::OpAssembleP ( DataAtIntegrationPts common_data,
BlockData data 
)
Examples
ElasticityMixedFormulation.hpp.

Definition at line 135 of file ElasticityMixedFormulation.hpp.

137  "P", "P", UserDataOperator::OPROWCOL),
138  commonData(common_data), dAta(data) {
139  sYmm = true;
140  }
DataAtIntegrationPts & commonData
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

PetscErrorCode OpAssembleP::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
DataForcesAndSourcesCore::EntData row_data,
DataForcesAndSourcesCore::EntData col_data 
)
Examples
ElasticityMixedFormulation.hpp.

Definition at line 142 of file ElasticityMixedFormulation.hpp.

145  {
146 
148  const int row_nb_dofs = row_data.getIndices().size();
149  if (!row_nb_dofs)
151  const int col_nb_dofs = col_data.getIndices().size();
152  if (!col_nb_dofs)
154 
155  if (dAta.tEts.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
156  dAta.tEts.end()) {
158  }
160  // Set size can clear local tangent matrix
161  locP.resize(row_nb_dofs, col_nb_dofs, false);
162  locP.clear();
163 
164  const int row_nb_gauss_pts = row_data.getN().size1();
165  if (!row_nb_gauss_pts)
167  const int row_nb_base_functions = row_data.getN().size2();
168  auto row_base_functions = row_data.getFTensor0N();
169 
170  // get data
173  const double lambda = commonData.lAmbda;
174 
175  double coefficient = commonData.pOisson == 0.5 ? 0. : 1 / lambda;
176 
177  // integration
178  if (coefficient != 0.) {
179  for (int gg = 0; gg != row_nb_gauss_pts; gg++) {
180 
181  // Get volume and integration weight
182  double w = getVolume() * getGaussPts()(3, gg);
183  if (getHoGaussPtsDetJac().size() > 0) {
184  w *= getHoGaussPtsDetJac()[gg]; ///< higher order geometry
185  }
186 
187  // INTEGRATION
188  int row_bb = 0;
189  for (; row_bb != row_nb_dofs; row_bb++) {
190  auto col_base_functions = col_data.getFTensor0N(gg, 0);
191  for (int col_bb = 0; col_bb != col_nb_dofs; col_bb++) {
192 
193  locP(row_bb, col_bb) -=
194  w * row_base_functions * col_base_functions * coefficient;
195 
196  ++col_base_functions;
197  }
198  ++row_base_functions;
199  }
200  for (; row_bb != row_nb_base_functions; row_bb++) {
201  ++row_base_functions;
202  }
203  }
204  }
205 
206 
208  getFEMethod()->ksp_B, row_nb_dofs, &*row_data.getIndices().begin(),
209  col_nb_dofs, &*col_data.getIndices().begin(), &locP(0, 0), ADD_VALUES);
210 
211  // is symmetric
212  if (row_side != col_side || row_type != col_type) {
213  translocP.resize(col_nb_dofs, row_nb_dofs, false);
214  noalias(translocP) = trans(locP);
215  CHKERR MatSetValues(getFEMethod()->ksp_B, col_nb_dofs,
216  &*col_data.getIndices().begin(), row_nb_dofs,
217  &*row_data.getIndices().begin(), &translocP(0, 0),
218  ADD_VALUES);
219  }
220 
222  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
DataAtIntegrationPts & commonData
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
static Index< 'i', 3 > i
static Index< 'j', 3 > j
#define CHKERR
Inline error check.
Definition: definitions.h:604
MoFEMErrorCode getBlockData(BlockData &data)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
double w(const double g, const double t)
Definition: ContactOps.hpp:169

Member Data Documentation

◆ commonData

DataAtIntegrationPts& OpAssembleP::commonData

◆ dAta

BlockData& OpAssembleP::dAta

◆ locP

MatrixDouble OpAssembleP::locP

◆ translocP

MatrixDouble OpAssembleP::translocP

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