v0.9.0
Public Member Functions | Private Member Functions | Private Attributes | List of all members
EshelbianPlasticity::CGGUserPolynomialBase Struct Reference
Inheritance diagram for EshelbianPlasticity::CGGUserPolynomialBase:
[legend]
Collaboration diagram for EshelbianPlasticity::CGGUserPolynomialBase:
[legend]

Public Member Functions

 CGGUserPolynomialBase ()
 
 ~CGGUserPolynomialBase ()
 
MoFEMErrorCode query_interface (const MOFEMuuid &uuid, BaseFunctionUnknownInterface **iface) const
 
MoFEMErrorCode getValue (MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::BaseFunction
 BaseFunction ()
 
virtual MoFEMErrorCode getValue (MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
virtual MoFEMErrorCode getValue (MatrixDouble &pts_x, MatrixDouble &pts_t, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::BaseFunctionUnknownInterface
virtual ~BaseFunctionUnknownInterface ()=default
 

Private Member Functions

MoFEMErrorCode getValueHdivForCGGBubble (MatrixDouble &pts)
 

Private Attributes

EntPolynomialBaseCtxcTx
 
MatrixDouble shapeFun
 

Detailed Description

Definition at line 550 of file EshelbianPlasticity.cpp.

Constructor & Destructor Documentation

◆ CGGUserPolynomialBase()

EshelbianPlasticity::CGGUserPolynomialBase::CGGUserPolynomialBase ( )

Definition at line 552 of file EshelbianPlasticity.cpp.

552 {}

◆ ~CGGUserPolynomialBase()

EshelbianPlasticity::CGGUserPolynomialBase::~CGGUserPolynomialBase ( )

Definition at line 553 of file EshelbianPlasticity.cpp.

553 {}

Member Function Documentation

◆ getValue()

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

Definition at line 569 of file EshelbianPlasticity.cpp.

570  {
572 
574  CHKERR ctx_ptr->query_interface(IDD_TET_BASE_FUNCTION, &iface);
575  cTx = reinterpret_cast<EntPolynomialBaseCtx *>(iface);
576 
577  int nb_gauss_pts = pts.size2();
578  if (!nb_gauss_pts) {
580  }
581 
582  if (pts.size1() < 3) {
583  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
584  "Wrong dimension of pts, should be at least 3 rows with "
585  "coordinates");
586  }
587 
588  switch (cTx->sPace) {
589  case HDIV:
591  break;
592  default:
593  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
594  }
595 
597  }
field with continuous normal traction
Definition: definitions.h:177
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:505
MoFEMErrorCode getValueHdivForCGGBubble(MatrixDouble &pts)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:512
static const MOFEMuuid IDD_TET_BASE_FUNCTION
#define CHKERR
Inline error check.
Definition: definitions.h:600

◆ getValueHdivForCGGBubble()

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

Definition at line 604 of file EshelbianPlasticity.cpp.

604  {
606 
607  // This should be used only in case USER_BASE is selected
608  if (cTx->bAse != USER_BASE) {
609  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
610  "Wrong base, should be USER_BASE");
611  }
612  // get access to data structures on element
614  // Get approximation order on element. Note that bubble functions are only
615  // on tetrahedron.
616  const int order = data.dataOnEntities[MBTET][0].getDataOrder();
617  /// number of integration points
618  const int nb_gauss_pts = pts.size2();
619  // number of base functions
620 
621  // calculate shape functions, i.e. barycentric coordinates
622  shapeFun.resize(nb_gauss_pts, 4, false);
623  CHKERR ShapeMBTET(&*shapeFun.data().begin(), &pts(0, 0), &pts(1, 0),
624  &pts(2, 0), nb_gauss_pts);
625  // direvatives of shape functions
626  double diff_shape_fun[12];
627  CHKERR ShapeDiffMBTET(diff_shape_fun);
628 
629  const int nb_base_functions = NBVOLUMETET_CCG_BUBBLE(order);
630  // get base functions and set size
631  MatrixDouble &phi = data.dataOnEntities[MBTET][0].getN(USER_BASE);
632  phi.resize(nb_gauss_pts, 9 * nb_base_functions, false);
633  // finally calculate base functions
635  &phi(0, 0), &phi(0, 1), &phi(0, 2),
636 
637  &phi(0, 3), &phi(0, 4), &phi(0, 5),
638 
639  &phi(0, 6), &phi(0, 7), &phi(0, 8));
640  CHKERR CGG_BubbleBase_MBTET(order, &shapeFun(0, 0), diff_shape_fun, t_phi,
641  nb_gauss_pts);
642 
644  }
user implemented approximation base
Definition: definitions.h:158
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
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:318
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:481
PetscErrorCode ShapeDiffMBTET(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:331
DataForcesAndSourcesCore & dAta
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.
#define CHKERR
Inline error check.
Definition: definitions.h:600
#define NBVOLUMETET_CCG_BUBBLE(P)
Bubble function for CGG H div space.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:411

◆ query_interface()

MoFEMErrorCode EshelbianPlasticity::CGGUserPolynomialBase::query_interface ( const MOFEMuuid uuid,
BaseFunctionUnknownInterface **  iface 
) const
virtual

Reimplemented from MoFEM::BaseFunction.

Definition at line 555 of file EshelbianPlasticity.cpp.

556  {
558  *iface = NULL;
559  if (uuid == IDD_TET_BASE_FUNCTION) {
560  *iface = const_cast<CGGUserPolynomialBase *>(this);
562  } else {
563  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
564  }
565  CHKERR BaseFunction::query_interface(uuid, iface);
567  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:481
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:512
static const MOFEMuuid IDD_TET_BASE_FUNCTION
#define CHKERR
Inline error check.
Definition: definitions.h:600
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:411

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* EshelbianPlasticity::CGGUserPolynomialBase::cTx
private

Definition at line 600 of file EshelbianPlasticity.cpp.

◆ shapeFun

MatrixDouble EshelbianPlasticity::CGGUserPolynomialBase::shapeFun
private

Definition at line 602 of file EshelbianPlasticity.cpp.


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