v0.15.0
Loading...
Searching...
No Matches
OpAssembleK Struct Reference

Assemble K *. More...

#include "tutorials/cor-7/src/ElasticityMixedFormulation.hpp"

Inheritance diagram for OpAssembleK:
[legend]
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, EntitiesFieldData::EntData &row_data, EntitiesFieldData::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
ElasticityMixedFormulation.hpp, and elasticity_mixed_formulation.cpp.

Definition at line 310 of file ElasticityMixedFormulation.hpp.

Constructor & Destructor Documentation

◆ OpAssembleK()

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

Definition at line 320 of file ElasticityMixedFormulation.hpp.

322 : VolumeElementForcesAndSourcesCore::UserDataOperator("U", "U", OPROWCOL,
323 symm),
324 commonData(common_data), dAta(data) {}
DataAtIntegrationPts & commonData

Member Function Documentation

◆ doWork()

PetscErrorCode OpAssembleK::doWork ( int row_side,
int col_side,
EntityType row_type,
EntityType col_type,
EntitiesFieldData::EntData & row_data,
EntitiesFieldData::EntData & col_data )
inline
Examples
ElasticityMixedFormulation.hpp.

Definition at line 326 of file ElasticityMixedFormulation.hpp.

329 {
331
332 auto get_tensor2 = [](MatrixDouble &m, const int r, const int c) {
334 &m(r + 0, c + 0), &m(r + 0, c + 1), &m(r + 0, c + 2),
335 &m(r + 1, c + 0), &m(r + 1, c + 1), &m(r + 1, c + 2),
336 &m(r + 2, c + 0), &m(r + 2, c + 1), &m(r + 2, c + 2));
337 };
338
339 const int row_nb_dofs = row_data.getIndices().size();
340 if (!row_nb_dofs)
342 const int col_nb_dofs = col_data.getIndices().size();
343 if (!col_nb_dofs)
345 if (dAta.tEts.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
346 dAta.tEts.end()) {
348 }
350
351 const bool diagonal_block =
352 (row_type == col_type) && (row_side == col_side);
353 // get number of integration points
354 // Set size can clear local tangent matrix
355 locK.resize(row_nb_dofs, col_nb_dofs, false);
356 locK.clear();
357
358 const int row_nb_gauss_pts = row_data.getN().size1();
359 const int row_nb_base_functions = row_data.getN().size2();
360
361 auto row_diff_base_functions = row_data.getFTensor1DiffN<3>();
362
363 FTensor::Index<'i', 3> i;
364 FTensor::Index<'j', 3> j;
365 FTensor::Index<'k', 3> k;
366 FTensor::Index<'l', 3> l;
367
368 // integrate local matrix for entity block
369 for (int gg = 0; gg != row_nb_gauss_pts; gg++) {
370
371 // Get volume and integration weight
372 double w = getVolume() * getGaussPts()(3, gg);
373
374 int row_bb = 0;
375 for (; row_bb != row_nb_dofs / 3; row_bb++) {
376
377 auto col_diff_base_functions = col_data.getFTensor1DiffN<3>(gg, 0);
378 const int final_bb = diagonal_block ? row_bb + 1 : col_nb_dofs / 3;
379 int col_bb = 0;
380 for (; col_bb != final_bb; col_bb++) {
381
382 auto t_assemble = get_tensor2(locK, 3 * row_bb, 3 * col_bb);
383
384 diffDiff(j, l) =
385 w * row_diff_base_functions(j) * col_diff_base_functions(l);
386
387 t_assemble(i, k) += diffDiff(j, l) * commonData.tD(i, j, k, l);
388 // Next base function for column
389 ++col_diff_base_functions;
390 }
391
392 ++row_diff_base_functions;
393 }
394 for (; row_bb != row_nb_base_functions; row_bb++) {
395 ++row_diff_base_functions;
396 }
397 }
398
399 if (diagonal_block) {
400 for (int row_bb = 0; row_bb != row_nb_dofs / 3; row_bb++) {
401 int col_bb = 0;
402 for (; col_bb != row_bb + 1; col_bb++) {
403 auto t_assemble = get_tensor2(locK, 3 * row_bb, 3 * col_bb);
404 auto t_off_side = get_tensor2(locK, 3 * col_bb, 3 * row_bb);
405 t_off_side(i, k) = t_assemble(k, i);
406 }
407 }
408 }
409
410 const int *row_ind = &*row_data.getIndices().begin();
411 const int *col_ind = &*col_data.getIndices().begin();
412 Mat B = getFEMethod()->ksp_B != PETSC_NULLPTR ? getFEMethod()->ksp_B
413 : getFEMethod()->ksp_B;
414 CHKERR MatSetValues(B, row_nb_dofs, row_ind, col_nb_dofs, col_ind,
415 &locK(0, 0), ADD_VALUES);
416
417 if (row_type != col_type || row_side != col_side) {
418 translocK.resize(col_nb_dofs, row_nb_dofs, false);
419 noalias(translocK) = trans(locK);
420 CHKERR MatSetValues(B, col_nb_dofs, col_ind, row_nb_dofs, row_ind,
421 &translocK(0, 0), ADD_VALUES);
422 }
423
425 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FTensor::Index< 'i', SPACE_DIM > i
const double c
speed of light (cm/ns)
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
int r
Definition sdf.py:8
FTensor::Index< 'm', 3 > m
MoFEMErrorCode getBlockData(BlockData &data)
FTensor::Ddg< double, 3, 3 > tD
FTensor::Tensor2< double, 3, 3 > diffDiff

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: