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

Calculate base functions on triangle. More...

#include <src/approximation/TriPolynomialBase.hpp>

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

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, BaseFunctionUnknownInterface **iface) const
 
 TriPolynomialBase ()
 
 ~TriPolynomialBase ()
 
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)
 
MoFEMErrorCode getValueH1AinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueH1BernsteinBezierBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2 (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdiv (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurl (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdivAinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurlAinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdivDemkowiczBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurlDemkowiczBase (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 30 of file TriPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ TriPolynomialBase()

TriPolynomialBase::TriPolynomialBase ( )

Definition at line 22 of file TriPolynomialBase.cpp.

22 {}

◆ ~TriPolynomialBase()

TriPolynomialBase::~TriPolynomialBase ( )

Definition at line 23 of file TriPolynomialBase.cpp.

23 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 614 of file TriPolynomialBase.cpp.

615  {
617 
619  CHKERR ctx_ptr->query_interface(IDD_TRI_BASE_FUNCTION, &iface);
620  cTx = reinterpret_cast<EntPolynomialBaseCtx *>(iface);
621 
622  int nb_gauss_pts = pts.size2();
623  if (!nb_gauss_pts) {
625  }
626 
627  if (pts.size1() < 2)
628  SETERRQ(
629  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
630  "Wrong dimension of pts, should be at least 3 rows with coordinates");
631 
632  const FieldApproximationBase base = cTx->bAse;
634 
635  if (base != AINSWORTH_BERNSTEIN_BEZIER_BASE) {
636  if (cTx->copyNodeBase == LASTBASE) {
637  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 3,
638  false);
640  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
641  &pts(0, 0), &pts(1, 0), nb_gauss_pts);
642  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(3, 2, false);
644  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
645  } else {
646  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
647  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
648  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(base) =
649  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(cTx->copyNodeBase);
650  }
651  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
652  (unsigned int)nb_gauss_pts) {
653  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
654  "Base functions or nodes has wrong number of integration points "
655  "for base %s",
656  ApproximationBaseNames[base]);
657  }
658  }
659  auto set_node_derivative_for_all_gauss_pts = [&]() {
661  // In linear geometry derivatives are constant,
662  // this in expense of efficiency makes implementation
663  // consistent between vertices and other types of entities
664  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 6,
665  false);
666  for (int gg = 0; gg != nb_gauss_pts; ++gg)
667  std::copy(Tools::diffShapeFunMBTRI.begin(),
669  &data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg, 0));
671  };
672 
673  CHKERR set_node_derivative_for_all_gauss_pts();
674 
675  switch (cTx->sPace) {
676  case H1:
677  CHKERR getValueH1(pts);
678  break;
679  case HDIV:
680  CHKERR getValueHdiv(pts);
681  break;
682  case HCURL:
683  CHKERR getValueHcurl(pts);
684  break;
685  case L2:
686  CHKERR getValueL2(pts);
687  break;
688  default:
689  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
690  }
691 
693 }
field with continuous normal traction
Definition: definitions.h:173
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
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
MoFEMErrorCode getValueH1(MatrixDouble &pts)
const FieldApproximationBase copyNodeBase
static const char *const ApproximationBaseNames[]
Definition: definitions.h:158
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getValueL2(MatrixDouble &pts)
PetscErrorCode ShapeDiffMBTRI(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:206
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:172
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition: Tools.hpp:87
#define CHKERR
Inline error check.
Definition: definitions.h:596
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition: fem_tools.c:194
static const MOFEMuuid IDD_TRI_BASE_FUNCTION
#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
field with C-1 continuity
Definition: definitions.h:174

◆ getValueH1()

MoFEMErrorCode TriPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 40 of file TriPolynomialBase.cpp.

40  {
42 
43  switch (cTx->bAse) {
47  break;
50  break;
51  default:
52  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
53  }
54 
56 }
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
EntPolynomialBaseCtx * cTx
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:146
#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
MoFEMErrorCode getValueH1BernsteinBezierBase(MatrixDouble &pts)
MoFEMErrorCode getValueH1AinsworthBase(MatrixDouble &pts)

◆ getValueH1AinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueH1AinsworthBase ( MatrixDouble pts)
private

Definition at line 198 of file TriPolynomialBase.cpp.

198  {
200 
202  const FieldApproximationBase base = cTx->bAse;
203  if (cTx->basePolynomialsType0 == NULL)
204  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
205  "Polynomial type not set");
206  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
207  double *diffL, const int dim) =
209 
210  int nb_gauss_pts = pts.size2();
211 
212  if (data.spacesOnEntities[MBEDGE].test(H1)) {
213  // edges
214  if (data.dataOnEntities[MBEDGE].size() != 3)
215  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
216 
217  int sense[3], order[3];
218  double *H1edgeN[3], *diffH1edgeN[3];
219  for (int ee = 0; ee < 3; ee++) {
220  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0)
221  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
222  "sense of the edge unknown");
223 
224  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
225  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
226  int nb_dofs = NBEDGE_H1(data.dataOnEntities[MBEDGE][ee].getDataOrder());
227  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
228  false);
229  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,
230  2 * nb_dofs, false);
231  H1edgeN[ee] = &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
232  diffH1edgeN[ee] =
233  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
234  }
236  sense, order,
237  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
238  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
239  H1edgeN, diffH1edgeN, nb_gauss_pts, base_polynomials);
240  }
241 
242  if (data.spacesOnEntities[MBTRI].test(H1)) {
243  // face
244  if (data.dataOnEntities[MBTRI].size() != 1) {
245  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
246  }
247  int nb_dofs = NBFACETRI_H1(data.dataOnEntities[MBTRI][0].getDataOrder());
248  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, nb_dofs,
249  false);
250  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
251  2 * nb_dofs, false);
252  const int face_nodes[] = {0, 1, 2};
254  face_nodes, data.dataOnEntities[MBTRI][0].getDataOrder(),
255  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
256  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
257  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
258  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
259  nb_gauss_pts, base_polynomials);
260  }
261 
263 }
#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,...
PetscErrorCode H1_EdgeShapeFunctions_MBTRI(int *sense, int *p, double *N, double *diffN, double *edgeN[3], double *diff_edgeN[3], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
H1_EdgeShapeFunctions_MBTRI.
Definition: h1.c:27
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
PetscErrorCode H1_FaceShapeFunctions_MBTRI(const int *face_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:162
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.
#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

◆ getValueH1BernsteinBezierBase()

MoFEMErrorCode TriPolynomialBase::getValueH1BernsteinBezierBase ( MatrixDouble pts)
private

Definition at line 59 of file TriPolynomialBase.cpp.

59  {
62  const std::string &field_name = cTx->fieldName;
63  int nb_gauss_pts = pts.size2();
64 
65  auto get_alpha = [field_name](auto &data) -> MatrixInt & {
66  auto &ptr = data.getBBAlphaIndicesSharedPtr(field_name);
67  if (!ptr)
68  ptr.reset(new MatrixInt());
69  return *ptr;
70  };
71 
72  auto get_base = [field_name](auto &data) -> MatrixDouble & {
73  auto &ptr = data.getBBNSharedPtr(field_name);
74  if (!ptr)
75  ptr.reset(new MatrixDouble());
76  return *ptr;
77  };
78 
79  auto get_diff_base = [field_name](auto &data) -> MatrixDouble & {
80  auto &ptr = data.getBBDiffNSharedPtr(field_name);
81  if (!ptr)
82  ptr.reset(new MatrixDouble());
83  return *ptr;
84  };
85 
86  auto &vert_get_n = get_base(data.dataOnEntities[MBVERTEX][0]);
87  auto &vert_get_diff_n = get_diff_base(data.dataOnEntities[MBVERTEX][0]);
88  vert_get_n.resize(nb_gauss_pts, 3, false);
89  vert_get_diff_n.resize(nb_gauss_pts, 6, false);
90 
91  if (data.dataOnEntities[MBVERTEX][0].getN(NOBASE).size1() !=
92  (unsigned int)nb_gauss_pts)
93  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
94  "Base functions or nodes has wrong number of integration points "
95  "for base %s",
97  auto &lambda = data.dataOnEntities[MBVERTEX][0].getN(NOBASE);
98 
99  auto &vertex_alpha = get_alpha(data.dataOnEntities[MBVERTEX][0]);
100  vertex_alpha.resize(3, 3, false);
101  vertex_alpha.clear();
102  for (int n = 0; n != 3; ++n)
103  vertex_alpha(n, n) = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[n];
104 
106  1, lambda.size1(), vertex_alpha.size1(), &vertex_alpha(0, 0),
107  &lambda(0, 0), Tools::diffShapeFunMBTRI.data(), &vert_get_n(0, 0),
108  &vert_get_diff_n(0, 0));
109  for (int n = 0; n != 3; ++n) {
110  const int f = boost::math::factorial<double>(
111  data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[n]);
112  for (int g = 0; g != nb_gauss_pts; ++g) {
113  vert_get_n(g, n) *= f;
114  for (int d = 0; d != 2; ++d)
115  vert_get_diff_n(g, 2 * n + d) *= f;
116  }
117  }
118 
119  // edges
120  if (data.spacesOnEntities[MBEDGE].test(H1)) {
121  if (data.dataOnEntities[MBEDGE].size() != 3)
122  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
123  "Wrong size of ent data");
124 
125  constexpr int edges_nodes[3][2] = {{0, 1}, {1, 2}, {2, 0}};
126  for (int ee = 0; ee != 3; ++ee) {
127  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0)
128  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
129  "Sense of the edge unknown");
130  const int sense = data.dataOnEntities[MBEDGE][ee].getSense();
131  const int order = data.dataOnEntities[MBEDGE][ee].getDataOrder();
132  const int nb_dofs =
133  NBEDGE_H1(data.dataOnEntities[MBEDGE][ee].getDataOrder());
134 
135  auto &get_n = get_base(data.dataOnEntities[MBEDGE][ee]);
136  auto &get_diff_n = get_diff_base(data.dataOnEntities[MBEDGE][ee]);
137  get_n.resize(nb_gauss_pts, nb_dofs, false);
138  get_diff_n.resize(nb_gauss_pts, 2 * nb_dofs, false);
139 
140  if (nb_dofs) {
141  auto &edge_alpha = get_alpha(data.dataOnEntities[MBEDGE][ee]);
142  edge_alpha.resize(nb_dofs, 3, false);
144  &edge_alpha(0, 0));
145  if (sense == -1)
146  for (int i = 0; i != edge_alpha.size1(); ++i) {
147  int a = edge_alpha(i, edges_nodes[ee][0]);
148  edge_alpha(i, edges_nodes[ee][0]) =
149  edge_alpha(i, edges_nodes[ee][1]);
150  edge_alpha(i, edges_nodes[ee][1]) = a;
151  }
153  order, lambda.size1(), edge_alpha.size1(), &edge_alpha(0, 0),
154  &lambda(0, 0), Tools::diffShapeFunMBTRI.data(), &get_n(0, 0),
155  &get_diff_n(0, 0));
156  }
157  }
158  } else {
159  for (int ee = 0; ee != 3; ++ee) {
160  auto &get_n = get_base(data.dataOnEntities[MBEDGE][ee]);
161  auto &get_diff_n = get_diff_base(data.dataOnEntities[MBEDGE][ee]);
162  get_n.resize(nb_gauss_pts, 0, false);
163  get_diff_n.resize(nb_gauss_pts, 0, false);
164  }
165  }
166 
167  // face
168  if (data.spacesOnEntities[MBTRI].test(H1)) {
169  if (data.dataOnEntities[MBTRI].size() != 1)
170  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
171  "Wrong size ent of ent data");
172 
173  const int order = data.dataOnEntities[MBTRI][0].getDataOrder();
174  const int nb_dofs = NBFACETRI_H1(order);
175  auto &get_n = get_base(data.dataOnEntities[MBTRI][0]);
176  auto &get_diff_n = get_diff_base(data.dataOnEntities[MBTRI][0]);
177  get_n.resize(nb_gauss_pts, nb_dofs, false);
178  get_diff_n.resize(nb_gauss_pts, 2 * nb_dofs, false);
179  if (nb_dofs) {
180  auto &face_alpha = get_alpha(data.dataOnEntities[MBTRI][0]);
181  face_alpha.resize(nb_dofs, 3, false);
182  CHKERR BernsteinBezier::generateIndicesTriTri(order, &face_alpha(0, 0));
184  order, lambda.size1(), face_alpha.size1(), &face_alpha(0, 0),
185  &lambda(0, 0), Tools::diffShapeFunMBTRI.data(), &get_n(0, 0),
186  &get_diff_n(0, 0));
187  }
188  } else {
189  auto &get_n = get_base(data.dataOnEntities[MBTRI][0]);
190  auto &get_diff_n = get_diff_base(data.dataOnEntities[MBTRI][0]);
191  get_n.resize(nb_gauss_pts, 0, false);
192  get_diff_n.resize(nb_gauss_pts, 0, false);
193  }
194 
196 }
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
static MoFEMErrorCode generateIndicesEdgeTri(const int N[], int *alpha[])
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
ublas::matrix< int, ublas::row_major, IntAllocator > MatrixInt
Definition: Types.hpp:73
static const char *const ApproximationBaseNames[]
Definition: definitions.h:158
EntPolynomialBaseCtx * cTx
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
DataForcesAndSourcesCore & dAta
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.
static MoFEMErrorCode baseFunctionsTri(const int N, const int gdim, const int n_alpha, const int *alpha, const double *lambda, const double *grad_lambda, double *base, double *grad_base)
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition: Tools.hpp:87
#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
static MoFEMErrorCode generateIndicesTriTri(const int N, int *alpha)

