v0.9.2
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 547 of file EshelbianPlasticity.cpp.

Constructor & Destructor Documentation

◆ CGGUserPolynomialBase()

EshelbianPlasticity::CGGUserPolynomialBase::CGGUserPolynomialBase ( )

Definition at line 549 of file EshelbianPlasticity.cpp.

549 {}

◆ ~CGGUserPolynomialBase()

EshelbianPlasticity::CGGUserPolynomialBase::~CGGUserPolynomialBase ( )

Definition at line 550 of file EshelbianPlasticity.cpp.

550 {}

Member Function Documentation

◆ getValue()

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

Definition at line 566 of file EshelbianPlasticity.cpp.

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

◆ getValueHdivForCGGBubble()

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

Definition at line 601 of file EshelbianPlasticity.cpp.

601  {
603 
604  // This should be used only in case USER_BASE is selected
605  if (cTx->bAse != USER_BASE) {
606  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
607  "Wrong base, should be USER_BASE");
608  }
609  // get access to data structures on element
611  // Get approximation order on element. Note that bubble functions are only
612  // on tetrahedron.
613  const int order = data.dataOnEntities[MBTET][0].getDataOrder();
614  /// number of integration points
615  const int nb_gauss_pts = pts.size2();
616  // number of base functions
617 
618  // calculate shape functions, i.e. barycentric coordinates
619  shapeFun.resize(nb_gauss_pts, 4, false);
620  CHKERR ShapeMBTET(&*shapeFun.data().begin(), &pts(0, 0), &pts(1, 0),
621  &pts(2, 0), nb_gauss_pts);
622  // direvatives of shape functions
623  double diff_shape_fun[12];
624  CHKERR ShapeDiffMBTET(diff_shape_fun);
625 
626  const int nb_base_functions = NBVOLUMETET_CCG_BUBBLE(order);
627  // get base functions and set size
628  MatrixDouble &phi = data.dataOnEntities[MBTET][0].getN(USER_BASE);
629  phi.resize(nb_gauss_pts, 9 * nb_base_functions, false);
630  // finally calculate base functions
632  &phi(0, 0), &phi(0, 1), &phi(0, 2),
633 
634  &phi(0, 3), &phi(0, 4), &phi(0, 5),
635 
636  &phi(0, 6), &phi(0, 7), &phi(0, 8));
637  CHKERR CGG_BubbleBase_MBTET(order, &shapeFun(0, 0), diff_shape_fun, t_phi,
638  nb_gauss_pts);
639 
641  }
user implemented approximation base
Definition: definitions.h:160
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:75
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
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:602
constexpr int order
#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:413

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 552 of file EshelbianPlasticity.cpp.

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

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* EshelbianPlasticity::CGGUserPolynomialBase::cTx
private

Definition at line 597 of file EshelbianPlasticity.cpp.

◆ shapeFun

MatrixDouble EshelbianPlasticity::CGGUserPolynomialBase::shapeFun
private

Definition at line 599 of file EshelbianPlasticity.cpp.


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