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

Calculate base functions on tetrahedral. More...

#include <src/approximation/EdgePolynomialBase.hpp>

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

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, BaseFunctionUnknownInterface **iface) const
 
 EdgePolynomialBase ()
 
 ~EdgePolynomialBase ()
 
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 getValueHcurlAinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurlDemkowiczBase (MatrixDouble &pts)
 

Private Attributes

EntPolynomialBaseCtxcTx
 
VectorDouble L
 
VectorDouble diffL
 

Detailed Description

Calculate base functions on tetrahedral.

Definition at line 33 of file EdgePolynomialBase.hpp.

Constructor & Destructor Documentation

◆ EdgePolynomialBase()

EdgePolynomialBase::EdgePolynomialBase ( )

Definition at line 38 of file EdgePolynomialBase.cpp.

38 {}

◆ ~EdgePolynomialBase()

EdgePolynomialBase::~EdgePolynomialBase ( )

Definition at line 37 of file EdgePolynomialBase.cpp.

37 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 41 of file EdgePolynomialBase.cpp.

42  {
44 
46  CHKERR ctx_ptr->query_interface(IDD_EDGE_BASE_FUNCTION, &iface);
47  cTx = reinterpret_cast<EntPolynomialBaseCtx *>(iface);
48 
49  int nb_gauss_pts = pts.size2();
50  if (!nb_gauss_pts)
52 
53  if (pts.size1() < 1)
54  SETERRQ(
55  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
56  "Wrong dimension of pts, should be at least 3 rows with coordinates");
57 
58  const FieldApproximationBase base = cTx->bAse;
60 
61  if (base != AINSWORTH_BERNSTEIN_BEZIER_BASE) {
62  if (cTx->copyNodeBase == LASTBASE) {
63  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 2,
64  false);
66  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
67  &pts(0, 0), nb_gauss_pts);
68  } else
69  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
70  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
71 
72  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
73  (unsigned int)nb_gauss_pts)
74  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
75  "Base functions or nodes has wrong number of integration points "
76  "for base %s",
78 
79  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(2, 1, false);
80  std::copy(Tools::diffShapeFunMBEDGE.begin(),
82  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
83  }
84 
85  switch (cTx->sPace) {
86  case H1:
87  CHKERR getValueH1(pts);
88  break;
89  case HDIV:
90  CHKERR getValueHdiv(pts);
91  break;
92  case HCURL:
93  CHKERR getValueHcurl(pts);
94  break;
95  case L2:
96  CHKERR getValueL2(pts);
97  break;
98  default:
99  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
100  }
101 
103 }
field with continuous normal traction
Definition: definitions.h:173
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
MoFEMErrorCode getValueH1(MatrixDouble &pts)
static constexpr std::array< double, 2 > diffShapeFunMBEDGE
Definition: Tools.hpp:71
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEMErrorCode getValueL2(MatrixDouble &pts)
#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
PetscErrorCode ShapeMBEDGE(double *N, const double *G_X, int DIM)
Definition: fem_tools.c:773
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:172
static const MOFEMuuid IDD_EDGE_BASE_FUNCTION
MoFEMErrorCode getValueHdiv(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
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
continuous field
Definition: definitions.h:171
EntPolynomialBaseCtx * cTx
field with C-1 continuity
Definition: definitions.h:174

◆ getValueH1()

MoFEMErrorCode EdgePolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 105 of file EdgePolynomialBase.cpp.

105  {
107 
108  switch (cTx->bAse) {
112  break;
115  break;
116  default:
117  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
118  }
119 
121 }
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEMErrorCode getValueH1BernsteinBezierBase(MatrixDouble &pts)
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:146
#define CHKERR
Inline error check.
Definition: definitions.h:596
MoFEMErrorCode getValueH1AinsworthBase(MatrixDouble &pts)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
EntPolynomialBaseCtx * cTx

◆ getValueH1AinsworthBase()

MoFEMErrorCode EdgePolynomialBase::getValueH1AinsworthBase ( MatrixDouble pts)
private

Definition at line 219 of file EdgePolynomialBase.cpp.

219  {
221 
223  const FieldApproximationBase base = cTx->bAse;
224  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
225  double *diffL, const int dim) =
227 
228  int nb_gauss_pts = pts.size2();
229 
230  const int side_number = 0;
231  int sense = data.dataOnEntities[MBEDGE][side_number].getSense();
232  int order = data.dataOnEntities[MBEDGE][side_number].getDataOrder();
233  data.dataOnEntities[MBEDGE][side_number].getN(base).resize(
234  nb_gauss_pts, NBEDGE_H1(order), false);
235  data.dataOnEntities[MBEDGE][side_number].getDiffN(base).resize(
236  nb_gauss_pts, NBEDGE_H1(order), false);
237 
238  data.dataOnEntities[MBEDGE][side_number].getN(base).clear();
239  data.dataOnEntities[MBEDGE][side_number].getDiffN(base).clear();
240 
241  L.resize(NBEDGE_H1(order), false);
242  diffL.resize(NBEDGE_H1(order), false);
243 
244  if (data.dataOnEntities[MBEDGE][side_number].getDataOrder() > 1) {
245 
246  double diff_s = 2.; // s = s(xi), ds/dxi = 2., because change of basis
247  for (int gg = 0; gg != nb_gauss_pts; gg++) {
248 
249  double s = 2 * pts(0, gg) - 1; // makes form -1..1
250  if (!sense) {
251  s *= -1;
252  diff_s *= -1;
253  }
254 
255  // calculate Legendre polynomials at integration points
256  CHKERR base_polynomials(NBEDGE_H1(order) - 1, s, &diff_s,
257  &*L.data().begin(), &*diffL.data().begin(), 1);
258 
259  for (unsigned int pp = 0;
260  pp < data.dataOnEntities[MBEDGE][side_number].getN(base).size2();
261  pp++) {
262 
263  // Calculate edge shape functions N0*N1*L(p), where N0 and N1 are nodal
264  // shape functions
265  double v = data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 0) *
266  data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 1);
267  data.dataOnEntities[MBEDGE][side_number].getN(base)(gg, pp) = v * L(pp);
268 
269  // Calculate derivative edge shape functions
270  // dN/dksi = dN0/dxi*N1*L + N0*dN1/ksi*L + N0*N1*dL/dxi
271  data.dataOnEntities[MBEDGE][side_number].getDiffN(base)(gg, pp) =
272  ((+1.) * data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 1) +
273  data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 0) * (-1.)) *
274  L(pp) +
275  v * diffL(pp);
276  }
277  }
278  }
279 
281 }
#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
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
const FieldApproximationBase bAse
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
#define CHKERR
Inline error check.
Definition: definitions.h:596
EntPolynomialBaseCtx * cTx

