v0.13.1
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 (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 EdgePolynomialBase ()
 
 ~EdgePolynomialBase ()
 
MoFEMErrorCode getValue (MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::BaseFunction
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, MoFEM::UnknownInterface **iface) const
 
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 MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const =0
 
virtual ~BaseFunctionUnknownInterface ()=default
 
- Public Member Functions inherited from MoFEM::UnknownInterface
virtual MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const =0
 
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 

Private Member Functions

MoFEMErrorCode getValueH1 (MatrixDouble &pts)
 
MoFEMErrorCode getValueH1AinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueH1DemkowiczBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueH1BernsteinBezierBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2 (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2DemkowiczBase (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
 

Additional Inherited Members

- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version. More...
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version. More...
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version. More...
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version. More...
 

Detailed Description

Calculate base functions on tetrahedral.

Definition at line 19 of file EdgePolynomialBase.hpp.

Constructor & Destructor Documentation

◆ EdgePolynomialBase()

EdgePolynomialBase::EdgePolynomialBase ( )

Definition at line 17 of file EdgePolynomialBase.cpp.

17{}

◆ ~EdgePolynomialBase()

EdgePolynomialBase::~EdgePolynomialBase ( )

Definition at line 16 of file EdgePolynomialBase.cpp.

16{}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 20 of file EdgePolynomialBase.cpp.

21 {
23
25
26 int nb_gauss_pts = pts.size2();
27 if (!nb_gauss_pts)
29
30 if (pts.size1() < 1)
31 SETERRQ(
32 PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
33 "Wrong dimension of pts, should be at least 3 rows with coordinates");
34
35 const FieldApproximationBase base = cTx->bAse;
36 EntitiesFieldData &data = cTx->dAta;
37
39 if (cTx->copyNodeBase == LASTBASE) {
40 data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 2,
41 false);
43 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
44 &pts(0, 0), nb_gauss_pts);
45 } else
46 data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
47 data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
48
49 if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
50 (unsigned int)nb_gauss_pts)
51 SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
52 "Base functions or nodes has wrong number of integration points "
53 "for base %s",
55
56 data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(2, 1, false);
57 std::copy(Tools::diffShapeFunMBEDGE.begin(),
59 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
60 }
61
62 switch (cTx->sPace) {
63 case H1:
64 CHKERR getValueH1(pts);
65 break;
66 case HDIV:
68 break;
69 case HCURL:
71 break;
72 case L2:
73 CHKERR getValueL2(pts);
74 break;
75 default:
76 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
77 }
78
80}
FieldApproximationBase
approximation base
Definition: definitions.h:58
@ LASTBASE
Definition: definitions.h:69
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:64
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ H1
continuous field
Definition: definitions.h:85
@ HCURL
field with continuous tangents
Definition: definitions.h:86
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
static const char *const ApproximationBaseNames[]
Definition: definitions.h:72
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
PetscErrorCode ShapeMBEDGE(double *N, const double *G_X, int DIM)
Definition: fem_tools.c:761
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
MoFEMErrorCode getValueL2(MatrixDouble &pts)
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getValueH1(MatrixDouble &pts)
Class used to pass element data to calculate base functions on tet,triangle,edge.
const FieldApproximationBase copyNodeBase
const FieldApproximationBase bAse
data structure for finite element entity
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
static constexpr std::array< double, 2 > diffShapeFunMBEDGE
Definition: Tools.hpp:68
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ getValueH1()

MoFEMErrorCode EdgePolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 82 of file EdgePolynomialBase.cpp.

82 {
84
85 switch (cTx->bAse) {
89 break;
92 break;
95 break;
96 default:
97 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
98 }
99
101}
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:62
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:66
MoFEMErrorCode getValueH1AinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueH1DemkowiczBase(MatrixDouble &pts)
MoFEMErrorCode getValueH1BernsteinBezierBase(MatrixDouble &pts)

◆ getValueH1AinsworthBase()