◆ getValueHcurl()

MoFEMErrorCode TriPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 595 of file TriPolynomialBase.cpp.

595  {
597 
598  switch (cTx->bAse) {
602  break;
605  break;
606  default:
607  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
608  }
609 
611 }
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:146
#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

◆ getValueHcurlAinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlAinsworthBase ( MatrixDouble pts)
private

Definition at line 415 of file TriPolynomialBase.cpp.

415  {
417 
419  const FieldApproximationBase base = cTx->bAse;
420  if (data.dataOnEntities[MBTRI].size() != 1)
421  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
422 
423  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
424  double *diffL, const int dim) =
426 
427  int nb_gauss_pts = pts.size2();
428 
429  // Calculation H-curl on triangle faces
430  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
431  if (data.dataOnEntities[MBEDGE].size() != 3) {
432  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
433  }
434  int sense[3], order[3];
435  double *hcurl_edge_n[3];
436  double *diff_hcurl_edge_n[3];
437  for (int ee = 0; ee < 3; ee++) {
438  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
439  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
440  "data inconsistency");
441  }
442  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
443  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
444  int nb_dofs = NBEDGE_AINSWORTH_HCURL(
445  data.dataOnEntities[MBEDGE][ee].getDataOrder());
446  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
447  3 * nb_dofs, false);
448  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
449  nb_gauss_pts, 2 * 3 * nb_dofs, false);
450  hcurl_edge_n[ee] =
451  &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
452  diff_hcurl_edge_n[ee] =
453  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
454  }
456  sense, order,
457  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
458  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
459  hcurl_edge_n, diff_hcurl_edge_n, nb_gauss_pts, base_polynomials);
460  } else {
461  for (int ee = 0; ee < 3; ee++) {
462  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, 0, false);
463  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts, 0,
464  false);
465  }
466  }
467 
468  if (data.spacesOnEntities[MBTRI].test(HCURL)) {
469 
470  // cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << endl;
471  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
472  //
473  // face
474  if (data.dataOnEntities[MBTRI].size() != 1) {
475  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
476  }
477  int order = data.dataOnEntities[MBTRI][0].getDataOrder();
478  int nb_dofs = NBFACETRI_AINSWORTH_HCURL(order);
479  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
480  false);
481  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
482  3 * 2 * nb_dofs, false);
483  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
484  int face_nodes[] = {0, 1, 2};
486  face_nodes, order,
487  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
488  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
489  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
490  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
491  nb_gauss_pts, base_polynomials);
492  // cerr << data.dataOnEntities[MBTRI][0].getN(base) << endl;
493  } else {
494  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
495  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
496  }
497 
499 }
MoFEMErrorCode Hcurl_Ainsworth_FaceFunctions_MBTET_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face H-curl functions.
Definition: Hcurl.cpp:1249
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
MoFEMErrorCode Hcurl_Ainsworth_EdgeBaseFunctions_MBTET_ON_FACE(int *sense, int *p, double *N, double *diffN, double *edgeN[], double *diff_edgeN[], int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Edge based H-curl base functions on face.
Definition: Hcurl.cpp:237
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
#define NBEDGE_AINSWORTH_HCURL(P)
field with continuous tangents
Definition: definitions.h:172
#define CHKERR
Inline error check.
Definition: definitions.h:596
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
#define NBFACETRI_AINSWORTH_HCURL(P)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 502 of file TriPolynomialBase.cpp.

502  {
504 
506  const FieldApproximationBase base = cTx->bAse;
507 
508  int nb_gauss_pts = pts.size2();
509 
510  // Calculation H-curl on triangle faces
511  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
512 
513  if (data.dataOnEntities[MBEDGE].size() != 3)
514  SETERRQ1(
515  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
516  "wrong number of data structures on edges, should be three but is %d",
517  data.dataOnEntities[MBEDGE].size());
518 
519  int sense[3], order[3];
520  double *hcurl_edge_n[3];
521  double *diff_hcurl_edge_n[3];
522 
523  for (int ee = 0; ee != 3; ++ee) {
524 
525  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0)
526  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
527  "orientation (sense) of edge is not set");
528 
529  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
530  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
531  int nb_dofs = NBEDGE_DEMKOWICZ_HCURL(
532  data.dataOnEntities[MBEDGE][ee].getDataOrder());
533  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
534  3 * nb_dofs, false);
535  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
536  nb_gauss_pts, 2 * 3 * nb_dofs, false);
537 
538  hcurl_edge_n[ee] =
539  &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
540  diff_hcurl_edge_n[ee] =
541  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
542  }
543 
545  sense, order,
546  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
547  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
548  hcurl_edge_n, diff_hcurl_edge_n, nb_gauss_pts);
549 
550  } else {
551 
552  // No DOFs on faces, resize base function matrices, indicating that no
553  // dofs on them.
554  for (int ee = 0; ee != 3; ++ee) {
555  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, 0, false);
556  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts, 0,
557  false);
558  }
559  }
560 
561  if (data.spacesOnEntities[MBTRI].test(HCURL)) {
562 
563  // face
564  if (data.dataOnEntities[MBTRI].size() != 1)
565  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
566  "No data struture to keep base functions on face");
567 
568  int order = data.dataOnEntities[MBTRI][0].getDataOrder();
569  int nb_dofs = NBFACETRI_DEMKOWICZ_HCURL(order);
570  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
571  false);
572  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
573  3 * 2 * nb_dofs, false);
574 
575  int face_nodes[] = {0, 1, 2};
577  face_nodes, order,
578  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
579  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
580  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
581  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
582  nb_gauss_pts);
583 
584  } else {
585 
586  // No DOFs on faces, resize base function matrices, indicating that no
587  // dofs on them.
588  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
589  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
590  }
591 
593 }
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
MoFEMErrorCode Hcurl_Demkowicz_FaceBaseFunctions_MBTRI(int *faces_nodes, int p, double *n, double *diff_n, double *phi, double *diff_phi, int nb_integration_pts)
Face base interior function.
Definition: Hcurl.cpp:2430
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:172
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define NBFACETRI_DEMKOWICZ_HCURL(P)
#define NBEDGE_DEMKOWICZ_HCURL(P)
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
MoFEMErrorCode Hcurl_Demkowicz_EdgeBaseFunctions_MBTRI(int *sense, int *p, double *n, double *diff_n, double *phi[], double *diff_phi[], int nb_integration_pts)
Edge based H-curl base functions on teriangle.
Definition: Hcurl.cpp:2105

◆ getValueHdiv()

MoFEMErrorCode TriPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 398 of file TriPolynomialBase.cpp.

398  {
400 
401  switch (cTx->bAse) {
404  return getValueHdivAinsworthBase(pts);
406  return getValueHdivDemkowiczBase(pts);
407  default:
408  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
409  }
410 
412 }
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
EntPolynomialBaseCtx * cTx
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:146
MoFEMErrorCode getValueHdivDemkowiczBase(MatrixDouble &pts)
MoFEMErrorCode getValueHdivAinsworthBase(MatrixDouble &pts)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ getValueHdivAinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueHdivAinsworthBase ( MatrixDouble pts)
private

Definition at line 297 of file TriPolynomialBase.cpp.

297  {
299 
301  const FieldApproximationBase base = cTx->bAse;
302  if (cTx->basePolynomialsType0 == NULL)
303  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
304  "Polynomial type not set");
305  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
306  double *diffL, const int dim) =
308 
309  int nb_gauss_pts = pts.size2();
310 
311  double *PHI_f_e[3];
312  double *PHI_f;
313 
314  N_face_edge.resize(1, 3, false);
315  N_face_bubble.resize(1, false);
316  int face_order = data.dataOnEntities[MBTRI][0].getDataOrder();
317  // three edges on face
318  for (int ee = 0; ee < 3; ee++) {
319  N_face_edge(0, ee).resize(
320  nb_gauss_pts, 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(face_order), false);
321  PHI_f_e[ee] = &((N_face_edge(0, ee))(0, 0));
322  }
323  N_face_bubble[0].resize(nb_gauss_pts,
324  3 * NBFACETRI_AINSWORTH_FACE_HDIV(face_order), false);
325  PHI_f = &*(N_face_bubble[0].data().begin());
326 
327  int face_nodes[3] = {0, 1, 2};
329  face_nodes, face_order,
330  &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0), NULL, PHI_f_e, NULL,
331  nb_gauss_pts, 3, base_polynomials);
333  face_nodes, face_order,
334  &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0), NULL, PHI_f, NULL,
335  nb_gauss_pts, 3, base_polynomials);
336 
337  // set shape functions into data structure
338  if (data.dataOnEntities[MBTRI].size() != 1) {
339  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
340  }
341  data.dataOnEntities[MBTRI][0].getN(base).resize(
342  nb_gauss_pts, 3 * NBFACETRI_AINSWORTH_HDIV(face_order), false);
343  int col = 0;
344  for (int oo = 0; oo < face_order; oo++) {
345  for (int ee = 0; ee < 3; ee++) {
346  for (int dd = 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(oo);
347  dd < 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(oo + 1); dd++, col++) {
348  for (int gg = 0; gg < nb_gauss_pts; gg++) {
349  data.dataOnEntities[MBTRI][0].getN(base)(gg, col) =
350  N_face_edge(0, ee)(gg, dd);
351  }
352  }
353  }
354  for (int dd = 3 * NBFACETRI_AINSWORTH_FACE_HDIV(oo);
355  dd < 3 * NBFACETRI_AINSWORTH_FACE_HDIV(oo + 1); dd++, col++) {
356  for (int gg = 0; gg < nb_gauss_pts; gg++) {
357  data.dataOnEntities[MBTRI][0].getN(base)(gg, col) =
358  N_face_bubble[0](gg, dd);
359  }
360  }
361  }
362 
364 }
#define NBFACETRI_AINSWORTH_EDGE_HDIV(P)
ublas::vector< MatrixDouble > N_face_bubble
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
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
EntPolynomialBaseCtx * cTx
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
ublas::matrix< MatrixDouble > N_face_edge
#define CHKERR
Inline error check.
Definition: definitions.h:596
MoFEMErrorCode Hdiv_Ainsworth_EdgeFaceShapeFunctions_MBTET_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f_e[3], double *diff_phi_f_e[3], int gdim, int nb, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Hdiv base functions, Edge-based face functions by Ainsworth .
Definition: Hdiv.cpp:37
MoFEMErrorCode Hdiv_Ainsworth_FaceBubbleShapeFunctions_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int gdim, int nb, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face bubble functions by Ainsworth .
Definition: Hdiv.cpp:160
#define NBFACETRI_AINSWORTH_FACE_HDIV(P)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
#define NBFACETRI_AINSWORTH_HDIV(P)

◆ getValueHdivDemkowiczBase()

MoFEMErrorCode TriPolynomialBase::getValueHdivDemkowiczBase ( MatrixDouble pts)
private

Definition at line 366 of file TriPolynomialBase.cpp.

366  {
368 
370  // set shape functions into data structure
371  if (data.dataOnEntities[MBTRI].size() != 1) {
372  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
373  }
374  const FieldApproximationBase base = cTx->bAse;
375  if (base != DEMKOWICZ_JACOBI_BASE) {
376  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
377  "This should be used only with DEMKOWICZ_JACOBI_BASE "
378  "but base is %s",
379  ApproximationBaseNames[base]);
380  }
381  int order = data.dataOnEntities[MBTRI][0].getDataOrder();
382  int nb_gauss_pts = pts.size2();
383  data.dataOnEntities[MBTRI][0].getN(base).resize(
384  nb_gauss_pts, 3 * NBFACETRI_DEMKOWICZ_HDIV(order), false);
385  double *phi_f = &*data.dataOnEntities[MBTRI][0].getN(base).data().begin();
386  if (NBFACETRI_DEMKOWICZ_HDIV(order) == 0)
388  int face_nodes[3] = {0, 1, 2};
390  face_nodes, order,
391  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
392  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(), phi_f,
393  NULL, nb_gauss_pts, 3);
394 
396 }
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
MoFEMErrorCode Hdiv_Demkowicz_Face_MBTET_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int gdim, int nb)
Definition: Hdiv.cpp:617
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
static const char *const ApproximationBaseNames[]
Definition: definitions.h:158
EntPolynomialBaseCtx * cTx
#define NBFACETRI_DEMKOWICZ_HDIV(P)
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
#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

◆ getValueL2()

MoFEMErrorCode TriPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 265 of file TriPolynomialBase.cpp.

265  {
267 
269  const FieldApproximationBase base = cTx->bAse;
270  if (cTx->basePolynomialsType0 == NULL)
271  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
272  "Polynomial type not set");
273  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
274  double *diffL, const int dim) =
276 
277  int nb_gauss_pts = pts.size2();
278 
279  data.dataOnEntities[MBTRI][0].getN(base).resize(
280  nb_gauss_pts, NBFACETRI_L2(data.dataOnEntities[MBTRI][0].getDataOrder()),
281  false);
282  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(
283  nb_gauss_pts,
284  2 * NBFACETRI_L2(data.dataOnEntities[MBTRI][0].getDataOrder()), false);
285 
287  data.dataOnEntities[MBTRI][0].getDataOrder(),
288  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
289  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
290  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
291  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
292  nb_gauss_pts, base_polynomials);
293 
295 }
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
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
#define NBFACETRI_L2(P)
Number of base functions on triangle for L2 space.
#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
PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTRI(int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Get base functions on triangle for L2 space.
Definition: l2.c:29

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 26 of file TriPolynomialBase.cpp.

27  {
29  *iface = NULL;
30  if (uuid == IDD_TET_BASE_FUNCTION) {
31  *iface = const_cast<TriPolynomialBase *>(this);
33  } else {
34  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
35  }
38 }
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
static const MOFEMuuid IDD_TET_BASE_FUNCTION
#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

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* MoFEM::TriPolynomialBase::cTx
private

Definition at line 42 of file TriPolynomialBase.hpp.

◆ diffN_face_bubble

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

Definition at line 53 of file TriPolynomialBase.hpp.

◆ diffN_face_edge

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

Definition at line 52 of file TriPolynomialBase.hpp.

◆ N_face_bubble

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

Definition at line 51 of file TriPolynomialBase.hpp.

◆ N_face_edge

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

Definition at line 50 of file TriPolynomialBase.hpp.


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