◆ getValueH1BernsteinBezierBase()

MoFEMErrorCode EdgePolynomialBase::getValueH1BernsteinBezierBase ( MatrixDouble pts)
private

Definition at line 124 of file EdgePolynomialBase.cpp.

124  {
126 
128  const std::string &field_name = cTx->fieldName;
129  const int nb_gauss_pts = pts.size2();
130 
131  auto get_alpha = [field_name](auto &data) -> MatrixInt & {
132  auto &ptr = data.getBBAlphaIndicesSharedPtr(field_name);
133  if (!ptr)
134  ptr.reset(new MatrixInt());
135  return *ptr;
136  };
137 
138  auto get_base = [field_name](auto &data) -> MatrixDouble & {
139  auto &ptr = data.getBBNSharedPtr(field_name);
140  if (!ptr)
141  ptr.reset(new MatrixDouble());
142  return *ptr;
143  };
144 
145  auto get_diff_base = [field_name](auto &data) -> MatrixDouble & {
146  auto &ptr = data.getBBDiffNSharedPtr(field_name);
147  if (!ptr)
148  ptr.reset(new MatrixDouble());
149  return *ptr;
150  };
151 
152  auto &get_n = get_base(data.dataOnEntities[MBVERTEX][0]);
153  auto &get_diff_n = get_diff_base(data.dataOnEntities[MBVERTEX][0]);
154  get_n.resize(nb_gauss_pts, 2, false);
155  get_diff_n.resize(nb_gauss_pts, 2, false);
156  get_n.clear();
157  get_diff_n.clear();
158 
159  if (data.dataOnEntities[MBVERTEX][0].getN(NOBASE).size1() !=
160  (unsigned int)nb_gauss_pts)
161  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
162  "Base functions or nodes has wrong number of integration points "
163  "for base %s",
165  auto &lambda = data.dataOnEntities[MBVERTEX][0].getN(NOBASE);
166 
167  auto &vertex_alpha = get_alpha(data.dataOnEntities[MBVERTEX][0]);
168  vertex_alpha.resize(2, 2, false);
169  vertex_alpha(0, 0) = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[0];
170  vertex_alpha(0, 1) = 0;
171  vertex_alpha(1, 0) = 0;
172  vertex_alpha(1, 1) = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[1];
173 
175  1, nb_gauss_pts, vertex_alpha.size1(), &vertex_alpha(0, 0), &lambda(0, 0),
176  Tools::diffShapeFunMBEDGE.data(), &get_n(0, 0), &get_diff_n(0, 0));
177  std::array<double, 2> f = {
178  boost::math::factorial<double>(
179  data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[0]),
180  boost::math::factorial<double>(
181  data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[1])};
182 
183  for (int g = 0; g != nb_gauss_pts; ++g)
184  for (int n = 0; n != 2; ++n) {
185  get_n(g, n) *= f[n];
186  get_diff_n(g, n) *= f[n];
187  }
188 
189  if (data.spacesOnEntities[MBEDGE].test(H1)) {
190  if (data.dataOnEntities[MBEDGE].size() != 1)
191  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
192  "Wrong size ent of ent data");
193 
194  int order = data.dataOnEntities[MBEDGE][0].getDataOrder();
195  const int nb_dofs = NBEDGE_H1(order);
196 
197  auto &get_n = get_base(data.dataOnEntities[MBEDGE][0]);
198  auto &get_diff_n = get_diff_base(data.dataOnEntities[MBEDGE][0]);
199  get_n.resize(nb_gauss_pts, nb_dofs, false);
200  get_diff_n.resize(nb_gauss_pts, nb_dofs, false);
201 
202  if (nb_dofs) {
203  auto &edge_alpha = get_alpha(data.dataOnEntities[MBEDGE][0]);
204  edge_alpha.resize(nb_dofs, 2);
205  CHKERR BernsteinBezier::generateIndicesEdgeEdge(order, &edge_alpha(0, 0));
207  order, nb_gauss_pts, edge_alpha.size1(), &edge_alpha(0, 0),
208  &lambda(0, 0), Tools::diffShapeFunMBEDGE.data(), &get_n(0, 0),
209  &get_diff_n(0, 0));
210  }
211  } else {
212  get_n.resize(nb_gauss_pts, 0, false);
213  get_diff_n.resize(nb_gauss_pts, 0, false);
214  }
215 
217 }
#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 constexpr std::array< double, 2 > diffShapeFunMBEDGE
Definition: Tools.hpp:71
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
static MoFEMErrorCode baseFunctionsEdge(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)
ublas::matrix< int, ublas::row_major, IntAllocator > MatrixInt
Definition: Types.hpp:73
static const char *const ApproximationBaseNames[]
Definition: definitions.h:158
DataForcesAndSourcesCore & dAta
static MoFEMErrorCode generateIndicesEdgeEdge(const int N, int *alpha)
#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

◆ getValueHcurl()

MoFEMErrorCode EdgePolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 373 of file EdgePolynomialBase.cpp.

373  {
375 
376  switch (cTx->bAse) {
380  break;
383  break;
384  default:
385  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
386  }
387 
389 }
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)
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
EntPolynomialBaseCtx * cTx

◆ getValueHcurlAinsworthBase()

MoFEMErrorCode EdgePolynomialBase::getValueHcurlAinsworthBase ( MatrixDouble pts)
private

Definition at line 296 of file EdgePolynomialBase.cpp.

296  {
298 
300  const FieldApproximationBase base = cTx->bAse;
301 
302  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
303  double *diffL, const int dim) =
305 
306  int nb_gauss_pts = pts.size2();
307  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
308 
309  if (data.dataOnEntities[MBEDGE].size() != 1)
310  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
311 
312  int sense = data.dataOnEntities[MBEDGE][0].getSense();
313  int order = data.dataOnEntities[MBEDGE][0].getDataOrder();
314  int nb_dofs =
315  NBEDGE_AINSWORTH_HCURL(data.dataOnEntities[MBEDGE][0].getDataOrder());
316  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
317  false);
318  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
319  false);
320 
322  sense, order, &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0),
323  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
324  &*data.dataOnEntities[MBEDGE][0].getN(base).data().begin(), nullptr,
325  nb_gauss_pts, base_polynomials);
326 
327  } else {
328  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, false);
329  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
330  false);
331  }
332 
334 }
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
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 MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
EntPolynomialBaseCtx * cTx
MoFEMErrorCode Hcurl_Ainsworth_EdgeBaseFunctions_MBTET_ON_EDGE(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 edge.
Definition: Hcurl.cpp:175

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode EdgePolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 337 of file EdgePolynomialBase.cpp.

337  {
339 
341  const FieldApproximationBase base = cTx->bAse;
342 
343  int nb_gauss_pts = pts.size2();
344  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
345 
346  if (data.dataOnEntities[MBEDGE].size() != 1)
347  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
348  "No data structure to store base functions");
349 
350  int sense = data.dataOnEntities[MBEDGE][0].getSense();
351  int order = data.dataOnEntities[MBEDGE][0].getDataOrder();
352  int nb_dofs =
353  NBEDGE_DEMKOWICZ_HCURL(data.dataOnEntities[MBEDGE][0].getDataOrder());
354  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
355  false);
356  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
357  false);
359  sense, order, &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0),
360  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
361  &*data.dataOnEntities[MBEDGE][0].getN(base).data().begin(), NULL,
362  nb_gauss_pts);
363  } else {
364 
365  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, false);
366  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
367  false);
368  }
369 
371 }
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
FieldApproximationBase
approximation base
Definition: definitions.h:144
DataForcesAndSourcesCore & dAta
MoFEMErrorCode Hcurl_Demkowicz_EdgeBaseFunctions_MBEDGE(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 edge.
Definition: Hcurl.cpp:2141
field with continuous tangents
Definition: definitions.h:172
#define CHKERR
Inline error check.
Definition: definitions.h:596
#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
EntPolynomialBaseCtx * cTx

◆ getValueHdiv()

MoFEMErrorCode EdgePolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 290 of file EdgePolynomialBase.cpp.

290  {
291  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
292  "Make no sense, unless problem is 2d (2d not implemented yet)");
293 }

◆ getValueL2()

MoFEMErrorCode EdgePolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 283 of file EdgePolynomialBase.cpp.

283  {
285  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
286  "Make no sense, unless problem is 1d (1d not implemented yet)");
288 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 21 of file EdgePolynomialBase.cpp.

22  {
23 
25  *iface = NULL;
26  if (uuid == IDD_EDGE_BASE_FUNCTION) {
27  *iface = const_cast<EdgePolynomialBase *>(this);
29  } else {
30  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
31  }
32  ierr = BaseFunction::query_interface(uuid, iface);
33  CHKERRG(ierr);
35 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, MoFEM::BaseFunctionUnknownInterface **iface) const
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
static const MOFEMuuid IDD_EDGE_BASE_FUNCTION

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* MoFEM::EdgePolynomialBase::cTx
private

Definition at line 45 of file EdgePolynomialBase.hpp.

◆ diffL

VectorDouble MoFEM::EdgePolynomialBase::diffL
private

Definition at line 47 of file EdgePolynomialBase.hpp.

◆ L

VectorDouble MoFEM::EdgePolynomialBase::L
private

Definition at line 47 of file EdgePolynomialBase.hpp.


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