MoFEMErrorCode EdgePolynomialBase::getValueH1AinsworthBase ( MatrixDouble pts)
private

Definition at line 199 of file EdgePolynomialBase.cpp.

199 {
201
202 EntitiesFieldData &data = cTx->dAta;
203 const FieldApproximationBase base = cTx->bAse;
204 PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
205 double *diffL, const int dim) =
207
208 int nb_gauss_pts = pts.size2();
209
210 const int side_number = 0;
211 int sense = data.dataOnEntities[MBEDGE][side_number].getSense();
212 int order = data.dataOnEntities[MBEDGE][side_number].getOrder();
213 data.dataOnEntities[MBEDGE][side_number].getN(base).resize(
214 nb_gauss_pts, NBEDGE_H1(order), false);
215 data.dataOnEntities[MBEDGE][side_number].getDiffN(base).resize(
216 nb_gauss_pts, NBEDGE_H1(order), false);
217
218 data.dataOnEntities[MBEDGE][side_number].getN(base).clear();
219 data.dataOnEntities[MBEDGE][side_number].getDiffN(base).clear();
220
221 L.resize(NBEDGE_H1(order), false);
222 diffL.resize(NBEDGE_H1(order), false);
223
224 if (data.dataOnEntities[MBEDGE][side_number].getOrder() > 1) {
225
226 double diff_s = 2.; // s = s(xi), ds/dxi = 2., because change of basis
227 for (int gg = 0; gg != nb_gauss_pts; gg++) {
228
229 double s = 2 * pts(0, gg) - 1; // makes form -1..1
230 if (!sense) {
231 s *= -1;
232 diff_s *= -1;
233 }
234
235 // calculate Legendre polynomials at integration points
236 CHKERR base_polynomials(NBEDGE_H1(order) - 1, s, &diff_s,
237 &*L.data().begin(), &*diffL.data().begin(), 1);
238
239 for (unsigned int pp = 0;
240 pp < data.dataOnEntities[MBEDGE][side_number].getN(base).size2();
241 pp++) {
242
243 // Calculate edge shape functions N0*N1*L(p), where N0 and N1 are nodal
244 // shape functions
245 double v = data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 0) *
246 data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 1);
247 data.dataOnEntities[MBEDGE][side_number].getN(base)(gg, pp) = v * L(pp);
248
249 // Calculate derivative edge shape functions
250 // dN/dksi = dN0/dxi*N1*L + N0*dN1/ksi*L + N0*N1*dL/dxi
251 data.dataOnEntities[MBEDGE][side_number].getDiffN(base)(gg, pp) =
252 ((+1.) * data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 1) +
253 data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 0) * (-1.)) *
254 L(pp) +
255 v * diffL(pp);
256 }
257 }
258 }
259
261}
static Index< 'p', 3 > p
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
const int dim
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
const double v
phase velocity of light in medium (cm/ns)
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)

◆ getValueH1BernsteinBezierBase()

MoFEMErrorCode EdgePolynomialBase::getValueH1BernsteinBezierBase ( MatrixDouble pts)
private

Definition at line 104 of file EdgePolynomialBase.cpp.

104 {
106
107 EntitiesFieldData &data = cTx->dAta;
108 const std::string &field_name = cTx->fieldName;
109 const int nb_gauss_pts = pts.size2();
110
111 auto get_alpha = [field_name](auto &data) -> MatrixInt & {
112 auto &ptr = data.getBBAlphaIndicesSharedPtr(field_name);
113 if (!ptr)
114 ptr.reset(new MatrixInt());
115 return *ptr;
116 };
117
118 auto get_base = [field_name](auto &data) -> MatrixDouble & {
119 auto &ptr = data.getBBNSharedPtr(field_name);
120 if (!ptr)
121 ptr.reset(new MatrixDouble());
122 return *ptr;
123 };
124
125 auto get_diff_base = [field_name](auto &data) -> MatrixDouble & {
126 auto &ptr = data.getBBDiffNSharedPtr(field_name);
127 if (!ptr)
128 ptr.reset(new MatrixDouble());
129 return *ptr;
130 };
131
132 auto &get_n = get_base(data.dataOnEntities[MBVERTEX][0]);
133 auto &get_diff_n = get_diff_base(data.dataOnEntities[MBVERTEX][0]);
134 get_n.resize(nb_gauss_pts, 2, false);
135 get_diff_n.resize(nb_gauss_pts, 2, false);
136 get_n.clear();
137 get_diff_n.clear();
138
139 if (data.dataOnEntities[MBVERTEX][0].getN(NOBASE).size1() !=
140 (unsigned int)nb_gauss_pts)
141 SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
142 "Base functions or nodes has wrong number of integration points "
143 "for base %s",
145 auto &lambda = data.dataOnEntities[MBVERTEX][0].getN(NOBASE);
146
147 auto &vertex_alpha = get_alpha(data.dataOnEntities[MBVERTEX][0]);
148 vertex_alpha.resize(2, 2, false);
149 vertex_alpha(0, 0) = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[0];
150 vertex_alpha(0, 1) = 0;
151 vertex_alpha(1, 0) = 0;
152 vertex_alpha(1, 1) = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[1];
153
155 1, nb_gauss_pts, vertex_alpha.size1(), &vertex_alpha(0, 0), &lambda(0, 0),
156 Tools::diffShapeFunMBEDGE.data(), &get_n(0, 0), &get_diff_n(0, 0));
157 std::array<double, 2> f = {
158 boost::math::factorial<double>(
159 data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[0]),
160 boost::math::factorial<double>(
161 data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[1])};
162
163 for (int g = 0; g != nb_gauss_pts; ++g)
164 for (int n = 0; n != 2; ++n) {
165 get_n(g, n) *= f[n];
166 get_diff_n(g, n) *= f[n];
167 }
168
169 if (data.spacesOnEntities[MBEDGE].test(H1)) {
170 if (data.dataOnEntities[MBEDGE].size() != 1)
171 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
172 "Wrong size ent of ent data");
173
174 int order = data.dataOnEntities[MBEDGE][0].getOrder();
175 const int nb_dofs = NBEDGE_H1(order);
176
177 auto &get_n = get_base(data.dataOnEntities[MBEDGE][0]);
178 auto &get_diff_n = get_diff_base(data.dataOnEntities[MBEDGE][0]);
179 get_n.resize(nb_gauss_pts, nb_dofs, false);
180 get_diff_n.resize(nb_gauss_pts, nb_dofs, false);
181
182 if (nb_dofs) {
183 auto &edge_alpha = get_alpha(data.dataOnEntities[MBEDGE][0]);
184 edge_alpha.resize(nb_dofs, 2);
187 order, nb_gauss_pts, edge_alpha.size1(), &edge_alpha(0, 0),
188 &lambda(0, 0), Tools::diffShapeFunMBEDGE.data(), &get_n(0, 0),
189 &get_diff_n(0, 0));
190 }
191 } else {
192 get_n.resize(nb_gauss_pts, 0, false);
193 get_diff_n.resize(nb_gauss_pts, 0, false);
194 }
195
197}
@ NOBASE
Definition: definitions.h:59
FTensor::Index< 'n', SPACE_DIM > n
constexpr double lambda
auto f
Definition: HenckyOps.hpp:5
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
UBlasMatrix< int > MatrixInt
Definition: Types.hpp:76
constexpr auto field_name
constexpr double g
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)
static MoFEMErrorCode generateIndicesEdgeEdge(const int N, int *alpha)
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types

◆ getValueH1DemkowiczBase()

MoFEMErrorCode EdgePolynomialBase::getValueH1DemkowiczBase ( MatrixDouble pts)
private

Definition at line 263 of file EdgePolynomialBase.cpp.

