v0.8.4
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 67 of file EdgePolynomialBase.cpp.

67 {}

◆ ~EdgePolynomialBase()

EdgePolynomialBase::~EdgePolynomialBase ( )

Definition at line 66 of file EdgePolynomialBase.cpp.

66 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 70 of file EdgePolynomialBase.cpp.

71  {
72 
74 
76  ierr = ctx_ptr->query_interface(IDD_EDGE_BASE_FUNCTION, &iface);
77  CHKERRG(ierr);
78  cTx = reinterpret_cast<EntPolynomialBaseCtx *>(iface);
79 
80  int nb_gauss_pts = pts.size2();
81  if (!nb_gauss_pts) {
83  }
84 
85  if (pts.size1() < 1) {
86  SETERRQ(
87  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
88  "Wrong dimension of pts, should be at least 3 rows with coordinates");
89  }
90 
91  const FieldApproximationBase base = cTx->bAse;
93  if (cTx->copyNodeBase == LASTBASE) {
94  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 2, false);
95  ierr = ShapeMBEDGE(
96  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
97  &pts(0, 0), nb_gauss_pts);
98  CHKERRG(ierr);
99  } else {
100  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
101  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
102  }
103  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
104  (unsigned int)nb_gauss_pts) {
105  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
106  "Base functions or nodes has wrong number of integration points "
107  "for base %s",
108  ApproximationBaseNames[base]);
109  }
110  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(2, 1, false);
112  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
113  CHKERRG(ierr);
114 
115  switch (cTx->sPace) {
116  case H1:
117  ierr = getValueH1(pts);
118  CHKERRG(ierr);
119  break;
120  case HDIV:
121  ierr = getValueHdiv(pts);
122  CHKERRG(ierr);
123  break;
124  case HCURL:
125  ierr = getValueHcurl(pts);
126  CHKERRG(ierr);
127  break;
128  case L2:
129  ierr = getValueL2(pts);
130  CHKERRG(ierr);
131  break;
132  default:
133  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
134  }
135 
137 }
field with continuous normal traction
Definition: definitions.h:179
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:522
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:565
MoFEMErrorCode getValueL2(MatrixDouble &pts)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
const FieldApproximationBase copyNodeBase
static const char *const ApproximationBaseNames[]
Definition: definitions.h:153
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
PetscErrorCode ShapeMBEDGE(double *N, const double *G_X, int DIM)
Definition: fem_tools.c:635
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:178
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:177
EntPolynomialBaseCtx * cTx
field with C-1 continuity
Definition: definitions.h:180

◆ getValueH1()

MoFEMErrorCode EdgePolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 139 of file EdgePolynomialBase.cpp.

139  {
140 
142 
144  const FieldApproximationBase base = cTx->bAse;
145  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
146  double *diffL, const int dim) =
148 
149  int nb_gauss_pts = pts.size2();
150 
151  // std::cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << std::endl;
152  // std::cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << std::endl;
153  //
154  // std::cerr << pts << std::endl;
155 
156  const int side_number = 0;
157  int sense = data.dataOnEntities[MBEDGE][side_number].getSense();
158  int order = data.dataOnEntities[MBEDGE][side_number].getDataOrder();
159  data.dataOnEntities[MBEDGE][side_number].getN(base).resize(
160  nb_gauss_pts, NBEDGE_H1(order), false);
161  data.dataOnEntities[MBEDGE][side_number].getDiffN(base).resize(
162  nb_gauss_pts, NBEDGE_H1(order), false);
163 
164  data.dataOnEntities[MBEDGE][side_number].getN(base).clear();
165  data.dataOnEntities[MBEDGE][side_number].getDiffN(base).clear();
166 
167  L.resize(NBEDGE_H1(order), false);
168  diffL.resize(NBEDGE_H1(order), false);
169 
170  // std::cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << std::endl;
171 
172  if (data.dataOnEntities[MBEDGE][side_number].getDataOrder() > 1) {
173 
174  double diff_s = 2.; // s = s(xi), ds/dxi = 2., because change of basis
175  for (int gg = 0; gg < nb_gauss_pts; gg++) {
176 
177  double s = 2 * pts(0, gg) - 1; // makes form -1..1
178  if (!sense) {
179  s *= -1;
180  diff_s *= -1;
181  }
182 
183  // calculate Legendre polynomials at integration points
184  ierr = base_polynomials(NBEDGE_H1(order) - 1, s, &diff_s,
185  &*L.data().begin(), &*diffL.data().begin(), 1);
186  CHKERRG(ierr);
187 
188  // std::cerr << "s " << s << " " << L << std::endl;
189 
190  for (unsigned int pp = 0;
191  pp < data.dataOnEntities[MBEDGE][side_number].getN(base).size2();
192  pp++) {
193 
194  // Calculate edge shape functions N0*N1*L(p), where N0 and N1 are nodal
195  // shape functions
196  double v = data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 0) *
197  data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 1);
198  data.dataOnEntities[MBEDGE][side_number].getN(base)(gg, pp) = v * L(pp);
199 
200  // Calculate derivative edge shape functions
201  // dN/dksi = dN0/dxi*N1*L + N0*dN1/ksi*L + N0*N1*dL/dxi
202  data.dataOnEntities[MBEDGE][side_number].getDiffN(base)(gg, pp) =
203  ((+1.) * data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 1) +
204  data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 0) * (-1.)) *
205  L(pp) +
206  v * diffL(pp);
207  }
208  }
209  }
210 
211  // std::cerr << data.dataOnEntities[MBEDGE][0].getN(base) << std::endl;
212  // std::cerr << data.dataOnEntities[MBEDGE][0].getDiffN(base) << std::endl;
213 
215 }
#define NBEDGE_H1(P)
Number of dofs 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:522
const FieldApproximationBase bAse
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:565
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
EntPolynomialBaseCtx * cTx

