v0.14.0
Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
EshelbianPlasticity::CGGUserPolynomialBase Struct Reference

#include <users_modules/eshelbian_plasticit/src/EshelbianPlasticity.hpp>

Inheritance diagram for EshelbianPlasticity::CGGUserPolynomialBase:
[legend]
Collaboration diagram for EshelbianPlasticity::CGGUserPolynomialBase:
[legend]

Public Types

using CachePhi = boost::tuple< int, int, MatrixDouble >
 

Public Member Functions

 CGGUserPolynomialBase (boost::shared_ptr< CachePhi > cache_phi=nullptr)
 
 ~CGGUserPolynomialBase ()=default
 
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, BaseFunctionUnknownInterface **iface) const
 
MoFEMErrorCode getValue (MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 

Private Member Functions

MoFEMErrorCode getValueHdivForCGGBubble (MatrixDouble &pts)
 

Private Attributes

MatrixDouble shapeFun
 
boost::shared_ptr< CachePhicachePhiPtr
 

Detailed Description

Definition at line 21 of file EshelbianPlasticity.hpp.

Member Typedef Documentation

◆ CachePhi

using EshelbianPlasticity::CGGUserPolynomialBase::CachePhi = boost::tuple<int, int, MatrixDouble>

Definition at line 23 of file EshelbianPlasticity.hpp.

Constructor & Destructor Documentation

◆ CGGUserPolynomialBase()

EshelbianPlasticity::CGGUserPolynomialBase::CGGUserPolynomialBase ( boost::shared_ptr< CachePhi cache_phi = nullptr)

Definition at line 1368 of file EshelbianPlasticity.cpp.

1370  : TetPolynomialBase(), cachePhiPtr(cache_phi_otr) {}

◆ ~CGGUserPolynomialBase()

EshelbianPlasticity::CGGUserPolynomialBase::~CGGUserPolynomialBase ( )
default

Member Function Documentation

◆ getValue()

MoFEMErrorCode EshelbianPlasticity::CGGUserPolynomialBase::getValue ( MatrixDouble &  pts,
boost::shared_ptr< BaseFunctionCtx ctx_ptr 
)

Definition at line 1380 of file EshelbianPlasticity.cpp.

1381  {
1383 
1384  this->cTx = ctx_ptr->getInterface<EntPolynomialBaseCtx>();
1385 
1386  int nb_gauss_pts = pts.size2();
1387  if (!nb_gauss_pts) {
1389  }
1390 
1391  if (pts.size1() < 3) {
1392  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1393  "Wrong dimension of pts, should be at least 3 rows with "
1394  "coordinates");
1395  }
1396 
1397  const auto base = this->cTx->bAse;
1398  EntitiesFieldData &data = this->cTx->dAta;
1399 
1400  switch (this->cTx->sPace) {
1401  case HDIV:
1403  break;
1404  case L2:
1405  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 4, false);
1406  CHKERR Tools::shapeFunMBTET(
1407  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
1408  &pts(0, 0), &pts(1, 0), &pts(2, 0), nb_gauss_pts);
1409  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(4, 3, false);
1410  std::copy(Tools::diffShapeFunMBTET.begin(), Tools::diffShapeFunMBTET.end(),
1411  data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
1412  this->cTx->basePolynomialsType0 = Legendre_polynomials;
1413  CHKERR getValueL2AinsworthBase(pts);
1414  break;
1415  default:
1416  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
1417  }
1418 
1420 }

◆ getValueHdivForCGGBubble()

MoFEMErrorCode EshelbianPlasticity::CGGUserPolynomialBase::getValueHdivForCGGBubble ( MatrixDouble &  pts)
private

number of integration points

Definition at line 1423 of file EshelbianPlasticity.cpp.

1423  {
1425 
1426  // This should be used only in case USER_BASE is selected
1427  if (this->cTx->bAse != USER_BASE) {
1428  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1429  "Wrong base, should be USER_BASE");
1430  }
1431  // get access to data structures on element
1432  EntitiesFieldData &data = this->cTx->dAta;
1433  // Get approximation order on element. Note that bubble functions are only
1434  // on tetrahedron.
1435  const int order = data.dataOnEntities[MBTET][0].getOrder();
1436  /// number of integration points
1437  const int nb_gauss_pts = pts.size2();
1438 
1439  auto check_cache = [this](int order, int nb_gauss_pts) -> bool {
1440  if (cachePhiPtr) {
1441  return cachePhiPtr->get<0>() == order &&
1442  cachePhiPtr->get<1>() == nb_gauss_pts;
1443  } else {
1444  return false;
1445  }
1446  };
1447 
1448  if (check_cache(order, nb_gauss_pts)) {
1449  auto &mat = cachePhiPtr->get<2>();
1450  auto &phi = data.dataOnEntities[MBTET][0].getN(USER_BASE);
1451  phi.resize(mat.size1(), mat.size2(), false);
1452  noalias(phi) = mat;
1453 
1454  } else {
1455  // calculate shape functions, i.e. barycentric coordinates
1456  shapeFun.resize(nb_gauss_pts, 4, false);
1457  CHKERR ShapeMBTET(&*shapeFun.data().begin(), &pts(0, 0), &pts(1, 0),
1458  &pts(2, 0), nb_gauss_pts);
1459  // derivatives of shape functions
1460  double diff_shape_fun[12];
1461  CHKERR ShapeDiffMBTET(diff_shape_fun);
1462 
1463  const int nb_base_functions = NBVOLUMETET_CCG_BUBBLE(order);
1464  // get base functions and set size
1465  MatrixDouble &phi = data.dataOnEntities[MBTET][0].getN(USER_BASE);
1466  phi.resize(nb_gauss_pts, 9 * nb_base_functions, false);
1467  // finally calculate base functions
1469  &phi(0, 0), &phi(0, 1), &phi(0, 2),
1470 
1471  &phi(0, 3), &phi(0, 4), &phi(0, 5),
1472 
1473  &phi(0, 6), &phi(0, 7), &phi(0, 8));
1474  CHKERR CGG_BubbleBase_MBTET(order, &shapeFun(0, 0), diff_shape_fun, t_phi,
1475  nb_gauss_pts);
1476 
1477  if (cachePhiPtr) {
1478  cachePhiPtr->get<0>() = order;
1479  cachePhiPtr->get<1>() = nb_gauss_pts;
1480  cachePhiPtr->get<2>().resize(phi.size1(), phi.size2(), false);
1481  noalias(cachePhiPtr->get<2>()) = phi;
1482  }
1483  }
1484 
1486 }