263 {
265
266 EntitiesFieldData &data = cTx->dAta;
267 const FieldApproximationBase base = cTx->bAse;
268
269 int nb_gauss_pts = pts.size2();
270
271 const int side_number = 0;
272 int order = data.dataOnEntities[MBEDGE][side_number].getOrder();
273
274 data.dataOnEntities[MBEDGE][side_number].getN(base).resize(
275 nb_gauss_pts, NBEDGE_H1(order), false);
276 data.dataOnEntities[MBEDGE][side_number].getDiffN(base).resize(
277 nb_gauss_pts, NBEDGE_H1(order), false);
278
279 data.dataOnEntities[MBEDGE][side_number].getN(base).clear();
280 data.dataOnEntities[MBEDGE][side_number].getDiffN(base).clear();
281
282 double diff_shape[] = {-1, 1};
283 MatrixDouble shape(nb_gauss_pts, 2);
284
285 if (data.dataOnEntities[MBEDGE][side_number].getOrder() > 1) {
286 for (int gg = 0; gg != nb_gauss_pts; gg++) {
287 const double mu0 = 1.0 - pts(0, gg); // pts ranges over [0, 1]
288 const double mu1 = pts(0, gg);
289 shape(gg, 0) = mu0;
290 shape(gg, 1) = mu1;
291 }
292
294 order, &*shape.data().begin(), diff_shape,
295 &*data.dataOnEntities[MBEDGE][side_number].getN(base).data().begin(),
296 &*data.dataOnEntities[MBEDGE][side_number]
297 .getDiffN(base)
298 .data()
299 .begin(),
300 nb_gauss_pts);
301 }
302
304}
MoFEMErrorCode H1_BubbleShapeFunctions_ONSEGMENT(int p, double *N, double *diffN, double *bubbleN, double *diff_bubbleN, int nb_integration_pts)

◆ getValueHcurl()

MoFEMErrorCode EdgePolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 444 of file EdgePolynomialBase.cpp.

444 {
446
447 switch (cTx->bAse) {
451 break;
454 break;
455 default:
456 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
457 }
458
460}
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)

◆ getValueHcurlAinsworthBase()

MoFEMErrorCode EdgePolynomialBase::getValueHcurlAinsworthBase ( MatrixDouble pts)
private

Definition at line 367 of file EdgePolynomialBase.cpp.

367 {
369
370 EntitiesFieldData &data = cTx->dAta;
371 const FieldApproximationBase base = cTx->bAse;
372
373 PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
374 double *diffL, const int dim) =
376
377 int nb_gauss_pts = pts.size2();
378 if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
379
380 if (data.dataOnEntities[MBEDGE].size() != 1)
381 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
382
383 int sense = data.dataOnEntities[MBEDGE][0].getSense();
384 int order = data.dataOnEntities[MBEDGE][0].getOrder();
385 int nb_dofs =
386 NBEDGE_AINSWORTH_HCURL(data.dataOnEntities[MBEDGE][0].getOrder());
387 data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
388 false);
389 data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
390 false);
391
393 sense, order, &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0),
394 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
395 &*data.dataOnEntities[MBEDGE][0].getN(base).data().begin(), nullptr,
396 nb_gauss_pts, base_polynomials);
397
398 } else {
399 data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, false);
400 data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
401 false);
402 }
403
405}
#define NBEDGE_AINSWORTH_HCURL(P)
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 408 of file EdgePolynomialBase.cpp.

408 {
410
411 EntitiesFieldData &data = cTx->dAta;
412 const FieldApproximationBase base = cTx->bAse;
413
414 int nb_gauss_pts = pts.size2();
415 if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
416
417 if (data.dataOnEntities[MBEDGE].size() != 1)
418 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
419 "No data structure to store base functions");
420
421 int sense = data.dataOnEntities[MBEDGE][0].getSense();
422 int order = data.dataOnEntities[MBEDGE][0].getOrder();
423 int nb_dofs =
424 NBEDGE_DEMKOWICZ_HCURL(data.dataOnEntities[MBEDGE][0].getOrder());
425 data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
426 false);
427 data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
428 false);
430 sense, order, &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0),
431 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
432 &*data.dataOnEntities[MBEDGE][0].getN(base).data().begin(), NULL,
433 nb_gauss_pts);
434 } else {
435
436 data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, false);
437 data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
438 false);
439 }
440
442}
#define NBEDGE_DEMKOWICZ_HCURL(P)
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:2144

◆ getValueHdiv()

MoFEMErrorCode EdgePolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 361 of file EdgePolynomialBase.cpp.

361 {
362 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
363 "Make no sense, unless problem is 2d (2d not implemented yet)");
364}

◆ getValueL2()

MoFEMErrorCode EdgePolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 306 of file EdgePolynomialBase.cpp.

306 {
308
309 switch (cTx->bAse) {
312 break;
313 default:
314 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
315 }
316
318}
MoFEMErrorCode getValueL2DemkowiczBase(MatrixDouble &pts)

◆ getValueL2DemkowiczBase()

MoFEMErrorCode EdgePolynomialBase::getValueL2DemkowiczBase ( MatrixDouble pts)
private

Definition at line 320 of file EdgePolynomialBase.cpp.

320 {
322
323 EntitiesFieldData &data = cTx->dAta;
324 const FieldApproximationBase base = cTx->bAse;
325
326 int nb_gauss_pts = pts.size2();
327
328 const int side_number = 0;
329 int order = data.dataOnEntities[MBEDGE][side_number].getOrder();
330
331 data.dataOnEntities[MBEDGE][side_number].getN(base).resize(
332 nb_gauss_pts, NBEDGE_L2(order), false);
333 data.dataOnEntities[MBEDGE][side_number].getDiffN(base).resize(
334 nb_gauss_pts, NBEDGE_L2(order), false);
335
336 data.dataOnEntities[MBEDGE][side_number].getN(base).clear();
337
338 double diff_n[] = {-1, 1};
339 MatrixDouble shape(nb_gauss_pts, 2);
340 if (data.dataOnEntities[MBEDGE][side_number].getOrder() > 1) {
341 for (int gg = 0; gg != nb_gauss_pts; gg++) {
342 const double mu0 = 1.0 - pts(0, gg); // pts ranges over [0, 1]
343 const double mu1 = pts(0, gg);
344 shape(gg, 0) = mu0;
345 shape(gg, 1) = mu1;
346 }
347
349 order, &*shape.data().begin(), diff_n,
350 &*data.dataOnEntities[MBEDGE][side_number].getN(base).data().begin(),
351 &*data.dataOnEntities[MBEDGE][side_number]
352 .getDiffN(base)
353 .data()
354 .begin(),
355 nb_gauss_pts);
356 }
357
359}
#define NBEDGE_L2(P)
Number of base functions on edge fro L2 space.
MoFEMErrorCode L2_ShapeFunctions_ONSEGMENT(int p, double *N, double *diffN, double *funN, double *funDiffN, int nb_integration_pts)

◆ query_interface()

MoFEMErrorCode EdgePolynomialBase::query_interface ( boost::typeindex::type_index  type_index,
UnknownInterface **  iface 
) const
virtual

Reimplemented from MoFEM::BaseFunction.

Definition at line 10 of file EdgePolynomialBase.cpp.

11 {
12 *iface = const_cast<EdgePolynomialBase *>(this);
13 return 0;
14}
Calculate base functions on tetrahedral.

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* MoFEM::EdgePolynomialBase::cTx
private

Definition at line 31 of file EdgePolynomialBase.hpp.

◆ diffL

VectorDouble MoFEM::EdgePolynomialBase::diffL
private

Definition at line 33 of file EdgePolynomialBase.hpp.

◆ L

VectorDouble MoFEM::EdgePolynomialBase::L
private

Definition at line 33 of file EdgePolynomialBase.hpp.


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