◆ getValueHcurl()

MoFEMErrorCode EdgePolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 306 of file EdgePolynomialBase.cpp.

306  {
308 
309  switch (cTx->bAse) {
313  break;
316  break;
317  default:
318  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
319  }
320 
322 }
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:142
#define CHKERR
Inline error check.
Definition: definitions.h:617
Like AINSWORTH_LEGENDRE_BASE but with Lobatto base instead Legendre .
Definition: definitions.h:143
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
EntPolynomialBaseCtx * cTx
Construction of base is by Demkowicz .
Definition: definitions.h:145

◆ getValueHcurlAinsworthBase()

MoFEMErrorCode EdgePolynomialBase::getValueHcurlAinsworthBase ( MatrixDouble pts)
private

Definition at line 229 of file EdgePolynomialBase.cpp.

229  {
231 
233  const FieldApproximationBase base = cTx->bAse;
234 
235  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
236  double *diffL, const int dim) =
238 
239  int nb_gauss_pts = pts.size2();
240  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
241  // edges
242  if (data.dataOnEntities[MBEDGE].size() != 1) {
243  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
244  }
245  int sense = data.dataOnEntities[MBEDGE][0].getSense();
246  int order = data.dataOnEntities[MBEDGE][0].getDataOrder();
247  int nb_dofs =
248  NBEDGE_AINSWORTH_HCURL(data.dataOnEntities[MBEDGE][0].getDataOrder());
249  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
250  false);
251  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
252  false);
253  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
255  sense, order, &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0),
256  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
257  &*data.dataOnEntities[MBEDGE][0].getN(base).data().begin(), NULL,
258  nb_gauss_pts, base_polynomials);
259  CHKERRG(ierr);
260  } else {
261  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, false);
262  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
263  false);
264  }
265 
267 }
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:522
const FieldApproximationBase bAse
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:565
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
#define NBEDGE_AINSWORTH_HCURL(P)
field with continuous tangents
Definition: definitions.h:178
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:181

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode EdgePolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 270 of file EdgePolynomialBase.cpp.

270  {
272 
274  const FieldApproximationBase base = cTx->bAse;
275 
276  int nb_gauss_pts = pts.size2();
277  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
278  // edges
279  if (data.dataOnEntities[MBEDGE].size() != 1) {
280  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
281  "No data structure to store base functions");
282  }
283  int sense = data.dataOnEntities[MBEDGE][0].getSense();
284  int order = data.dataOnEntities[MBEDGE][0].getDataOrder();
285  int nb_dofs =
286  NBEDGE_DEMKOWICZ_HCURL(data.dataOnEntities[MBEDGE][0].getDataOrder());
287  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
288  false);
289  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
290  false);
291  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
293  sense, order, &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0),
294  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
295  &*data.dataOnEntities[MBEDGE][0].getN(base).data().begin(), NULL,
296  nb_gauss_pts);
297  } else {
298  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, false);
299  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
300  false);
301  }
302 
304 }
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:498
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
FieldApproximationBase
approximation base
Definition: definitions.h:140
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:2284
field with continuous tangents
Definition: definitions.h:178
#define CHKERR
Inline error check.
Definition: definitions.h:617
#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:443
EntPolynomialBaseCtx * cTx

◆ getValueHdiv()

MoFEMErrorCode EdgePolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 224 of file EdgePolynomialBase.cpp.

224  {
225  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
226  "Make no sense, unless problem is 2d (2d not implemented yet)");
227 }

◆ getValueL2()

MoFEMErrorCode EdgePolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 217 of file EdgePolynomialBase.cpp.

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

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 50 of file EdgePolynomialBase.cpp.

51  {
52 
54  *iface = NULL;
55  if (uuid == IDD_EDGE_BASE_FUNCTION) {
56  *iface = const_cast<EdgePolynomialBase *>(this);
58  } else {
59  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
60  }
61  ierr = BaseFunction::query_interface(uuid, iface);
62  CHKERRG(ierr);
64 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:522
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:565
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
Calculate base functions on tetrahedral.
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
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: