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

Static Public Member Functions

static MoFEMErrorCode setDofsSideMap (const FieldSpace space, const FieldContinuity continuity, const FieldApproximationBase base, DofsSideMap &)
 Set map of dof to side number.
 
- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version.
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version.
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version.
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version.
 

Private Member Functions

MoFEMErrorCode getValueH1 (MatrixDouble &pts)
 
MoFEMErrorCode getValueH1AinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueH1BernsteinBezierBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2 (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2AinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2BernsteinBezierBase (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)
 
MoFEMErrorCode getValueHcurlAinsworthBrokenBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurlDemkowiczBrokenBase (MatrixDouble &pts)
 

Static Private Member Functions

static MoFEMErrorCode setDofsSideMapHcurl (const FieldSpace space, const FieldContinuity continuity, const FieldApproximationBase base, DofsSideMap &dofs_side_map)
 

Private Attributes

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

Additional Inherited Members

- Public Types inherited from MoFEM::BaseFunction
using DofsSideMap = multi_index_container< DofsSideMapData, indexed_by< ordered_non_unique< tag< TypeSide_mi_tag >, composite_key< DofsSideMapData, member< DofsSideMapData, EntityType, &DofsSideMapData::type >, member< DofsSideMapData, int, &DofsSideMapData::side > > >, ordered_unique< tag< EntDofIdx_mi_tag >, member< DofsSideMapData, int, &DofsSideMapData::dof > > > >
 Map entity stype and side to element/entity dof index.
 

Detailed Description

Calculate base functions on triangle.

Definition at line 16 of file TriPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ TriPolynomialBase()

MoFEM::TriPolynomialBase::TriPolynomialBase ( )
default

◆ ~TriPolynomialBase()

virtual MoFEM::TriPolynomialBase::~TriPolynomialBase ( )
virtualdefault

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 1127 of file TriPolynomialBase.cpp.

1128 {
1130
1132
1133 int nb_gauss_pts = pts.size2();
1134 if (!nb_gauss_pts) {
1136 }
1137
1138 if (pts.size1() < 2)
1139 SETERRQ(
1140 PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1141 "Wrong dimension of pts, should be at least 3 rows with coordinates");
1142
1143 const FieldApproximationBase base = cTx->bAse;
1144 EntitiesFieldData &data = cTx->dAta;
1145
1146 if (base != AINSWORTH_BERNSTEIN_BEZIER_BASE) {
1147 if (cTx->copyNodeBase == LASTBASE) {
1148 data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 3,
1149 false);
1151 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
1152 &pts(0, 0), &pts(1, 0), nb_gauss_pts);
1153 data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(3, 2, false);
1155 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
1156 } else {
1157 data.dataOnEntities[MBVERTEX][0].getN(base) =
1158 data.dataOnEntities[MBVERTEX][0].getN(cTx->copyNodeBase);
1159 data.dataOnEntities[MBVERTEX][0].getDiffN(base) =
1160 data.dataOnEntities[MBVERTEX][0].getDiffN(cTx->copyNodeBase);
1161 }
1162 if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
1163 (unsigned int)nb_gauss_pts) {
1164 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1165 "Base functions or nodes has wrong number of integration points "
1166 "for base %s",
1168 }
1169 }
1170
1171 auto set_node_derivative_for_all_gauss_pts = [&]() {
1173 // In linear geometry derivatives are constant,
1174 // this in expense of efficiency makes implementation
1175 // consistent between vertices and other types of entities
1176 data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 6,
1177 false);
1178 for (int gg = 0; gg != nb_gauss_pts; ++gg)
1179 std::copy(Tools::diffShapeFunMBTRI.begin(),
1181 &data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg, 0));
1183 };
1184
1185 CHKERR set_node_derivative_for_all_gauss_pts();
1186
1187 switch (cTx->spaceContinuity) {
1188 case CONTINUOUS:
1189
1190 switch (cTx->sPace) {
1191 case H1:
1193 case HDIV:
1195 case HCURL:
1197 case L2:
1199 default:
1200 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
1201 }
1202
1203 case DISCONTINUOUS:
1204 switch (cTx->sPace) {
1205 case HCURL:
1207 default:
1208 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
1209 }
1210
1211 default:
1212 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Unknown continuity");
1213 }
1214
1216}
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()
@ 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
@ CONTINUOUS
Regular field.
@ DISCONTINUOUS
Broken continuity (No effect on L2 space)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ 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()
#define CHKERR
Inline error check.
PetscErrorCode ShapeDiffMBTRI(double *diffN)
calculate derivatives of shape functions
Definition fem_tools.c:194
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition fem_tools.c:182
Class used to pass element data to calculate base functions on tet,triangle,edge.
const FieldApproximationBase copyNodeBase
const FieldContinuity spaceContinuity
const FieldApproximationBase bAse
data structure for finite element entity
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition Tools.hpp:104
MoFEMErrorCode getValueH1(MatrixDouble &pts)
MoFEMErrorCode getValueL2(MatrixDouble &pts)
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.

◆ getValueH1()

MoFEMErrorCode TriPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 16 of file TriPolynomialBase.cpp.

16 {
18
19 switch (cTx->bAse) {
25 default:
26 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
27 }
28
30}
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition definitions.h:60
@ AINSWORTH_LOBATTO_BASE
Definition definitions.h:62
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
MoFEMErrorCode getValueH1BernsteinBezierBase(MatrixDouble &pts)
MoFEMErrorCode getValueH1AinsworthBase(MatrixDouble &pts)

◆ getValueH1AinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueH1AinsworthBase ( MatrixDouble pts)
private

Definition at line 246 of file TriPolynomialBase.cpp.

246 {
248
249 EntitiesFieldData &data = cTx->dAta;
250 const FieldApproximationBase base = cTx->bAse;
251
252 if (cTx->basePolynomialsType0 == NULL)
253 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
254 "Polynomial type not set");
255
256 PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
257 double *diffL, const int dim) =
259
260 int nb_gauss_pts = pts.size2();
261
262 if (data.spacesOnEntities[MBEDGE].test(H1)) {
263 // edges
264 if (data.dataOnEntities[MBEDGE].size() != 3)
265 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
266 "expected size of data.dataOnEntities[MBEDGE] is 3");
267
268 int sense[3], order[3];
269 double *H1edgeN[3], *diffH1edgeN[3];
270 for (int ee = 0; ee < 3; ee++) {
271
272 if (data.dataOnEntities[MBEDGE][ee].getSense() == 0)
273 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
274 "sense of the edge unknown");
275
276 sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
277 order[ee] = data.dataOnEntities[MBEDGE][ee].getOrder();
278 int nb_dofs = NBEDGE_H1(data.dataOnEntities[MBEDGE][ee].getOrder());
279 data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
280 false);
281 data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,
282 2 * nb_dofs, false);
283 H1edgeN[ee] = &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
284 diffH1edgeN[ee] =
285 &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
286 }
288 sense, order,
289 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
290 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
291 H1edgeN, diffH1edgeN, nb_gauss_pts, base_polynomials);
292 }
293
294 if (data.spacesOnEntities[MBTRI].test(H1)) {
295 // face
296 if (data.dataOnEntities[MBTRI].size() != 1) {
297 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
298 "expected that size data.dataOnEntities[MBTRI] is one");
299 }
300
301 int nb_dofs = NBFACETRI_H1(data.dataOnEntities[MBTRI][0].getOrder());
302 data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, nb_dofs,
303 false);
304 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
305 2 * nb_dofs, false);
306 const int face_nodes[] = {0, 1, 2};
308 face_nodes, data.dataOnEntities[MBTRI][0].getOrder(),
309 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
310 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
311 &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
312 &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
313 nb_gauss_pts, base_polynomials);
314 }
315
317}
constexpr int order
#define NBEDGE_H1(P)
Number of base function on edge for H1 space.
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:17
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:191
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types

◆ getValueH1BernsteinBezierBase()

MoFEMErrorCode TriPolynomialBase::getValueH1BernsteinBezierBase ( MatrixDouble pts)
private

Definition at line 33 of file TriPolynomialBase.cpp.

33 {
35 EntitiesFieldData &data = cTx->dAta;
36 const std::string &field_name = cTx->fieldName;
37 int nb_gauss_pts = pts.size2();
38
39 auto get_alpha = [field_name](auto &data) -> MatrixInt & {
40 auto &ptr = data.getBBAlphaIndicesSharedPtr(field_name);
41 if (!ptr)
42 ptr.reset(new MatrixInt());
43 return *ptr;
44 };
45
46 auto get_base = [field_name](auto &data) -> MatrixDouble & {
47 auto &ptr = data.getBBNSharedPtr(field_name);
48 if (!ptr)
49 ptr.reset(new MatrixDouble());
50 return *ptr;
51 };
52
53 auto get_diff_base = [field_name](auto &data) -> MatrixDouble & {
54 auto &ptr = data.getBBDiffNSharedPtr(field_name);
55 if (!ptr)
56 ptr.reset(new MatrixDouble());
57 return *ptr;
58 };
59
60 auto get_alpha_by_name_ptr =
61 [](auto &data,
62 const std::string &field_name) -> boost::shared_ptr<MatrixInt> & {
63 return data.getBBAlphaIndicesSharedPtr(field_name);
64 };
65
66 auto get_base_by_name_ptr =
67 [](auto &data,
68 const std::string &field_name) -> boost::shared_ptr<MatrixDouble> & {
69 return data.getBBNSharedPtr(field_name);
70 };
71
72 auto get_diff_base_by_name_ptr =
73 [](auto &data,
74 const std::string &field_name) -> boost::shared_ptr<MatrixDouble> & {
75 return data.getBBDiffNSharedPtr(field_name);
76 };
77
78 auto get_alpha_by_order_ptr =
79 [](auto &data, const size_t o) -> boost::shared_ptr<MatrixInt> & {
80 return data.getBBAlphaIndicesByOrderSharedPtr(o);
81 };
82
83 auto get_base_by_order_ptr =
84 [](auto &data, const size_t o) -> boost::shared_ptr<MatrixDouble> & {
85 return data.getBBNByOrderSharedPtr(o);
86 };
87
88 auto get_diff_base_by_order_ptr =
89 [](auto &data, const size_t o) -> boost::shared_ptr<MatrixDouble> & {
90 return data.getBBDiffNByOrderSharedPtr(o);
91 };
92
93 auto &vert_get_n = get_base(data.dataOnEntities[MBVERTEX][0]);
94 auto &vert_get_diff_n = get_diff_base(data.dataOnEntities[MBVERTEX][0]);
95 vert_get_n.resize(nb_gauss_pts, 3, false);
96 vert_get_diff_n.resize(nb_gauss_pts, 6, false);
97
98 if (data.dataOnEntities[MBVERTEX][0].getN(NOBASE).size1() !=
99 (unsigned int)nb_gauss_pts)
100 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
101 "Base functions or nodes has wrong number of integration points "
102 "for base %s",
104 auto &lambda = data.dataOnEntities[MBVERTEX][0].getN(NOBASE);
105
106 auto &vertex_alpha = get_alpha(data.dataOnEntities[MBVERTEX][0]);
107 vertex_alpha.resize(3, 3, false);
108 vertex_alpha.clear();
109 for (int n = 0; n != 3; ++n)
110 vertex_alpha(n, n) = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[n];
111
113 1, lambda.size1(), vertex_alpha.size1(), &vertex_alpha(0, 0),
114 &lambda(0, 0), Tools::diffShapeFunMBTRI.data(), &vert_get_n(0, 0),
115 &vert_get_diff_n(0, 0));
116 for (int n = 0; n != 3; ++n) {
117 const int f = boost::math::factorial<double>(
118 data.dataOnEntities[MBVERTEX][0].getBBNodeOrder()[n]);
119 for (int g = 0; g != nb_gauss_pts; ++g) {
120 vert_get_n(g, n) *= f;
121 for (int d = 0; d != 2; ++d)
122 vert_get_diff_n(g, 2 * n + d) *= f;
123 }
124 }
125
126 // edges
127 if (data.spacesOnEntities[MBEDGE].test(H1)) {
128 if (data.dataOnEntities[MBEDGE].size() != 3)
129 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
130 "Wrong size of ent data");
131
132 constexpr int edges_nodes[3][2] = {{0, 1}, {1, 2}, {2, 0}};
133 for (int ee = 0; ee != 3; ++ee) {
134 auto &ent_data = data.dataOnEntities[MBEDGE][ee];
135
136 if (ent_data.getSense() == 0)
137 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
138 "Sense of the edge unknown");
139
140 const int sense = ent_data.getSense();
141 const int order = ent_data.getOrder();
142 const int nb_dofs = NBEDGE_H1(ent_data.getOrder());
143
144 if (nb_dofs) {
145 if (get_alpha_by_order_ptr(ent_data, order)) {
146 get_alpha_by_name_ptr(ent_data, field_name) =
147 get_alpha_by_order_ptr(ent_data, order);
148 get_base_by_name_ptr(ent_data, field_name) =
149 get_base_by_order_ptr(ent_data, order);
150 get_diff_base_by_name_ptr(ent_data, field_name) =
151 get_diff_base_by_order_ptr(ent_data, order);
152 } else {
153 auto &get_n = get_base(ent_data);
154 auto &get_diff_n = get_diff_base(ent_data);
155 get_n.resize(nb_gauss_pts, nb_dofs, false);
156 get_diff_n.resize(nb_gauss_pts, 2 * nb_dofs, false);
157
158 auto &edge_alpha = get_alpha(data.dataOnEntities[MBEDGE][ee]);
159 edge_alpha.resize(nb_dofs, 3, false);
161 &edge_alpha(0, 0));
162 if (sense == -1) {
163 for (int i = 0; i != edge_alpha.size1(); ++i) {
164 int a = edge_alpha(i, edges_nodes[ee][0]);
165 edge_alpha(i, edges_nodes[ee][0]) =
166 edge_alpha(i, edges_nodes[ee][1]);
167 edge_alpha(i, edges_nodes[ee][1]) = a;
168 }
169 }
171 order, lambda.size1(), edge_alpha.size1(), &edge_alpha(0, 0),
172 &lambda(0, 0), Tools::diffShapeFunMBTRI.data(), &get_n(0, 0),
173 &get_diff_n(0, 0));
174
175 get_alpha_by_order_ptr(ent_data, order) =
176 get_alpha_by_name_ptr(ent_data, field_name);
177 get_base_by_order_ptr(ent_data, order) =
178 get_base_by_name_ptr(ent_data, field_name);
179 get_diff_base_by_order_ptr(ent_data, order) =
180 get_diff_base_by_name_ptr(ent_data, field_name);
181 }
182 }
183 }
184 } else {
185 for (int ee = 0; ee != 3; ++ee) {
186 auto &ent_data = data.dataOnEntities[MBEDGE][ee];
187 ent_data.getBBAlphaIndicesSharedPtr(field_name).reset();
188 auto &get_n = get_base(ent_data);
189 auto &get_diff_n = get_diff_base(ent_data);
190 get_n.resize(nb_gauss_pts, 0, false);
191 get_diff_n.resize(nb_gauss_pts, 0, false);
192 }
193 }
194
195 if (data.spacesOnEntities[MBTRI].test(H1)) {
196 if (data.dataOnEntities[MBTRI].size() != 1)
197 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
198 "Wrong size ent of ent data");
199
200 auto &ent_data = data.dataOnEntities[MBTRI][0];
201 const int order = ent_data.getOrder();
202 const int nb_dofs = NBFACETRI_H1(order);
203 if (get_alpha_by_order_ptr(ent_data, order)) {
204 get_alpha_by_name_ptr(ent_data, field_name) =
205 get_alpha_by_order_ptr(ent_data, order);
206 get_base_by_name_ptr(ent_data, field_name) =
207 get_base_by_order_ptr(ent_data, order);
208 get_diff_base_by_name_ptr(ent_data, field_name) =
209 get_diff_base_by_order_ptr(ent_data, order);
210 } else {
211
212 auto &get_n = get_base(ent_data);
213 auto &get_diff_n = get_diff_base(ent_data);
214 get_n.resize(nb_gauss_pts, nb_dofs, false);
215 get_diff_n.resize(nb_gauss_pts, 2 * nb_dofs, false);
216 if (nb_dofs) {
217 auto &tri_alpha = get_alpha(ent_data);
218 tri_alpha.resize(nb_dofs, 3, false);
219
222 order, lambda.size1(), tri_alpha.size1(), &tri_alpha(0, 0),
223 &lambda(0, 0), Tools::diffShapeFunMBTRI.data(), &get_n(0, 0),
224 &get_diff_n(0, 0));
225
226 get_alpha_by_order_ptr(ent_data, order) =
227 get_alpha_by_name_ptr(ent_data, field_name);
228 get_base_by_order_ptr(ent_data, order) =
229 get_base_by_name_ptr(ent_data, field_name);
230 get_diff_base_by_order_ptr(ent_data, order) =
231 get_diff_base_by_name_ptr(ent_data, field_name);
232 }
233 }
234 } else {
235 auto &ent_data = data.dataOnEntities[MBTRI][0];
236 ent_data.getBBAlphaIndicesSharedPtr(field_name).reset();
237 auto &get_n = get_base(ent_data);
238 auto &get_diff_n = get_diff_base(ent_data);
239 get_n.resize(nb_gauss_pts, 0, false);
240 get_diff_n.resize(nb_gauss_pts, 0, false);
241 }
242
244}
constexpr double a
@ NOBASE
Definition definitions.h:59
FTensor::Index< 'i', SPACE_DIM > i
static double lambda
const double n
refractive index of diffusive medium
UBlasMatrix< double > MatrixDouble
Definition Types.hpp:77
UBlasMatrix< int > MatrixInt
Definition Types.hpp:76
constexpr auto field_name
constexpr double g
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 MoFEMErrorCode generateIndicesTriTri(const int N, int *alpha)
static MoFEMErrorCode generateIndicesEdgeTri(const int N[], int *alpha[])

