v0.9.0
Public Member Functions | Private Member Functions | Private Attributes | List of all members
MoFEM::QuadPolynomialBase Struct Reference

Calculate base functions on triangle. More...

#include <src/approximation/QuadPolynomialBase.hpp>

Inheritance diagram for MoFEM::QuadPolynomialBase:
[legend]
Collaboration diagram for MoFEM::QuadPolynomialBase:
[legend]

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, BaseFunctionUnknownInterface **iface) const
 
 QuadPolynomialBase ()
 
 ~QuadPolynomialBase ()
 
MoFEMErrorCode getValue (MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::BaseFunction
 BaseFunction ()
 
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 getValueH1 (MatrixDouble &pts)
 

Private Attributes

EntPolynomialBaseCtxcTx
 
ublas::matrix< MatrixDoubleN_face_edge
 
ublas::vector< MatrixDoubleN_face_bubble
 
ublas::matrix< MatrixDoublediffN_face_edge
 
ublas::vector< MatrixDoublediffN_face_bubble
 

Detailed Description

Calculate base functions on triangle.

Definition at line 33 of file QuadPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ QuadPolynomialBase()

QuadPolynomialBase::QuadPolynomialBase ( )

Definition at line 22 of file QuadPolynomialBase.cpp.

22 {}

◆ ~QuadPolynomialBase()

QuadPolynomialBase::~QuadPolynomialBase ( )

Definition at line 23 of file QuadPolynomialBase.cpp.

23 {}

Member Function Documentation

◆ getValue()

MoFEMErrorCode QuadPolynomialBase::getValue ( MatrixDouble pts,
boost::shared_ptr< BaseFunctionCtx ctx_ptr 
)
virtual

Reimplemented from MoFEM::BaseFunction.

Definition at line 106 of file QuadPolynomialBase.cpp.

107  {
109 
111  CHKERR ctx_ptr->query_interface(IDD_QUAD_BASE_FUNCTION, &iface);
112  cTx = reinterpret_cast<EntPolynomialBaseCtx *>(iface);
113 
114  int nb_gauss_pts = pts.size2();
115  if (!nb_gauss_pts)
117 
118  if (pts.size1() < 2)
119  SETERRQ(
120  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
121  "Wrong dimension of pts, should be at least 3 rows with coordinates");
122 
123  const FieldApproximationBase base = cTx->bAse;
125  if (cTx->copyNodeBase != NOBASE)
126  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
127  "Shape base has to be on NOBASE", ApproximationBaseNames[base]);
128 
129  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
130  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
131  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(base) =
132  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(cTx->copyNodeBase);
133 
134  auto &base_shape = data.dataOnEntities[MBVERTEX][0].getN(base);
135  auto &diff_base = data.dataOnEntities[MBVERTEX][0].getDiffN(base);
136 
137  if (base_shape.size1() != pts.size2())
138  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
139  "Number of base functions integration points not equal number of "
140  "set integration point");
141  if (base_shape.size2() != 4)
142  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
143  "Number of shape functions should be four");
144  if (diff_base.size1() != pts.size2())
145  SETERRQ(
146  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
147  "Number of diff base functions integration points not equal number of "
148  "set integration point");
149  if (diff_base.size2() != 8)
150  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
151  "Number of shape functions should be four");
152 
153  switch (cTx->sPace) {
154  case H1:
155  CHKERR getValueH1(pts);
156  break;
157  default:
158  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
159  }
160 
162 }
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
const FieldApproximationBase copyNodeBase
static const char *const ApproximationBaseNames[]
Definition: definitions.h:158
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
MoFEMErrorCode getValueH1(MatrixDouble &pts)
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
continuous field
Definition: definitions.h:171
EntPolynomialBaseCtx * cTx
static const MOFEMuuid IDD_QUAD_BASE_FUNCTION

◆ getValueH1()

MoFEMErrorCode QuadPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 39 of file QuadPolynomialBase.cpp.

