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

Assemble G *

\[ {\bf{G}} = - \int\limits_\Omega {{{\bf{B}}^T}{\bf m}{{\bf{N}}_p}d\Omega } \]

. More...

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

Inherits UserDataOperator.

Collaboration diagram for OpAssembleG:
[legend]

Public Member Functions

 OpAssembleG (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 locG
 
BlockDatadAta
 

Detailed Description

Assemble G *

\[ {\bf{G}} = - \int\limits_\Omega {{{\bf{B}}^T}{\bf m}{{\bf{N}}_p}d\Omega } \]

.

Examples
elasticity_mixed_formulation.cpp, and ElasticityMixedFormulation.hpp.

Definition at line 231 of file ElasticityMixedFormulation.hpp.

Constructor & Destructor Documentation

◆ OpAssembleG()

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

Definition at line 238 of file ElasticityMixedFormulation.hpp.

240  "U", "P", UserDataOperator::OPROWCOL),
241  commonData(common_data), dAta(data) {
242  sYmm = false;
243  }
DataAtIntegrationPts & commonData
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

PetscErrorCode OpAssembleG::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 245 of file ElasticityMixedFormulation.hpp.

248  {
249 
251 
252  const int row_nb_dofs = row_data.getIndices().size();
253  if (!row_nb_dofs)
255  const int col_nb_dofs = col_data.getIndices().size();
256  if (!col_nb_dofs)
258 
259  if (dAta.tEts.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
260  dAta.tEts.end()) {
262  }
264 
265  // Set size can clear local tangent matrix
266  locG.resize(row_nb_dofs, col_nb_dofs, false);
267  locG.clear();
268  const int row_nb_gauss_pts = row_data.getN().size1();
269  if (!row_nb_gauss_pts)
271  const int row_nb_base_functions = row_data.getN().size2();
272  auto row_diff_base_functions = row_data.getFTensor1DiffN<3>();
273 
277 
278  // INTEGRATION
279  for (int gg = 0; gg != row_nb_gauss_pts; gg++) {
280 
281  // Get volume and integration weight
282  double w = getVolume() * getGaussPts()(3, gg);
283  if (getHoGaussPtsDetJac().size() > 0) {
284  w *= getHoGaussPtsDetJac()[gg]; ///< higher order geometry
285  }
286 
287  int row_bb = 0;
288  for (; row_bb != row_nb_dofs / 3; row_bb++) {
289 
290  t1(i) = w * row_diff_base_functions(i);
291 
292  auto base_functions = col_data.getFTensor0N(gg, 0);
293  for (int col_bb = 0; col_bb != col_nb_dofs; col_bb++) {
294 
295  FTensor::Tensor1<double *, 3> k(&locG(3 * row_bb + 0, col_bb),
296  &locG(3 * row_bb + 1, col_bb),
297  &locG(3 * row_bb + 2, col_bb));
298 
299  k(i) += t1(i) * base_functions;
300  ++base_functions;
301  }
302  ++row_diff_base_functions;
303  }
304  for (; row_bb != row_nb_base_functions; row_bb++) {
305  ++row_diff_base_functions;
306  }
307  }
308 
309  CHKERR MatSetValues(getFEMethod()->ksp_B, row_nb_dofs,
310  &*row_data.getIndices().begin(), col_nb_dofs,
311  &*col_data.getIndices().begin(), &*locG.data().begin(),
312  ADD_VALUES);
313 
314  // ASSEMBLE THE TRANSPOSE
315  locG = trans(locG);
316  CHKERR MatSetValues(getFEMethod()->ksp_B, col_nb_dofs,
317  &*col_data.getIndices().begin(), row_nb_dofs,
318  &*row_data.getIndices().begin(), &*locG.data().begin(),
319  ADD_VALUES);
321  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
DataAtIntegrationPts & commonData
static Index< 'i', 3 > i
static Index< 'j', 3 > j
#define CHKERR
Inline error check.
Definition: definitions.h:604
MoFEMErrorCode getBlockData(BlockData &data)
static Index< 'k', 3 > k
#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& OpAssembleG::commonData

◆ dAta

BlockData& OpAssembleG::dAta

◆ locG

MatrixDouble OpAssembleG::locG

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