◆ getValueHcurl()

MoFEMErrorCode TriPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 1091 of file TriPolynomialBase.cpp.

1091 {
1093 switch (cTx->spaceContinuity) {
1094 case CONTINUOUS:
1095
1096 switch (cTx->bAse) {
1102 default:
1103 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
1104 }
1105
1106 break;
1107
1108 case DISCONTINUOUS:
1109 switch (cTx->bAse) {
1115 default:
1116 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
1117 }
1118
1119 default:
1120 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Unknown continuity");
1121 }
1122
1124}
@ DEMKOWICZ_JACOBI_BASE
Definition definitions.h:66
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlAinsworthBrokenBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlDemkowiczBrokenBase(MatrixDouble &pts)

◆ getValueHcurlAinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlAinsworthBase ( MatrixDouble pts)
private

Definition at line 663 of file TriPolynomialBase.cpp.

663 {
665
666 EntitiesFieldData &data = cTx->dAta;
667 const FieldApproximationBase base = cTx->bAse;
668 if (data.dataOnEntities[MBTRI].size() != 1)
669 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
670
671 PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
672 double *diffL, const int dim) =
674
675 int nb_gauss_pts = pts.size2();
676
677 // Calculation H-curl on triangle faces
678 if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
679 if (data.dataOnEntities[MBEDGE].size() != 3) {
680 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
681 }
682 int sense[3], order[3];
683 double *hcurl_edge_n[3];
684 double *diff_hcurl_edge_n[3];
685 for (int ee = 0; ee < 3; ee++) {
686 if (data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
687 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
688 "data inconsistency");
689 }
690 sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
691 order[ee] = data.dataOnEntities[MBEDGE][ee].getOrder();
692 int nb_dofs =
693 NBEDGE_AINSWORTH_HCURL(data.dataOnEntities[MBEDGE][ee].getOrder());
694 data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
695 3 * nb_dofs, false);
696 data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
697 nb_gauss_pts, 2 * 3 * nb_dofs, false);
698 hcurl_edge_n[ee] =
699 &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
700 diff_hcurl_edge_n[ee] =
701 &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
702 }
704 sense, order,
705 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
706 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
707 hcurl_edge_n, diff_hcurl_edge_n, nb_gauss_pts, base_polynomials);
708 } else {
709 for (int ee = 0; ee < 3; ee++) {
710 data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, 0, false);
711 data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts, 0,
712 false);
713 }
714 }
715
716 if (data.spacesOnEntities[MBTRI].test(HCURL)) {
717
718 // cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << endl;
719 // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
720 //
721 // face
722 if (data.dataOnEntities[MBTRI].size() != 1) {
723 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
724 }
725 int order = data.dataOnEntities[MBTRI][0].getOrder();
726 int nb_dofs = NBFACETRI_AINSWORTH_HCURL(order);
727 data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
728 false);
729 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
730 3 * 2 * nb_dofs, false);
731 // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
732 int face_nodes[] = {0, 1, 2};
734 face_nodes, order,
735 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
736 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
737 &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
738 &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
739 nb_gauss_pts, base_polynomials);
740 // cerr << data.dataOnEntities[MBTRI][0].getN(base) << endl;
741 } else {
742 data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
743 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
744 }
745
747}
#define NBFACETRI_AINSWORTH_HCURL(P)
#define NBEDGE_AINSWORTH_HCURL(P)
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
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

◆ getValueHcurlAinsworthBrokenBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlAinsworthBrokenBase ( MatrixDouble pts)
private

Definition at line 750 of file TriPolynomialBase.cpp.

750 {
752
753 EntitiesFieldData &data = cTx->dAta;
754 const FieldApproximationBase base = cTx->bAse;
755 if (data.dataOnEntities[MBTRI].size() != 1)
756 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
757
758 PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
759 double *diffL, const int dim) =
761
762 int nb_gauss_pts = pts.size2();
763
764 if (data.spacesOnEntities[MBTRI].test(HCURL)) {
765
766 // face
767 if (data.dataOnEntities[MBTRI].size() != 1) {
768 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
769 }
770
771 int order = data.dataOnEntities[MBTRI][0].getOrder();
772 auto nb_edge_dofs = NBEDGE_AINSWORTH_HCURL(order);
773 int nb_dofs_face = NBFACETRI_AINSWORTH_HCURL(order);
774
775 // three faces on triangle
776 int nb_dofs = 3 * nb_edge_dofs + nb_dofs_face;
777 data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
778 false);
779 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
780 3 * 2 * nb_dofs, false);
781
782 MatrixDouble edge_bases[] = {
783 MatrixDouble(nb_gauss_pts, 3 * nb_edge_dofs, false),
784 MatrixDouble(nb_gauss_pts, 3 * nb_edge_dofs, false),
785 MatrixDouble(nb_gauss_pts, 3 * nb_edge_dofs, false)};
786 MatrixDouble diff_edge_bases[] = {
787 MatrixDouble(nb_gauss_pts, 3 * 2 * nb_edge_dofs, false),
788 MatrixDouble(nb_gauss_pts, 3 * 2 * nb_edge_dofs, false),
789 MatrixDouble(nb_gauss_pts, 3 * 2 * nb_edge_dofs, false)};
790
791 std::array<double *, 3> phi{&*edge_bases[0].data().begin(),
792 &*edge_bases[1].data().begin(),
793 &*edge_bases[2].data().begin()};
794 std::array<double *, 3> diff_phi{&*diff_edge_bases[0].data().begin(),
795 &*diff_edge_bases[1].data().begin(),
796 &*diff_edge_bases[2].data().begin()};
797
798 std::array<int, 3> edge_order{order, order, order};
799 std::array<int, 3> edge_sense{1, 1, 1};
801 edge_sense.data(), edge_order.data(),
802 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
803 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
804 phi.data(), diff_phi.data(), nb_gauss_pts, base_polynomials);
805
806 MatrixDouble face_bases(nb_gauss_pts, 3 * nb_dofs_face, false);
807 MatrixDouble diff_face_bases(nb_gauss_pts, 3 * 2 * nb_dofs_face, false);
808 int face_nodes[] = {0, 1, 2};
810 face_nodes, order,
811 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
812 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
813 &*face_bases.data().begin(), &*diff_face_bases.data().begin(),
814 nb_gauss_pts, base_polynomials);
815
816 // edges
817 using Tensor1Ptr3 =
818 decltype(getFTensor1FromPtr<3>(static_cast<double *>(nullptr)));
819 Tensor1Ptr3 t_edge_base[] = {
820 getFTensor1FromPtr<3>(&*edge_bases[0].data().begin()),
821 getFTensor1FromPtr<3>(&*edge_bases[1].data().begin()),
822 getFTensor1FromPtr<3>(&*edge_bases[2].data().begin())};
823 FTensor::Tensor2<FTensor::PackPtr<double *, 6>, 3, 2> t_edge_diff_base[] = {
824 getFTensor2HVecFromPtr<3, 2>(&*diff_edge_bases[0].data().begin()),
825 getFTensor2HVecFromPtr<3, 2>(&*diff_edge_bases[1].data().begin()),
826 getFTensor2HVecFromPtr<3, 2>(&*diff_edge_bases[2].data().begin())};
827
828 // face
829 auto t_vol_base = getFTensor1FromPtr<3>(&*face_bases.data().begin());
830 auto t_vol_diff_base =
831 getFTensor2HVecFromPtr<3, 2>(&*diff_face_bases.data().begin());
832
833 auto t_base = getFTensor1FromPtr<3>(
834 &*data.dataOnEntities[MBTRI][0].getN(base).data().begin());
835 auto t_diff_base = getFTensor2HVecFromPtr<3, 2>(
836 &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin());
837
838 FTENSOR_INDEX(3, i);
839 FTENSOR_INDEX(2, j);
840
841 for (int gg = 0; gg != nb_gauss_pts; gg++) {
842 for (int oo = 0; oo < order; oo++) {
843 // edges
844 for (int dd = NBEDGE_AINSWORTH_HCURL(oo);
845 dd != NBEDGE_AINSWORTH_HCURL(oo + 1); ++dd) {
846 for (int ee = 0; ee != 3; ++ee) {
847 t_base(i) = t_edge_base[ee](i);
848 ++t_base;
849 ++t_edge_base[ee];
850 t_diff_base(i, j) = t_edge_diff_base[ee](i, j);
851 ++t_diff_base;
852 ++t_edge_diff_base[ee];
853 }
854 }
855 // face
856 for (int dd = NBFACETRI_AINSWORTH_HCURL(oo);
857 dd != NBFACETRI_AINSWORTH_HCURL(oo + 1); ++dd) {
858 t_base(i) = t_vol_base(i);
859 ++t_base;
860 ++t_vol_base;
861 t_diff_base(i, j) = t_vol_diff_base(i, j);
862 ++t_diff_base;
863 ++t_vol_diff_base;
864 }
865 }
866 }
867
868 } else {
869 data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
870 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
871 }
872
874}
#define FTENSOR_INDEX(DIM, I)
FTensor::Index< 'j', 3 > j
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
auto getFTensor2HVecFromPtr< 3, 2 >(double *ptr)
static double phi

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 877 of file TriPolynomialBase.cpp.

877 {
879
880 EntitiesFieldData &data = cTx->dAta;
881 const FieldApproximationBase base = cTx->bAse;
882
883 int nb_gauss_pts = pts.size2();
884
885 // Calculation H-curl on triangle faces
886 if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
887
888 if (data.dataOnEntities[MBEDGE].size() != 3)
889 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
890 "wrong number of data structures on edges, should be three but "
891 "is %zu",
892 data.dataOnEntities[MBEDGE].size());
893
894 int sense[3], order[3];
895 double *hcurl_edge_n[3];
896 double *diff_hcurl_edge_n[3];
897
898 for (int ee = 0; ee != 3; ++ee) {
899
900 if (data.dataOnEntities[MBEDGE][ee].getSense() == 0)
901 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
902 "orientation (sense) of edge is not set");
903
904 sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
905 order[ee] = data.dataOnEntities[MBEDGE][ee].getOrder();
906 int nb_dofs =
907 NBEDGE_DEMKOWICZ_HCURL(data.dataOnEntities[MBEDGE][ee].getOrder());
908 data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
909 3 * nb_dofs, false);
910 data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
911 nb_gauss_pts, 2 * 3 * nb_dofs, false);
912
913 hcurl_edge_n[ee] =
914 &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
915 diff_hcurl_edge_n[ee] =
916 &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
917 }
918
920 sense, order,
921 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
922 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
923 hcurl_edge_n, diff_hcurl_edge_n, nb_gauss_pts);
924
925 } else {
926
927 // No DOFs on faces, resize base function matrices, indicating that no
928 // dofs on them.
929 for (int ee = 0; ee != 3; ++ee) {
930 data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, 0, false);
931 data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts, 0,
932 false);
933 }
934 }
935
936 if (data.spacesOnEntities[MBTRI].test(HCURL)) {
937
938 // face
939 if (data.dataOnEntities[MBTRI].size() != 1)
940 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
941 "No data struture to keep base functions on face");
942
943 int order = data.dataOnEntities[MBTRI][0].getOrder();
944 int nb_dofs = NBFACETRI_DEMKOWICZ_HCURL(order);
945 data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
946 false);
947 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
948 3 * 2 * nb_dofs, false);
949
950 int face_nodes[] = {0, 1, 2};
952 face_nodes, order,
953 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
954 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
955 &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
956 &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
957 nb_gauss_pts);
958
959 } else {
960
961 // No DOFs on faces, resize base function matrices, indicating that no
962 // dofs on them.
963 data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
964 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
965 }
966
968}
#define NBEDGE_DEMKOWICZ_HCURL(P)
#define NBFACETRI_DEMKOWICZ_HCURL(P)
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:2116
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:2440

◆ getValueHcurlDemkowiczBrokenBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlDemkowiczBrokenBase ( MatrixDouble pts)
private

Definition at line 971 of file TriPolynomialBase.cpp.

971 {
973
974 EntitiesFieldData &data = cTx->dAta;
975 const FieldApproximationBase base = cTx->bAse;
976
977 int nb_gauss_pts = pts.size2();
978
979 if (data.spacesOnEntities[MBTRI].test(HCURL)) {
980
981 // face
982 if (data.dataOnEntities[MBTRI].size() != 1)
983 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
984 "No data struture to keep base functions on face");
985
986 int order = data.dataOnEntities[MBTRI][0].getOrder();
987 int nb_edge_dofs = NBEDGE_DEMKOWICZ_HCURL(order);
988 int nb_volume_dofs = NBFACETRI_DEMKOWICZ_HCURL(order);
989 int nb_dofs = 3 * nb_edge_dofs + nb_volume_dofs;
990 data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
991 false);
992 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
993 3 * 2 * nb_dofs, false);
994
995 MatrixDouble edge_bases[] = {
996 MatrixDouble(nb_gauss_pts, 3 * nb_edge_dofs, false),
997 MatrixDouble(nb_gauss_pts, 3 * nb_edge_dofs, false),
998 MatrixDouble(nb_gauss_pts, 3 * nb_edge_dofs, false)};
999 MatrixDouble diff_edge_bases[] = {
1000 MatrixDouble(nb_gauss_pts, 3 * 2 * nb_edge_dofs, false),
1001 MatrixDouble(nb_gauss_pts, 3 * 2 * nb_edge_dofs, false),
1002 MatrixDouble(nb_gauss_pts, 3 * 2 * nb_edge_dofs, false)};
1003
1004 std::array<int, 3> edge_order{order, order, order};
1005 std::array<int, 3> edge_sense{1, 1, 1};
1006 std::array<double *, 3> phi{&*edge_bases[0].data().begin(),
1007 &*edge_bases[1].data().begin(),
1008 &*edge_bases[2].data().begin()};
1009 std::array<double *, 3> diff_phi{&*diff_edge_bases[0].data().begin(),
1010 &*diff_edge_bases[1].data().begin(),
1011 &*diff_edge_bases[2].data().begin()};
1013 edge_sense.data(), edge_order.data(),
1014 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
1015 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
1016 phi.data(), diff_phi.data(), nb_gauss_pts);
1017
1018 MatrixDouble face_bases(nb_gauss_pts, 3 * nb_volume_dofs, false);
1019 MatrixDouble diff_face_bases(nb_gauss_pts, 3 * 2 * nb_volume_dofs, false);
1020 int face_nodes[] = {0, 1, 2};
1022 face_nodes, order,
1023 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
1024 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
1025 &*face_bases.data().begin(), &*diff_face_bases.data().begin(),
1026 nb_gauss_pts);
1027
1028 // edges
1029 using Tensor1Ptr3 =
1030 decltype(getFTensor1FromPtr<3>(static_cast<double *>(nullptr)));
1031 Tensor1Ptr3 t_edge_base[] = {
1032 getFTensor1FromPtr<3>(&*edge_bases[0].data().begin()),
1033 getFTensor1FromPtr<3>(&*edge_bases[1].data().begin()),
1034 getFTensor1FromPtr<3>(&*edge_bases[2].data().begin())};
1035 FTensor::Tensor2<FTensor::PackPtr<double *, 6>, 3, 2> t_edge_diff_base[] = {
1036 getFTensor2HVecFromPtr<3, 2>(&*diff_edge_bases[0].data().begin()),
1037 getFTensor2HVecFromPtr<3, 2>(&*diff_edge_bases[1].data().begin()),
1038 getFTensor2HVecFromPtr<3, 2>(&*diff_edge_bases[2].data().begin())};
1039 // face
1040 auto t_vol_base = getFTensor1FromPtr<3>(&*face_bases.data().begin());
1041 auto t_vol_diff_base =
1042 getFTensor2HVecFromPtr<3, 2>(&*diff_face_bases.data().begin());
1043
1044 auto t_base = getFTensor1FromPtr<3>(
1045 &*data.dataOnEntities[MBTRI][0].getN(base).data().begin());
1046 auto t_diff_base = getFTensor2HVecFromPtr<3, 2>(
1047 &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin());
1048
1049 FTENSOR_INDEX(3, i);
1050 FTENSOR_INDEX(2, j);
1051
1052 // Fill the base functions and their derivatives
1053 for (int gg = 0; gg != nb_gauss_pts; gg++) {
1054 for (int oo = 0; oo < order; oo++) {
1055 // edges
1056 for (int dd = NBEDGE_DEMKOWICZ_HCURL(oo);
1057 dd != NBEDGE_DEMKOWICZ_HCURL(oo + 1); ++dd) {
1058 for (int ee = 0; ee != 3; ++ee) {
1059 t_base(i) = t_edge_base[ee](i);
1060 ++t_base;
1061 ++t_edge_base[ee];
1062 t_diff_base(i, j) = t_edge_diff_base[ee](i, j);
1063 ++t_diff_base;
1064 ++t_edge_diff_base[ee];
1065 }
1066 }
1067 // faces
1068 for (int dd = NBFACETRI_DEMKOWICZ_HCURL(oo);
1069 dd != NBFACETRI_DEMKOWICZ_HCURL(oo + 1); ++dd) {
1070 t_base(i) = t_vol_base(i);
1071 ++t_base;
1072 ++t_vol_base;
1073 t_diff_base(i, j) = t_vol_diff_base(i, j);
1074 ++t_diff_base;
1075 ++t_vol_diff_base;
1076 }
1077 }
1078 }
1079
1080 } else {
1081
1082 // No DOFs on faces, resize base function matrices, indicating that no
1083 // dofs on them.
1084 data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
1085 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
1086 }
1087
1089}

◆ getValueHdiv()

MoFEMErrorCode TriPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 646 of file TriPolynomialBase.cpp.

646 {
648
649 switch (cTx->bAse) {
655 default:
656 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
657 }
658
660}
MoFEMErrorCode getValueHdivAinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueHdivDemkowiczBase(MatrixDouble &pts)

◆ getValueHdivAinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueHdivAinsworthBase ( MatrixDouble pts)
private

Definition at line 524 of file TriPolynomialBase.cpp.

524 {
526
527 EntitiesFieldData &data = cTx->dAta;
528 const FieldApproximationBase base = cTx->bAse;
529 if (cTx->basePolynomialsType0 == NULL)
530 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
531 "Polynomial type not set");
532 PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
533 double *diffL, const int dim) =
535
536 int nb_gauss_pts = pts.size2();
537
538 double *phi_f_e[3];
539 double *phi_f;
540
541 N_face_edge.resize(1, 3, false);
542 N_face_bubble.resize(1, false);
543 int face_order = data.dataOnEntities[MBTRI][0].getOrder();
544 auto nb_dofs_on_face =
549 if (nb_dofs_on_face) {
550 auto face_edge_dofs = NBFACETRI_AINSWORTH_EDGE_HDIV(
552 // three edges on face
553 for (int ee = 0; ee < 3; ee++) {
554 N_face_edge(0, ee).resize(nb_gauss_pts, 3 * face_edge_dofs, false);
555 phi_f_e[ee] = &((N_face_edge(0, ee))(0, 0));
556 }
557 auto face_bubble_dofs = NBFACETRI_AINSWORTH_FACE_HDIV(
559 N_face_bubble[0].resize(nb_gauss_pts, 3 * face_bubble_dofs, false);
560 phi_f = &*(N_face_bubble[0].data().begin());
561
562 int face_nodes[3] = {0, 1, 2};
563 if (face_edge_dofs)
565 face_nodes,
567 &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0), NULL, phi_f_e,
568 NULL, nb_gauss_pts, 3, base_polynomials);
569 if (face_bubble_dofs)
571 face_nodes,
573 &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0), NULL, phi_f, NULL,
574 nb_gauss_pts, 3, base_polynomials);
575
576 // set shape functions into data structure
577 if (data.dataOnEntities[MBTRI].size() != 1) {
578 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
579 }
580 data.dataOnEntities[MBTRI][0].getN(base).resize(
581 nb_gauss_pts, 3 * NBFACETRI_AINSWORTH_HDIV(face_order), false);
582
583 auto max_face_order =
584 std::max(face_order,
586
587 int col = 0;
588 for (int oo = 0; oo < max_face_order; oo++) {
590 for (int ee = 0; ee < 3; ee++) {
591 for (int dd = 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(oo);
592 dd < 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(oo + 1); dd++, col++) {
593 for (int gg = 0; gg < nb_gauss_pts; gg++) {
594 data.dataOnEntities[MBTRI][0].getN(base)(gg, col) =
595 N_face_edge(0, ee)(gg, dd);
596 }
597 }
598 }
599
601 for (int dd = 3 * NBFACETRI_AINSWORTH_FACE_HDIV(oo);
602 dd < 3 * NBFACETRI_AINSWORTH_FACE_HDIV(oo + 1); dd++, col++) {
603 for (int gg = 0; gg < nb_gauss_pts; gg++) {
604 data.dataOnEntities[MBTRI][0].getN(base)(gg, col) =
605 N_face_bubble[0](gg, dd);
606 }
607 }
608 }
609 }
610
612}
#define NBFACETRI_AINSWORTH_FACE_HDIV(P)
#define NBFACETRI_AINSWORTH_EDGE_HDIV(P)
#define NBFACETRI_AINSWORTH_HDIV(P)
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:97
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:224
static boost::function< int(int)> broken_nbfacetri_face_hdiv
Definition Hdiv.hpp:26
static boost::function< int(int)> broken_nbfacetri_edge_hdiv
Definition Hdiv.hpp:25
ublas::matrix< MatrixDouble > N_face_edge
ublas::vector< MatrixDouble > N_face_bubble

◆ getValueHdivDemkowiczBase()

MoFEMErrorCode TriPolynomialBase::getValueHdivDemkowiczBase ( MatrixDouble pts)
private

Definition at line 614 of file TriPolynomialBase.cpp.

614 {
616
617 EntitiesFieldData &data = cTx->dAta;
618 // set shape functions into data structure
619 if (data.dataOnEntities[MBTRI].size() != 1) {
620 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
621 }
622 const FieldApproximationBase base = cTx->bAse;
623 if (base != DEMKOWICZ_JACOBI_BASE) {
624 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
625 "This should be used only with DEMKOWICZ_JACOBI_BASE "
626 "but base is %s",
628 }
629 int order = data.dataOnEntities[MBTRI][0].getOrder();
630 int nb_gauss_pts = pts.size2();
631 data.dataOnEntities[MBTRI][0].getN(base).resize(
632 nb_gauss_pts, 3 * NBFACETRI_DEMKOWICZ_HDIV(order), false);
633 double *phi_f = &*data.dataOnEntities[MBTRI][0].getN(base).data().begin();
636 int face_nodes[3] = {0, 1, 2};
638 face_nodes, order,
639 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
640 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(), phi_f,
641 NULL, nb_gauss_pts, 3);
642
644}
#define NBFACETRI_DEMKOWICZ_HDIV(P)
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:683

◆ getValueL2()

MoFEMErrorCode TriPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 319 of file TriPolynomialBase.cpp.

319 {
321
322 switch (cTx->bAse) {
327 break;
330 break;
331 default:
332 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
333 }
334
336}
MoFEMErrorCode getValueL2BernsteinBezierBase(MatrixDouble &pts)
MoFEMErrorCode getValueL2AinsworthBase(MatrixDouble &pts)

◆ getValueL2AinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueL2AinsworthBase ( MatrixDouble pts)
private

Definition at line 338 of file TriPolynomialBase.cpp.

338 {
340
341 EntitiesFieldData &data = cTx->dAta;
342 const FieldApproximationBase base = cTx->bAse;
343 if (cTx->basePolynomialsType0 == NULL)
344 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
345 "Polynomial type not set");
346 PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
347 double *diffL, const int dim) =
349
350 int nb_gauss_pts = pts.size2();
351
352 data.dataOnEntities[MBTRI][0].getN(base).resize(
353 nb_gauss_pts, NBFACETRI_L2(data.dataOnEntities[MBTRI][0].getOrder()),
354 false);
355 data.dataOnEntities[MBTRI][0].getDiffN(base).resize(
356 nb_gauss_pts, 2 * NBFACETRI_L2(data.dataOnEntities[MBTRI][0].getOrder()),
357 false);
358
360 data.dataOnEntities[MBTRI][0].getOrder(),
361 &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
362 &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
363 &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
364 &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
365 nb_gauss_pts, base_polynomials);
366
368}
#define NBFACETRI_L2(P)
Number of base functions on triangle for L2 space.
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:19

◆ getValueL2BernsteinBezierBase()

MoFEMErrorCode TriPolynomialBase::getValueL2BernsteinBezierBase ( MatrixDouble pts)
private

Definition at line 371 of file TriPolynomialBase.cpp.

371 {
373
374 EntitiesFieldData &data = cTx->dAta;
375 const std::string &field_name = cTx->fieldName;
376 int nb_gauss_pts = pts.size2();
377
378 auto get_alpha = [field_name](auto &data) -> MatrixInt & {
379 auto &ptr = data.getBBAlphaIndicesSharedPtr(field_name);
380 if (!ptr)
381 ptr.reset(new MatrixInt());
382 return *ptr;
383 };
384
385 auto get_base = [field_name](auto &data) -> MatrixDouble & {
386 auto &ptr = data.getBBNSharedPtr(field_name);
387 if (!ptr)
388 ptr.reset(new MatrixDouble());
389 return *ptr;
390 };
391
392 auto get_diff_base = [field_name](auto &data) -> MatrixDouble & {
393 auto &ptr = data.getBBDiffNSharedPtr(field_name);
394 if (!ptr)
395 ptr.reset(new MatrixDouble());
396 return *ptr;
397 };
398
399 auto get_alpha_by_name_ptr =
400 [](auto &data,
401 const std::string &field_name) -> boost::shared_ptr<MatrixInt> & {
402 return data.getBBAlphaIndicesSharedPtr(field_name);
403 };
404
405 auto get_base_by_name_ptr =
406 [](auto &data,
407 const std::string &field_name) -> boost::shared_ptr<MatrixDouble> & {
408 return data.getBBNSharedPtr(field_name);
409 };
410
411 auto get_diff_base_by_name_ptr =
412 [](auto &data,
413 const std::string &field_name) -> boost::shared_ptr<MatrixDouble> & {
414 return data.getBBDiffNSharedPtr(field_name);
415 };
416
417 auto get_alpha_by_order_ptr =
418 [](auto &data, const size_t o) -> boost::shared_ptr<MatrixInt> & {
419 return data.getBBAlphaIndicesByOrderSharedPtr(o);
420 };
421
422 auto get_base_by_order_ptr =
423 [](auto &data, const size_t o) -> boost::shared_ptr<MatrixDouble> & {
424 return data.getBBNByOrderSharedPtr(o);
425 };
426
427 auto get_diff_base_by_order_ptr =
428 [](auto &data, const size_t o) -> boost::shared_ptr<MatrixDouble> & {
429 return data.getBBDiffNByOrderSharedPtr(o);
430 };
431
432 if (data.spacesOnEntities[MBTRI].test(L2)) {
433 if (data.dataOnEntities[MBTRI].size() != 1)
434 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
435 "Wrong size ent of ent data");
436
437 if (data.dataOnEntities[MBVERTEX][0].getN(NOBASE).size1() !=
438 (unsigned int)nb_gauss_pts)
439 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
440 "Base functions or nodes has wrong number of integration points "
441 "for base %s",
443 auto &lambda = data.dataOnEntities[MBVERTEX][0].getN(NOBASE);
444
445 auto &ent_data = data.dataOnEntities[MBTRI][0];
446 const int order = ent_data.getOrder();
447 const int nb_dofs = NBFACETRI_L2(order);
448
449 if (get_alpha_by_order_ptr(ent_data, order)) {
450 get_alpha_by_name_ptr(ent_data, field_name) =
451 get_alpha_by_order_ptr(ent_data, order);
452 get_base_by_name_ptr(ent_data, field_name) =
453 get_base_by_order_ptr(ent_data, order);
454 get_diff_base_by_name_ptr(ent_data, field_name) =
455 get_diff_base_by_order_ptr(ent_data, order);
456 } else {
457
458 auto &get_n = get_base(ent_data);
459 auto &get_diff_n = get_diff_base(ent_data);
460 get_n.resize(nb_gauss_pts, nb_dofs, false);
461 get_diff_n.resize(nb_gauss_pts, 2 * nb_dofs, false);
462 if (nb_dofs) {
463
464 if (order == 0) {
465
466 if (nb_dofs != 1)
467 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
468 "Inconsistent number of DOFs");
469
470 auto &tri_alpha = get_alpha(ent_data);
471 tri_alpha.clear();
472 get_n(0, 0) = 1;
473 get_diff_n.clear();
474
475 } else {
476
477 if (nb_dofs != 3 + 3 * NBEDGE_H1(order) + NBFACETRI_H1(order))
478 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
479 "Inconsistent number of DOFs");
480
481 auto &tri_alpha = get_alpha(ent_data);
482 tri_alpha.resize(nb_dofs, 3, false);
483
485 &tri_alpha(0, 0));
486
487 if (order > 1) {
488 std::array<int, 3> face_n{order, order, order};
489 std::array<int *, 3> face_edge_ptr{
490 &tri_alpha(3, 0), &tri_alpha(3 + NBEDGE_H1(order), 0),
491 &tri_alpha(3 + 2 * NBEDGE_H1(order), 0)};
493 face_n.data(), face_edge_ptr.data());
494 if (order > 2)
496 order, &tri_alpha(3 + 3 * NBEDGE_H1(order), 0));
497 }
499 order, lambda.size1(), tri_alpha.size1(), &tri_alpha(0, 0),
500 &lambda(0, 0), Tools::diffShapeFunMBTRI.data(), &get_n(0, 0),
501 &get_diff_n(0, 0));
502 }
503
504 get_alpha_by_order_ptr(ent_data, order) =
505 get_alpha_by_name_ptr(ent_data, field_name);
506 get_base_by_order_ptr(ent_data, order) =
507 get_base_by_name_ptr(ent_data, field_name);
508 get_diff_base_by_order_ptr(ent_data, order) =
509 get_diff_base_by_name_ptr(ent_data, field_name);
510 }
511 }
512 } else {
513 auto &ent_data = data.dataOnEntities[MBTRI][0];
514 ent_data.getBBAlphaIndicesSharedPtr(field_name).reset();
515 auto &get_n = get_base(ent_data);
516 auto &get_diff_n = get_diff_base(ent_data);
517 get_n.resize(nb_gauss_pts, 0, false);
518 get_diff_n.resize(nb_gauss_pts, 0, false);
519 }
520
522}
static MoFEMErrorCode generateIndicesVertexTri(const int N, int *alpha)

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 9 of file TriPolynomialBase.cpp.

10 {
12 *iface = const_cast<TriPolynomialBase *>(this);
14}
Calculate base functions on triangle.

◆ setDofsSideMap()

MoFEMErrorCode TriPolynomialBase::setDofsSideMap ( const FieldSpace  space,
const FieldContinuity  continuity,
const FieldApproximationBase  base,
DofsSideMap dofs_side_map 
)
static

Set map of dof to side number.

That is used for broken space to establish connection between dofs in the interior of element/entity and side of element/entity to which that dof is associated. That depends on implementation of the base for given space, and has to be implemented while implementing base function for given space.

Parameters
space
continuity
base
DofsSideMap
Returns
MoFEMErrorCode

Definition at line 1218 of file TriPolynomialBase.cpp.

1223 {
1225
1226 switch (continuity) {
1227 case DISCONTINUOUS:
1228
1229 switch (space) {
1230 case HCURL:
1232 setDofsSideMapHcurl(space, continuity, base, dofs_side_map));
1233 default:
1234 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Unknown space %s",
1235 FieldSpaceNames[space]);
1236 }
1237 break;
1238
1239 default:
1240 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1241 "Unknown (or not implemented) continuity");
1242 }
1243
1245}
static const char *const FieldSpaceNames[]
Definition definitions.h:92
static MoFEMErrorCode setDofsSideMapHcurl(const FieldSpace space, const FieldContinuity continuity, const FieldApproximationBase base, DofsSideMap &dofs_side_map)

◆ setDofsSideMapHcurl()

MoFEMErrorCode TriPolynomialBase::setDofsSideMapHcurl ( const FieldSpace  space,
const FieldContinuity  continuity,
const FieldApproximationBase  base,
DofsSideMap dofs_side_map 
)
staticprivate

Definition at line 1247 of file TriPolynomialBase.cpp.

1249 {
1251
1252 // That has to be consistent with implementation of getValueHdiv for
1253 // particular base functions.
1254
1255 auto set_ainsworth = [&dofs_side_map]() {
1257
1258 dofs_side_map.clear();
1259
1260 int dof = 0;
1261 for (int oo = 0; oo < Field::maxBrokenDofsOrder; oo++) {
1262
1263 // edges
1264 for (int dd = NBEDGE_AINSWORTH_HCURL(oo);
1265 dd != NBEDGE_AINSWORTH_HCURL(oo + 1); ++dd) {
1266 for (int ee = 0; ee != 3; ++ee) {
1267 dofs_side_map.insert(DofsSideMapData{MBEDGE, ee, dof});
1268 ++dof;
1269 }
1270 }
1271 // face
1272 for (int dd = NBFACETRI_AINSWORTH_HCURL(oo);
1273 dd != NBFACETRI_AINSWORTH_HCURL(oo + 1); ++dd) {
1274 dofs_side_map.insert(DofsSideMapData{MBTRI, 0, dof});
1275 ++dof;
1276 }
1277 }
1278
1280 };
1281
1282 auto set_demkowicz = [&dofs_side_map]() {
1284
1285 dofs_side_map.clear();
1286
1287 int dof = 0;
1288 for (int oo = 0; oo < Field::maxBrokenDofsOrder; oo++) {
1289
1290 // edges
1291 for (int dd = NBEDGE_DEMKOWICZ_HCURL(oo);
1292 dd != NBEDGE_DEMKOWICZ_HCURL(oo + 1); ++dd) {
1293 for (int ee = 0; ee != 3; ++ee) {
1294 dofs_side_map.insert(DofsSideMapData{MBEDGE, ee, dof});
1295 ++dof;
1296 }
1297 }
1298 // faces
1299 for (int dd = NBFACETRI_DEMKOWICZ_HCURL(oo);
1300 dd != NBFACETRI_DEMKOWICZ_HCURL(oo + 1); ++dd) {
1301 dofs_side_map.insert(DofsSideMapData{MBTRI, 0, dof});
1302 ++dof;
1303 }
1304 }
1305
1307 };
1308
1309 switch (continuity) {
1310 case DISCONTINUOUS:
1311 switch (base) {
1314 MoFEMFunctionReturnHot(set_ainsworth());
1316 MoFEMFunctionReturnHot(set_demkowicz());
1317 default:
1318 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
1319 }
1320 break;
1321
1322 default:
1323 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1324 "Unknown (or not implemented) continuity");
1325 }
1326
1328}
static constexpr int maxBrokenDofsOrder
Maximum order for broken space DOFs.

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* MoFEM::TriPolynomialBase::cTx
private

Definition at line 47 of file TriPolynomialBase.hpp.

◆ diffN_face_bubble

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

Definition at line 60 of file TriPolynomialBase.hpp.

◆ diffN_face_edge

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

Definition at line 59 of file TriPolynomialBase.hpp.

◆ N_face_bubble

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

Definition at line 58 of file TriPolynomialBase.hpp.

◆ N_face_edge

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

Definition at line 57 of file TriPolynomialBase.hpp.


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