39  {
41 
43  const FieldApproximationBase base = cTx->bAse;
44  if (cTx->basePolynomialsType0 == NULL)
45  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
46  "Polynomial type not set");
47  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
48  double *diffL, const int dim) =
50 
51  int nb_gauss_pts = pts.size2();
52  auto &vert_dat = data.dataOnEntities[MBVERTEX][0];
53 
54  if (data.spacesOnEntities[MBEDGE].test(H1)) {
55  // edges
56  if (data.dataOnEntities[MBEDGE].size() != 4)
57  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
58  "should be four edges on quad");
59 
60  int sense[4], order[4];
61  double *H1edgeN[4], *diffH1edgeN[4];
62  for (int ee = 0; ee != 4; ++ee) {
63  auto &ent_dat = data.dataOnEntities[MBEDGE][ee];
64  if (ent_dat.getSense() == 0)
65  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "sense not set");
66 
67  sense[ee] = ent_dat.getSense();
68  order[ee] = ent_dat.getDataOrder();
69  int nb_dofs = NBEDGE_H1(ent_dat.getDataOrder());
70  ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
71  ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
72  H1edgeN[ee] = &*ent_dat.getN(base).data().begin();
73  diffH1edgeN[ee] = &*ent_dat.getDiffN(base).data().begin();
74  }
76  sense, order, &*vert_dat.getN(base).data().begin(),
77  &*vert_dat.getDiffN(base).data().begin(), H1edgeN, diffH1edgeN,
78  nb_gauss_pts, base_polynomials);
79  }
80 
81  if (data.spacesOnEntities[MBQUAD].test(H1)) {
82 
83  // face
84  if (data.dataOnEntities[MBQUAD].size() != 1)
85  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
86  "should be one quad to store bubble base on quad");
87 
88  auto &ent_dat = data.dataOnEntities[MBQUAD][0];
89  int nb_dofs = NBFACEQUAD_H1(ent_dat.getDataOrder());
90  ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
91  ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
92  int face_nodes[] = {0, 1, 2, 3};
94  face_nodes, ent_dat.getDataOrder(),
95  &*vert_dat.getN(base).data().begin(),
96  &*vert_dat.getDiffN(base).data().begin(),
97  &*ent_dat.getN(base).data().begin(),
98  &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts,
99  base_polynomials);
100  }
101 
103 }
PetscErrorCode H1_EdgeShapeFunctions_MBQUAD(int *sense, int *p, double *N, double *diffN, double *edgeN[4], double *diff_edgeN[4], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:937
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
const FieldApproximationBase bAse
PetscErrorCode H1_QuadShapeFunctions_MBQUAD(int *faces_nodes, int p, double *N, double *diffN, double *faceN, double *diff_faceN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:827
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
FieldApproximationBase
approximation base
Definition: definitions.h:144
#define NBFACEQUAD_H1(P)
Number of base functions on quad for H1 space.
DataForcesAndSourcesCore & dAta
#define CHKERR
Inline error check.
Definition: definitions.h:596
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
continuous field
Definition: definitions.h:171
EntPolynomialBaseCtx * cTx

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 25 of file QuadPolynomialBase.cpp.

26  {
28  *iface = NULL;
29  if (uuid == IDD_QUAD_BASE_FUNCTION) {
30  *iface = const_cast<QuadPolynomialBase *>(this);
32  } else
33  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
34 
37 }
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, MoFEM::BaseFunctionUnknownInterface **iface) const
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
static const MOFEMuuid IDD_QUAD_BASE_FUNCTION

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* MoFEM::QuadPolynomialBase::cTx
private

Definition at line 45 of file QuadPolynomialBase.hpp.

◆ diffN_face_bubble

ublas::vector<MatrixDouble> MoFEM::QuadPolynomialBase::diffN_face_bubble
private

Definition at line 52 of file QuadPolynomialBase.hpp.

◆ diffN_face_edge

ublas::matrix<MatrixDouble> MoFEM::QuadPolynomialBase::diffN_face_edge
private

Definition at line 51 of file QuadPolynomialBase.hpp.

◆ N_face_bubble

ublas::vector<MatrixDouble> MoFEM::QuadPolynomialBase::N_face_bubble
private

Definition at line 50 of file QuadPolynomialBase.hpp.

◆ N_face_edge

ublas::matrix<MatrixDouble> MoFEM::QuadPolynomialBase::N_face_edge
private

Definition at line 49 of file QuadPolynomialBase.hpp.


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