◆ query_interface()

MoFEMErrorCode EshelbianPlasticity::CGGUserPolynomialBase::query_interface ( boost::typeindex::type_index  type_index,
BaseFunctionUnknownInterface **  iface 
) const

Definition at line 1372 of file EshelbianPlasticity.cpp.

1374  {
1375  *iface = const_cast<CGGUserPolynomialBase *>(this);
1376  return 0;
1377 }

Member Data Documentation

◆ cachePhiPtr

boost::shared_ptr<CachePhi> EshelbianPlasticity::CGGUserPolynomialBase::cachePhiPtr
private

Definition at line 36 of file EshelbianPlasticity.hpp.

◆ shapeFun

MatrixDouble EshelbianPlasticity::CGGUserPolynomialBase::shapeFun
private

Definition at line 35 of file EshelbianPlasticity.hpp.


The documentation for this struct was generated from the following files:
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::EntPolynomialBaseCtx
Class used to pass element data to calculate base functions on tet,triangle,edge.
Definition: EntPolynomialBaseCtx.hpp:22
L2
@ L2
field with C-1 continuity
Definition: definitions.h:88
phi
static double phi
Definition: poisson_2d_dis_galerkin.cpp:30
MoFEM::Types::MatrixDouble
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
FTensor::Tensor2
Definition: Tensor2_value.hpp:16
EshelbianPlasticity::CGGUserPolynomialBase::shapeFun
MatrixDouble shapeFun
Definition: EshelbianPlasticity.hpp:35
USER_BASE
@ USER_BASE
user implemented approximation base
Definition: definitions.h:68
order
constexpr int order
Definition: dg_projection.cpp:18
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
EshelbianPlasticity::CGGUserPolynomialBase::getValueHdivForCGGBubble
MoFEMErrorCode getValueHdivForCGGBubble(MatrixDouble &pts)
Definition: EshelbianPlasticity.cpp:1423
EshelbianPlasticity::CGG_BubbleBase_MBTET
MoFEMErrorCode CGG_BubbleBase_MBTET(const int p, const double *N, const double *diffN, FTensor::Tensor2< FTensor::PackPtr< double *, 9 >, 3, 3 > &phi, const int gdim)
Calculate CGGT tonsorial bubble base.
Definition: CGGTonsorialBubbleBase.cpp:20
Legendre_polynomials
PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Legendre approximation basis.
Definition: base_functions.c:15
EshelbianPlasticity::CGGUserPolynomialBase::cachePhiPtr
boost::shared_ptr< CachePhi > cachePhiPtr
Definition: EshelbianPlasticity.hpp:36
NBVOLUMETET_CCG_BUBBLE
#define NBVOLUMETET_CCG_BUBBLE(P)
Bubble function for CGG H div space.
Definition: CGGTonsorialBubbleBase.hpp:19
ShapeDiffMBTET
PetscErrorCode ShapeDiffMBTET(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:319
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::EntitiesFieldData::dataOnEntities
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
Definition: EntitiesFieldData.hpp:57
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
MoFEM::EntitiesFieldData
data structure for finite element entity
Definition: EntitiesFieldData.hpp:40
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
HDIV
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
EshelbianPlasticity::CGGUserPolynomialBase::CGGUserPolynomialBase
CGGUserPolynomialBase(boost::shared_ptr< CachePhi > cache_phi=nullptr)
Definition: EshelbianPlasticity.cpp:1368
MoFEM::TetPolynomialBase
Calculate base functions on tetrahedral.
Definition: TetPolynomialBase.hpp:17
ShapeMBTET
PetscErrorCode ShapeMBTET(double *N, const double *G_X, const double *G_Y, const double *G_Z, int DIM)
calculate shape functions
Definition: fem_tools.c:306