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

Assemble K *

\[ {\bf{K}} = \int\limits_\Omega {{{\bf{B}}^T}{{\bf{D}}_d}{\bf{B}}d\Omega } \]

. More...

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

Inherits UserDataOperator.

Collaboration diagram for OpAssembleK:
[legend]

Public Member Functions

 OpAssembleK (DataAtIntegrationPts &common_data, BlockData &data, bool symm=true)
 
PetscErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 

Public Attributes

MatrixDouble locK
 
MatrixDouble translocK
 
FTensor::Tensor2< double, 3, 3 > diffDiff
 
DataAtIntegrationPtscommonData
 
BlockDatadAta
 

Detailed Description

Assemble K *

\[ {\bf{K}} = \int\limits_\Omega {{{\bf{B}}^T}{{\bf{D}}_d}{\bf{B}}d\Omega } \]

.

Examples
elasticity_mixed_formulation.cpp, and ElasticityMixedFormulation.hpp.

Definition at line 330 of file ElasticityMixedFormulation.hpp.

Constructor & Destructor Documentation

◆ OpAssembleK()

OpAssembleK::OpAssembleK ( DataAtIntegrationPts common_data,
BlockData data,
bool  symm = true 
)
Examples
ElasticityMixedFormulation.hpp.

Definition at line 340 of file ElasticityMixedFormulation.hpp.

343  symm),
344  commonData(common_data), dAta(data) {}
DataAtIntegrationPts & commonData
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

PetscErrorCode OpAssembleK::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 346 of file ElasticityMixedFormulation.hpp.

349  {
351 
352  auto get_tensor2 = [](MatrixDouble &m, const int r, const int c) {
354  &m(r + 0, c + 0), &m(r + 0, c + 1), &m(r + 0, c + 2),
355  &m(r + 1, c + 0), &m(r + 1, c + 1), &m(r + 1, c + 2),
356  &m(r + 2, c + 0), &m(r + 2, c + 1), &m(r + 2, c + 2));
357  };
358 
359  const int row_nb_dofs = row_data.getIndices().size();
360  if (!row_nb_dofs)
362  const int col_nb_dofs = col_data.getIndices().size();
363  if (!col_nb_dofs)
365  if (dAta.tEts.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
366  dAta.tEts.end()) {
368  }
370 
371  const bool diagonal_block =
372  (row_type == col_type) && (row_side == col_side);
373  // get number of integration points
374  // Set size can clear local tangent matrix
375  locK.resize(row_nb_dofs, col_nb_dofs, false);
376  locK.clear();
377 
378  const int row_nb_gauss_pts = row_data.getN().size1();
379  const int row_nb_base_functions = row_data.getN().size2();
380 
381  auto row_diff_base_functions = row_data.getFTensor1DiffN<3>();
382 
387 
388  // integrate local matrix for entity block
389  for (int gg = 0; gg != row_nb_gauss_pts; gg++) {
390 
391  // Get volume and integration weight
392  double w = getVolume() * getGaussPts()(3, gg);
393  if (getHoGaussPtsDetJac().size() > 0) {
394  w *= getHoGaussPtsDetJac()[gg]; ///< higher order geometry
395  }
396 
397  int row_bb = 0;
398  for (; row_bb != row_nb_dofs / 3; row_bb++) {
399 
400  auto col_diff_base_functions = col_data.getFTensor1DiffN<3>(gg, 0);
401  const int final_bb = diagonal_block ? row_bb + 1 : col_nb_dofs / 3;
402  int col_bb = 0;
403  for (; col_bb != final_bb; col_bb++) {
404 
405  auto t_assemble = get_tensor2(locK, 3 * row_bb, 3 * col_bb);
406 
407  diffDiff(j, l) =
408  w * row_diff_base_functions(j) * col_diff_base_functions(l);
409 
410  t_assemble(i, k) += diffDiff(j, l) * commonData.tD(i, j, k, l);
411  // Next base function for column
412  ++col_diff_base_functions;
413  }
414 
415  ++row_diff_base_functions;
416  }
417  for (; row_bb != row_nb_base_functions; row_bb++) {
418  ++row_diff_base_functions;
419  }
420  }
421 
422  if (diagonal_block) {
423  for (int row_bb = 0; row_bb != row_nb_dofs / 3; row_bb++) {
424  int col_bb = 0;
425  for (; col_bb != row_bb + 1; col_bb++) {
426  auto t_assemble = get_tensor2(locK, 3 * row_bb, 3 * col_bb);
427  auto t_off_side = get_tensor2(locK, 3 * col_bb, 3 * row_bb);
428  t_off_side(i, k) = t_assemble(k, i);
429  }
430  }
431  }
432 
433  const int *row_ind = &*row_data.getIndices().begin();
434  const int *col_ind = &*col_data.getIndices().begin();
435  Mat B = getFEMethod()->ksp_B != PETSC_NULL ? getFEMethod()->ksp_B
436  : getFEMethod()->ksp_B;
437  CHKERR MatSetValues(B, row_nb_dofs, row_ind, col_nb_dofs, col_ind,
438  &locK(0, 0), ADD_VALUES);
439 
440  if (row_type != col_type || row_side != col_side) {
441  translocK.resize(col_nb_dofs, row_nb_dofs, false);
442  noalias(translocK) = trans(locK);
443  CHKERR MatSetValues(B, col_nb_dofs, col_ind, row_nb_dofs, row_ind,
444  &translocK(0, 0), ADD_VALUES);
445  }
446 
448  }
DataAtIntegrationPts & commonData
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
FTensor::Ddg< double, 3, 3 > tD
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
FTensor::Tensor2< double, 3, 3 > diffDiff
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
static Index< 'l', 3 > l
static Index< 'm', 3 > m
static Index< 'i', 3 > i
static Index< 'j', 3 > j
#define CHKERR
Inline error check.
Definition: definitions.h:604
const double r
rate factor
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& OpAssembleK::commonData

◆ dAta

BlockData& OpAssembleK::dAta

◆ diffDiff

FTensor::Tensor2<double, 3, 3> OpAssembleK::diffDiff

◆ locK

MatrixDouble OpAssembleK::locK

◆ translocK

MatrixDouble OpAssembleK::translocK

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