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

◆ getValueH1()

MoFEMErrorCode EdgePolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 102 of file EdgePolynomialBase.cpp.

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

◆ getValueHcurl()

MoFEMErrorCode EdgePolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 269 of file EdgePolynomialBase.cpp.

269  {
271 
272  switch (cTx->bAse) {
276  break;
279  break;
280  default:
281  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
282  }
283 
285 }
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:145
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
EntPolynomialBaseCtx * cTx

◆ getValueHcurlAinsworthBase()

MoFEMErrorCode EdgePolynomialBase::getValueHcurlAinsworthBase ( MatrixDouble pts)
private

Definition at line 192 of file EdgePolynomialBase.cpp.

192  {
194 
196  const FieldApproximationBase base = cTx->bAse;
197 
198  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
199  double *diffL, const int dim) =
201 
202  int nb_gauss_pts = pts.size2();
203  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
204 
205  if (data.dataOnEntities[MBEDGE].size() != 1)
206  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
207 
208  int sense = data.dataOnEntities[MBEDGE][0].getSense();
209  int order = data.dataOnEntities[MBEDGE][0].getDataOrder();
210  int nb_dofs =
211  NBEDGE_AINSWORTH_HCURL(data.dataOnEntities[MBEDGE][0].getDataOrder());
212  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
213  false);
214  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
215  false);
216 
218  sense, order, &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0),
219  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
220  &*data.dataOnEntities[MBEDGE][0].getN(base).data().begin(), nullptr,
221  nb_gauss_pts, base_polynomials);
222 
223  } else {
224  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, false);
225  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
226  false);
227  }
228 
230 }
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,...
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
FieldApproximationBase
approximation base
Definition: definitions.h:143
DataForcesAndSourcesCore & dAta
#define NBEDGE_AINSWORTH_HCURL(P)
field with continuous tangents
Definition: definitions.h:172
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
EntPolynomialBaseCtx * cTx
const int order
approximation order
MoFEMErrorCode Hcurl_Ainsworth_EdgeBaseFunctions_MBTET_ON_EDGE(int sense, int p, double *N, double *diffN, double *edgeN, double *diff_edgeN, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Edge based H-curl base functions on edge.
Definition: Hcurl.cpp:175

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode EdgePolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 233 of file EdgePolynomialBase.cpp.

233  {
235 
237  const FieldApproximationBase base = cTx->bAse;
238 
239  int nb_gauss_pts = pts.size2();
240  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
241 
242  if (data.dataOnEntities[MBEDGE].size() != 1)
243  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
244  "No data structure to store base functions");
245 
246  int sense = data.dataOnEntities[MBEDGE][0].getSense();
247  int order = data.dataOnEntities[MBEDGE][0].getDataOrder();
248  int nb_dofs =
249  NBEDGE_DEMKOWICZ_HCURL(data.dataOnEntities[MBEDGE][0].getDataOrder());
250  data.dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
251  false);
252  data.dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
253  false);
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);
259  } else {
260 
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
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:477
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
FieldApproximationBase
approximation base
Definition: definitions.h:143
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:2143
field with continuous tangents
Definition: definitions.h:172
#define CHKERR
Inline error check.
Definition: definitions.h:596
#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:407
EntPolynomialBaseCtx * cTx
const int order
approximation order

◆ getValueHdiv()

MoFEMErrorCode EdgePolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 186 of file EdgePolynomialBase.cpp.

186  {
187  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
188  "Make no sense, unless problem is 2d (2d not implemented yet)");
189 }

◆ getValueL2()

MoFEMErrorCode EdgePolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 179 of file EdgePolynomialBase.cpp.

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

◆ query_interface()

MoFEMErrorCode EdgePolynomialBase::query_interface ( const MOFEMuuid uuid,
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:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
static const MOFEMuuid IDD_EDGE_BASE_FUNCTION
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: