v0.8.20
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, UnknownInterface **iface) const
 
 EdgePolynomialBase ()
 
 ~EdgePolynomialBase ()
 
MoFEMErrorCode getValue (MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::BaseFunction
 BaseFunction ()
 
 ~BaseFunction ()
 
virtual MoFEMErrorCode getValue (MatrixDouble &pts_x, MatrixDouble &pts_t, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of 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 ()
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 

Private Member Functions

MoFEMErrorCode getValueH1 (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
 

Additional Inherited Members

- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

Calculate base functions on tetrahedral.

Definition at line 33 of file EdgePolynomialBase.hpp.

Constructor & Destructor Documentation

◆ EdgePolynomialBase()

EdgePolynomialBase::EdgePolynomialBase ( )

Definition at line 39 of file EdgePolynomialBase.cpp.

39 {}

◆ ~EdgePolynomialBase()

EdgePolynomialBase::~EdgePolynomialBase ( )

Definition at line 38 of file EdgePolynomialBase.cpp.

38 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 42 of file EdgePolynomialBase.cpp.

43  {
44 
46 
48  ierr = ctx_ptr->query_interface(IDD_EDGE_BASE_FUNCTION, &iface);
49  CHKERRG(ierr);
50  cTx = reinterpret_cast<EntPolynomialBaseCtx *>(iface);
51 
52  int nb_gauss_pts = pts.size2();
53  if (!nb_gauss_pts) {
55  }
56 
57  if (pts.size1() < 1) {
58  SETERRQ(
59  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
60  "Wrong dimension of pts, should be at least 3 rows with coordinates");
61  }
62 
63  const FieldApproximationBase base = cTx->bAse;
65  if (cTx->copyNodeBase == LASTBASE) {
66  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 2, false);
67  ierr = ShapeMBEDGE(
68  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
69  &pts(0, 0), nb_gauss_pts);
70  CHKERRG(ierr);
71  } else {
72  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
73  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
74  }
75  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
76  (unsigned int)nb_gauss_pts) {
77  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
78  "Base functions or nodes has wrong number of integration points "
79  "for base %s",
81  }
82  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(2, 1, false);
84  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
85  CHKERRG(ierr);
86 
87  switch (cTx->sPace) {
88  case H1:
89  ierr = getValueH1(pts);
90  CHKERRG(ierr);
91  break;
92  case HDIV:
93  ierr = getValueHdiv(pts);
94  CHKERRG(ierr);
95  break;
96  case HCURL:
97  ierr = getValueHcurl(pts);
98  CHKERRG(ierr);
99  break;
100  case L2:
101  ierr = getValueL2(pts);
102  CHKERRG(ierr);
103  break;
104  default:
105  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
106  }
107 
109 }
field with continuous normal traction
Definition: definitions.h:171
data structure for finite element entityIt keeps that about indices of degrees of freedom...
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
MoFEMErrorCode getValueH1(MatrixDouble &pts)
Class used to pass element data to calculate base functions on tet,triangle,edge. ...
const FieldApproximationBase bAse
base class for all interface classes
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
MoFEMErrorCode getValueL2(MatrixDouble &pts)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
const FieldApproximationBase copyNodeBase
static const char *const ApproximationBaseNames[]
Definition: definitions.h:156
PetscErrorCode ShapeMBEDGE(double *N, const double *G_X, int DIM)
Definition: fem_tools.c:635
FieldApproximationBase
approximation base
Definition: definitions.h:141
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:170
static const MOFEMuuid IDD_EDGE_BASE_FUNCTION
PetscErrorCode ShapeDiffMBEDGE(double *diffN)
Definition: fem_tools.c:645
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
continuous field
Definition: definitions.h:169
EntPolynomialBaseCtx * cTx
field with C-1 continuity
Definition: definitions.h:172

◆ getValueH1()

MoFEMErrorCode EdgePolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 111 of file EdgePolynomialBase.cpp.

111  {
112 
114 
116  const FieldApproximationBase base = cTx->bAse;
117  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
118  double *diffL, const int dim) =
120 
121  int nb_gauss_pts = pts.size2();
122 
123  // std::cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << std::endl;
124  // std::cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << std::endl;
125  //
126  // std::cerr << pts << std::endl;
127 
128  const int side_number = 0;
129  int sense = data.dataOnEntities[MBEDGE][side_number].getSense();
130  int order = data.dataOnEntities[MBEDGE][side_number].getDataOrder();
131  data.dataOnEntities[MBEDGE][side_number].getN(base).resize(
132  nb_gauss_pts, NBEDGE_H1(order), false);
133  data.dataOnEntities[MBEDGE][side_number].getDiffN(base).resize(
134  nb_gauss_pts, NBEDGE_H1(order), false);
135 
136  data.dataOnEntities[MBEDGE][side_number].getN(base).clear();
137  data.dataOnEntities[MBEDGE][side_number].getDiffN(base).clear();
138 
139  L.resize(NBEDGE_H1(order), false);
140  diffL.resize(NBEDGE_H1(order), false);
141 
142  // std::cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << std::endl;
143 
144  if (data.dataOnEntities[MBEDGE][side_number].getDataOrder() > 1) {
145 
146  double diff_s = 2.; // s = s(xi), ds/dxi = 2., because change of basis
147  for (int gg = 0; gg < nb_gauss_pts; gg++) {
148 
149  double s = 2 * pts(0, gg) - 1; // makes form -1..1
150  if (!sense) {
151  s *= -1;
152  diff_s *= -1;
153  }
154 
155  // calculate Legendre polynomials at integration points
156  ierr = base_polynomials(NBEDGE_H1(order) - 1, s, &diff_s,
157  &*L.data().begin(), &*diffL.data().begin(), 1);
158  CHKERRG(ierr);
159 
160  // std::cerr << "s " << s << " " << L << std::endl;
161 
162  for (unsigned int pp = 0;
163  pp < data.dataOnEntities[MBEDGE][side_number].getN(base).size2();
164  pp++) {
165 
166  // Calculate edge shape functions N0*N1*L(p), where N0 and N1 are nodal
167  // shape functions
168  double v = data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 0) *
169  data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 1);
170  data.dataOnEntities[MBEDGE][side_number].getN(base)(gg, pp) = v * L(pp);
171 
172  // Calculate derivative edge shape functions
173  // dN/dksi = dN0/dxi*N1*L + N0*dN1/ksi*L + N0*N1*dL/dxi
174  data.dataOnEntities[MBEDGE][side_number].getDiffN(base)(gg, pp) =
175  ((+1.) * data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 1) +
176  data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 0) * (-1.)) *
177  L(pp) +
178  v * diffL(pp);
179  }
180  }
181  }
182 
183  // std::cerr << data.dataOnEntities[MBEDGE][0].getN(base) << std::endl;
184  // std::cerr << data.dataOnEntities[MBEDGE][0].getDiffN(base) << std::endl;
185 
187 }
#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...
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
const FieldApproximationBase bAse
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
FieldApproximationBase
approximation base
Definition: definitions.h:141
DataForcesAndSourcesCore & dAta
EntPolynomialBaseCtx * cTx

◆ getValueHcurl()

MoFEMErrorCode EdgePolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 278 of file EdgePolynomialBase.cpp.

278  {
280 
281  switch (cTx->bAse) {
285  break;
288  break;
289  default:
290  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
291  }
292 
294 }
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:143
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
EntPolynomialBaseCtx * cTx

◆ getValueHcurlAinsworthBase()

MoFEMErrorCode EdgePolynomialBase::getValueHcurlAinsworthBase ( MatrixDouble pts)
private

Definition at line 201 of file EdgePolynomialBase.cpp.

201  {
203 
205  const FieldApproximationBase base = cTx->bAse;
206 
207  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
208  double *diffL, const int dim) =
210 
211  int nb_gauss_pts = pts.size2();
212  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
213  // edges
214  if (data.dataOnEntities[MBEDGE].size() != 1) {
215  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
216  }
217  int sense = data.dataOnEntities[MBEDGE][0].getSense();
218  int order = data.dataOnEntities[MBEDGE][0].getDataOrder();
219  int nb_dofs =
220  NBEDGE_AINSWORTH_HCURL(data.dataOnEntities[MBEDGE][0].getDataOrder());
221  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
222  false);
223  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
224  false);
225  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
227  sense, order, &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0),
228  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
229  &*data.dataOnEntities[MBEDGE][0].getN(base).data().begin(), NULL,
230  nb_gauss_pts, base_polynomials);
231  CHKERRG(ierr);
232  } else {
233  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, false);
234  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
235  false);
236  }
237 
239 }
std::bitset< LASTSPACE > spacesOnEntities[MBMAXTYPE]
spaces on entity types
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...
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
const FieldApproximationBase bAse
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
FieldApproximationBase
approximation base
Definition: definitions.h:141
DataForcesAndSourcesCore & dAta
#define NBEDGE_AINSWORTH_HCURL(P)
field with continuous tangents
Definition: definitions.h:170
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:174

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode EdgePolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 242 of file EdgePolynomialBase.cpp.

242  {
244 
246  const FieldApproximationBase base = cTx->bAse;
247 
248  int nb_gauss_pts = pts.size2();
249  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
250  // edges
251  if (data.dataOnEntities[MBEDGE].size() != 1) {
252  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
253  "No data structure to store base functions");
254  }
255  int sense = data.dataOnEntities[MBEDGE][0].getSense();
256  int order = data.dataOnEntities[MBEDGE][0].getDataOrder();
257  int nb_dofs =
258  NBEDGE_DEMKOWICZ_HCURL(data.dataOnEntities[MBEDGE][0].getDataOrder());
259  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
260  false);
261  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
262  false);
263  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
265  sense, order, &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0),
266  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
267  &*data.dataOnEntities[MBEDGE][0].getN(base).data().begin(), NULL,
268  nb_gauss_pts);
269  } else {
270  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, false);
271  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
272  false);
273  }
274 
276 }
std::bitset< LASTSPACE > spacesOnEntities[MBMAXTYPE]
spaces on entity types
data structure for finite element entityIt keeps that about indices of degrees of freedom...
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
FieldApproximationBase
approximation base
Definition: definitions.h:141
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:2278
field with continuous tangents
Definition: definitions.h:170
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define NBEDGE_DEMKOWICZ_HCURL(P)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
EntPolynomialBaseCtx * cTx

◆ getValueHdiv()

MoFEMErrorCode EdgePolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 196 of file EdgePolynomialBase.cpp.

196  {
197  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
198  "Make no sense, unless problem is 2d (2d not implemented yet)");
199 }

◆ getValueL2()

MoFEMErrorCode EdgePolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 189 of file EdgePolynomialBase.cpp.

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

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 22 of file EdgePolynomialBase.cpp.

23  {
24 
26  *iface = NULL;
27  if (uuid == IDD_EDGE_BASE_FUNCTION) {
28  *iface = const_cast<EdgePolynomialBase *>(this);
30  } else {
31  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
32  }
33  ierr = BaseFunction::query_interface(uuid, iface);
34  CHKERRG(ierr);
36 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
Calculate base functions on tetrahedral.
static const MOFEMuuid IDD_EDGE_BASE_FUNCTION
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, MoFEM::UnknownInterface **iface) const

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: