v0.13.0
Public Member Functions | Private Member Functions | Private Attributes | List of all members
MoFEM::HexPolynomialBase Struct Reference

Calculate base functions on tetrahedral. More...

#include <src/approximation/HexPolynomialBase.hpp>

Inheritance diagram for MoFEM::HexPolynomialBase:
[legend]
Collaboration diagram for MoFEM::HexPolynomialBase:
[legend]

Public Member Functions

MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 HexPolynomialBase ()=default
 
 ~HexPolynomialBase ()=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. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 

Private Member Functions

MoFEMErrorCode getValueH1 (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2 (MatrixDouble &pts)
 Get base functions for L2 space. More...
 
MoFEMErrorCode getValueHdiv (MatrixDouble &pts)
 Get base functions for Hdiv space. More...
 
MoFEMErrorCode getValueHcurl (MatrixDouble &pts)
 Get base functions for Hcurl space. More...
 
MoFEMErrorCode getValueH1DemkowiczBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2DemkowiczBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdivDemkowiczBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurlDemkowiczBase (MatrixDouble &pts)
 

Private Attributes

EntPolynomialBaseCtxcTx
 
std::array< MatrixDouble, 6 > faceFamily
 
std::array< MatrixDouble, 6 > diffFaceFamily
 
MatrixDouble volFamily
 
MatrixDouble diffVolFamily
 

Additional Inherited Members

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

Detailed Description

Calculate base functions on tetrahedral.

Definition at line 31 of file HexPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ HexPolynomialBase()

MoFEM::HexPolynomialBase::HexPolynomialBase ( )
default

◆ ~HexPolynomialBase()

MoFEM::HexPolynomialBase::~HexPolynomialBase ( )
default

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 723 of file HexPolynomialBase.cpp.

724  {
726 
727  cTx = ctx_ptr->getInterface<EntPolynomialBaseCtx>();
728 
729  int nb_gauss_pts = pts.size2();
730  if (!nb_gauss_pts)
732 
733  if (pts.size1() < 3)
734  SETERRQ(
735  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
736  "Wrong dimension of pts, should be at least 3 rows with coordinates");
737 
738  switch (cTx->sPace) {
739  case H1:
740  CHKERR getValueH1(pts);
741  break;
742  case HDIV:
743  CHKERR getValueHdiv(pts);
744  break;
745  case HCURL:
746  CHKERR getValueHcurl(pts);
747  break;
748  case L2:
749  CHKERR getValueL2(pts);
750  break;
751  default:
752  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Unknown space");
753  }
754 
756 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
@ L2
field with C-1 continuity
Definition: definitions.h:101
@ H1
continuous field
Definition: definitions.h:98
@ HCURL
field with continuous tangents
Definition: definitions.h:99
@ HDIV
field with continuous normal traction
Definition: definitions.h:100
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
Class used to pass element data to calculate base functions on tet,triangle,edge.
MoFEMErrorCode getValueL2(MatrixDouble &pts)
Get base functions for L2 space.
MoFEMErrorCode getValueH1(MatrixDouble &pts)
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
Get base functions for Hcurl space.
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
Get base functions for Hdiv space.
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ getValueH1()

MoFEMErrorCode HexPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 32 of file HexPolynomialBase.cpp.

32  {
34 
35  switch (cTx->bAse) {
38  break;
39  default:
40  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
41  }
42 
44 }
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:79
const FieldApproximationBase bAse
MoFEMErrorCode getValueH1DemkowiczBase(MatrixDouble &pts)

◆ getValueH1DemkowiczBase()

MoFEMErrorCode HexPolynomialBase::getValueH1DemkowiczBase ( MatrixDouble pts)
private

Definition at line 46 of file HexPolynomialBase.cpp.

46  {
48 
49  auto &data = cTx->dAta;
50  const auto base = cTx->bAse;
51  const auto copy_base = cTx->copyNodeBase;
52  int nb_gauss_pts = pts.size2();
53 
54  auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
55  auto &copy_diff_base_fun =
56  data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
57  if (copy_base_fun.size1() != nb_gauss_pts)
58  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
59  "Inconsistent number of integration pts");
60 
61  auto add_base_on_verts = [&] {
63  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts,
64  copy_base_fun.size2());
65  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(
66  nb_gauss_pts, copy_diff_base_fun.size2());
67  noalias(data.dataOnEntities[MBVERTEX][0].getN(base)) = copy_base_fun;
68  noalias(data.dataOnEntities[MBVERTEX][0].getDiffN(base)) =
69  copy_diff_base_fun;
71  };
72 
73  // Edges
74  auto add_base_on_edges = [&] {
76  std::array<int, 12> sense;
77  std::array<int, 12> order;
78  if (data.spacesOnEntities[MBEDGE].test(H1)) {
79  if (data.dataOnEntities[MBEDGE].size() != 12)
80  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
81  "Expected 12 data on entities");
82 
83  std::array<double *, 12> h1_edge_n;
84  std::array<double *, 12> diff_h1_egde_n;
85  bool nb_dofs_sum = false;
86  for (int ee = 0; ee != 12; ++ee) {
87  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0)
88  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
89  "Sense of entity not set");
90 
91  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
92  order[ee] = data.dataOnEntities[MBEDGE][ee].getOrder();
93 
94  int nb_dofs = NBEDGE_H1(data.dataOnEntities[MBEDGE][ee].getOrder());
95  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
96  false);
97  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
98  nb_gauss_pts, 3 * nb_dofs, false);
99  h1_edge_n[ee] =
100  &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
101  diff_h1_egde_n[ee] =
102  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
103 
104  nb_dofs_sum |= (nb_dofs > 0);
105  }
106  if (nb_dofs_sum) {
108  sense.data(), order.data(), &*copy_base_fun.data().begin(),
109  &*copy_diff_base_fun.data().begin(), h1_edge_n.data(),
110  diff_h1_egde_n.data(), nb_gauss_pts);
111  }
112  } else {
113  for (int ee = 0; ee != 12; ++ee) {
114  data.dataOnEntities[MBEDGE][ee].getN(base).resize(0, 0, false);
115  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(0, 0, false);
116  }
117  }
119  };
120 
121  // Face
122  auto add_base_on_quads = [&]() {
124  std::array<int, 6> order;
125  if (data.spacesOnEntities[MBQUAD].test(H1)) {
126  // faces
127  if (data.dataOnEntities[MBQUAD].size() != 6)
128  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
129  "Expected six faces on hex");
130 
131  std::array<double *, 6> h1_face_n;
132  std::array<double *, 6> diff_h1_face_n;
133  bool nb_dofs_sum = false;
134  for (int ff = 0; ff != 6; ++ff) {
135 
136  order[ff] = data.dataOnEntities[MBQUAD][ff].getOrder();
137  const int nb_dofs = NBFACEQUAD_H1(order[ff]);
138 
139  data.dataOnEntities[MBQUAD][ff].getN(base).resize(nb_gauss_pts, nb_dofs,
140  false);
141  data.dataOnEntities[MBQUAD][ff].getDiffN(base).resize(
142  nb_gauss_pts, 3 * nb_dofs, false);
143 
144  h1_face_n[ff] =
145  &*data.dataOnEntities[MBQUAD][ff].getN(base).data().begin();
146  diff_h1_face_n[ff] =
147  &*data.dataOnEntities[MBQUAD][ff].getDiffN(base).data().begin();
148 
149  nb_dofs_sum |= (nb_dofs > 0);
150  }
151  if (data.facesNodes.size1() != 6)
152  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
153  "Expected six face nodes");
154  if (data.facesNodes.size2() != 4)
155  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
156  "Expected four nodes on face");
157 
158  if (nb_dofs_sum) {
160  &*data.facesNodes.data().begin(),
161  &*data.facesNodesOrder.data().begin(), order.data(),
162  &*copy_base_fun.data().begin(), &*copy_diff_base_fun.data().begin(),
163  h1_face_n.data(), diff_h1_face_n.data(), nb_gauss_pts);
164  }
165 
166  } else {
167  for (int ff = 0; ff != 6; ++ff) {
168  data.dataOnEntities[MBQUAD][ff].getN(base).resize(0, false);
169  data.dataOnEntities[MBQUAD][ff].getDiffN(base).resize(0, 0, false);
170  }
171  }
172 
174  };
175 
176  // Face
177  auto add_base_on_volume = [&]() {
179 
180  if (data.spacesOnEntities[MBHEX].test(H1)) {
181  // volume
182  int order = data.dataOnEntities[MBHEX][0].getOrder();
183  int nb_vol_dofs = NBVOLUMEHEX_H1(order);
184  data.dataOnEntities[MBHEX][0].getN(base).resize(nb_gauss_pts, nb_vol_dofs,
185  false);
186  data.dataOnEntities[MBHEX][0].getDiffN(base).resize(
187  nb_gauss_pts, 3 * nb_vol_dofs, false);
188 
189  if (nb_vol_dofs) {
190  const std::array<int, 3> p{order, order, order};
192  p.data(), &*copy_base_fun.data().begin(),
193  &*copy_diff_base_fun.data().begin(),
194  &*data.dataOnEntities[MBHEX][0].getN(base).data().begin(),
195  &*data.dataOnEntities[MBHEX][0].getDiffN(base).data().begin(),
196  nb_gauss_pts);
197  }
198  } else {
199  data.dataOnEntities[MBHEX][0].getN(base).resize(0, 0, false);
200  data.dataOnEntities[MBHEX][0].getDiffN(base).resize(0, 0, false);
201  }
202 
204  };
205 
206  CHKERR add_base_on_verts();
207  CHKERR add_base_on_edges();
208  CHKERR add_base_on_quads();
209  CHKERR add_base_on_volume();
210 
212 }
static Index< 'p', 3 > p
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
#define NBFACEQUAD_H1(P)
Number of base functions on quad for H1 space.
#define NBVOLUMEHEX_H1(P)
Number of base functions on hex for H1 space.
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
MoFEMErrorCode H1_EdgeShapeFunctions_ONHEX(int *sense, int *p, double *N, double *N_diff, double *edgeN[12], double *diff_edgeN[12], int nb_integration_pts)
MoFEMErrorCode H1_FaceShapeFunctions_ONHEX(int *face_nodes, int *face_nodes_order, int *p, double *N, double *N_diff, double *faceN[6], double *diff_faceN[6], int nb_integration_pts)
MoFEMErrorCode H1_InteriorShapeFunctions_ONHEX(const int *p, double *N, double *N_diff, double *faceN, double *diff_faceN, int nb_integration_pts)
const FieldApproximationBase copyNodeBase

◆ getValueHcurl()

MoFEMErrorCode HexPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Get base functions for Hcurl space.

Parameters
ptsmatrix of intergation pts
Returns
MoFEMErrorCode
Note
matrix of integration points on rows has local coordinates of finite element on columns are integration pts.

Definition at line 448 of file HexPolynomialBase.cpp.

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

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode HexPolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 463 of file HexPolynomialBase.cpp.

463  {
465 
466  auto &data = cTx->dAta;
467  const auto base = cTx->bAse;
468  const auto copy_base = cTx->copyNodeBase;
469  const int nb_gauss_pts = pts.size2();
470 
471  auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
472  auto &copy_diff_base_fun =
473  data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
474 
475  if (nb_gauss_pts != copy_base_fun.size1())
476  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
477  "Wrong number of gauss pts");
478 
479  if (nb_gauss_pts != copy_diff_base_fun.size1())
480  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
481  "Wrong number of gauss pts");
482 
483  // edges
484  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
485  std::array<int, 12> sense;
486  std::array<int, 12> order;
487  if (data.dataOnEntities[MBEDGE].size() != 12)
488  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
489  "Expected 12 edges data structures on Hex");
490 
491  std::array<double *, 12> hcurl_edge_n;
492  std::array<double *, 12> diff_hcurl_edge_n;
493  bool sum_nb_dofs = false;
494 
495  for (int ee = 0; ee != 12; ++ee) {
496  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0)
497  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
498  "Sense on edge <%d> on Hex not set", ee);
499 
500  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
501  order[ee] = data.dataOnEntities[MBEDGE][ee].getOrder();
502  const int nb_dofs = NBEDGE_DEMKOWICZ_HCURL(order[ee]);
503  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
504  3 * nb_dofs, false);
505  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,
506  9 * nb_dofs, false);
507  hcurl_edge_n[ee] =
508  &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
509  diff_hcurl_edge_n[ee] =
510  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
511 
512  sum_nb_dofs |= (nb_dofs > 0);
513  }
514 
515  if (sum_nb_dofs) {
517  sense.data(), order.data(), &*copy_base_fun.data().begin(),
518  &*copy_diff_base_fun.data().begin(), hcurl_edge_n.data(),
519  diff_hcurl_edge_n.data(), nb_gauss_pts);
520  }
521 
522  } else {
523  for (int ee = 0; ee != 12; ee++) {
524  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, 0, false);
525  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts, 0,
526  false);
527  }
528  }
529 
530  // Quad
531  if (data.spacesOnEntities[MBQUAD].test(HCURL)) {
532 
533  if (data.dataOnEntities[MBQUAD].size() != 6)
534  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
535  "Expected six data structures on Hex");
536 
537  std::array<int, 6> order;
538  double *face_family_ptr[6][2];
539  double *diff_face_family_ptr[6][2];
540 
541  bool sum_nb_dofs = false;
542  for (int ff = 0; ff != 6; ff++) {
543  if (data.dataOnEntities[MBQUAD][ff].getSense() == 0)
544  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
545  "Sense pn quad <%d> not set", ff);
546 
547  order[ff] = data.dataOnEntities[MBQUAD][ff].getOrder();
548  if (data.facesNodes.size1() != 6)
549  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
550  "Expected six faces");
551  if (data.facesNodes.size2() != 4)
552  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
553  "Expected four nodes on face");
554 
555  const int nb_family_dofs =
557  faceFamily[ff].resize(2, 3 * nb_family_dofs * nb_gauss_pts, false);
558  diffFaceFamily[ff].resize(2, 9 * nb_family_dofs * nb_gauss_pts, false);
559 
560  if (nb_family_dofs) {
561  face_family_ptr[ff][0] = &(faceFamily[ff](0, 0));
562  face_family_ptr[ff][1] = &(faceFamily[ff](1, 0));
563  diff_face_family_ptr[ff][0] = &(diffFaceFamily[ff](0, 0));
564  diff_face_family_ptr[ff][1] = &(diffFaceFamily[ff](1, 0));
565  }
566 
567  sum_nb_dofs |= (nb_family_dofs > 0);
568  }
569 
570  if (sum_nb_dofs) {
572  &*data.facesNodes.data().begin(),
573  &*data.facesNodesOrder.data().begin(), order.data(),
574  &*copy_base_fun.data().begin(), &*copy_diff_base_fun.data().begin(),
575  face_family_ptr, diff_face_family_ptr, nb_gauss_pts);
576 
577  for (int ff = 0; ff != 6; ++ff) {
578 
579  // put family back
580 
581  int nb_dofs = NBFACEQUAD_DEMKOWICZ_HCURL(order[ff]);
582  if (nb_dofs) {
583  auto &face_n = data.dataOnEntities[MBQUAD][ff].getN(base);
584  auto &diff_face_n = data.dataOnEntities[MBQUAD][ff].getDiffN(base);
585  face_n.resize(nb_gauss_pts, 3 * nb_dofs, false);
586  diff_face_n.resize(nb_gauss_pts, 9 * nb_dofs, false);
587 
588  auto ptr_f0 = &(faceFamily[ff](0, 0));
589  auto ptr_f1 = &(faceFamily[ff](1, 0));
590  double *ptr = &face_n(0, 0);
591  for (int n = 0; n != faceFamily[ff].size2() / 3; ++n) {
592  for (int j = 0; j != 3; ++j) {
593  *ptr = *ptr_f0;
594  ++ptr;
595  ++ptr_f0;
596  }
597  for (int j = 0; j != 3; ++j) {
598  *ptr = *ptr_f1;
599  ++ptr;
600  ++ptr_f1;
601  }
602  }
603 
604  auto diff_ptr_f0 = &(diffFaceFamily[ff](0, 0));
605  auto diff_ptr_f1 = &(diffFaceFamily[ff](1, 0));
606  double *diff_ptr = &diff_face_n(0, 0);
607  for (int n = 0; n != diffFaceFamily[ff].size2() / 9; ++n) {
608  for (int j = 0; j != 9; ++j) {
609  *diff_ptr = *diff_ptr_f0;
610  ++diff_ptr;
611  ++diff_ptr_f0;
612  }
613  for (int j = 0; j != 9; ++j) {
614  *diff_ptr = *diff_ptr_f1;
615  ++diff_ptr;
616  ++diff_ptr_f1;
617  }
618  }
619  }
620  }
621  } else {
622  for (int ff = 0; ff != 6; ff++) {
623  data.dataOnEntities[MBQUAD][ff].getN(base).resize(nb_gauss_pts, 0,
624  false);
625  data.dataOnEntities[MBQUAD][ff].getDiffN(base).resize(nb_gauss_pts, 0,
626  false);
627  }
628  }
629 
630  } else {
631  for (int ff = 0; ff != 6; ff++) {
632  data.dataOnEntities[MBQUAD][ff].getN(base).resize(nb_gauss_pts, 0, false);
633  data.dataOnEntities[MBQUAD][ff].getDiffN(base).resize(nb_gauss_pts, 0,
634  false);
635  }
636  }
637 
638  // Hex
639  if (data.spacesOnEntities[MBHEX].test(HCURL)) {
640 
641  const int order = data.dataOnEntities[MBHEX][0].getOrder();
642  const int nb_dofs = NBVOLUMEHEX_DEMKOWICZ_FAMILY_HCURL(order, order, order);
643 
644  volFamily.resize(3, 3 * nb_dofs * nb_gauss_pts);
645  diffVolFamily.resize(3, 9 * nb_dofs * nb_gauss_pts);
646  if (nb_dofs) {
647 
648  std::array<double *, 3> family_ptr = {&volFamily(0, 0), &volFamily(1, 0),
649  &volFamily(2, 0)};
650  std::array<double *, 3> diff_family_ptr = {
651  &diffVolFamily(0, 0), &diffVolFamily(1, 0), &diffVolFamily(2, 0)};
652 
653  std::array<int, 3> p{order, order, order};
655  p.data(), &*copy_base_fun.data().begin(),
656  &*copy_diff_base_fun.data().begin(), family_ptr.data(),
657  diff_family_ptr.data(), nb_gauss_pts);
658 
659  const int nb_vol_dofs = NBVOLUMEHEX_DEMKOWICZ_HCURL(order);
660  auto &face_n = data.dataOnEntities[MBHEX][0].getN(base);
661  auto &diff_face_n = data.dataOnEntities[MBHEX][0].getDiffN(base);
662  face_n.resize(nb_gauss_pts, 3 * nb_vol_dofs, false);
663  diff_face_n.resize(nb_gauss_pts, 9 * nb_vol_dofs, false);
664 
665  auto ptr_f0 = &(volFamily(0, 0));
666  auto ptr_f1 = &(volFamily(1, 0));
667  auto ptr_f2 = &(volFamily(2, 0));
668  double *ptr = &face_n(0, 0);
669  for (int n = 0; n != volFamily.size2() / 3; ++n) {
670  for (int j = 0; j != 3; ++j) {
671  *ptr = *ptr_f0;
672  ++ptr;
673  ++ptr_f0;
674  }
675  for (int j = 0; j != 3; ++j) {
676  *ptr = *ptr_f1;
677  ++ptr;
678  ++ptr_f1;
679  }
680  for (int j = 0; j != 3; ++j) {
681  *ptr = *ptr_f2;
682  ++ptr;
683  ++ptr_f2;
684  }
685  }
686 
687  auto diff_ptr_f0 = &(diffVolFamily(0, 0));
688  auto diff_ptr_f1 = &(diffVolFamily(1, 0));
689  auto diff_ptr_f2 = &(diffVolFamily(2, 0));
690  double *diff_ptr = &diff_face_n(0, 0);
691  for (int n = 0; n != diffVolFamily.size2() / 9; ++n) {
692  for (int j = 0; j != 9; ++j) {
693  *diff_ptr = *diff_ptr_f0;
694  ++diff_ptr;
695  ++diff_ptr_f0;
696  }
697  for (int j = 0; j != 9; ++j) {
698  *diff_ptr = *diff_ptr_f1;
699  ++diff_ptr;
700  ++diff_ptr_f1;
701  }
702  for (int j = 0; j != 9; ++j) {
703  *diff_ptr = *diff_ptr_f2;
704  ++diff_ptr;
705  ++diff_ptr_f2;
706  }
707  }
708  } else {
709  data.dataOnEntities[MBHEX][0].getN(base).resize(nb_gauss_pts, 0, false);
710  data.dataOnEntities[MBHEX][0].getDiffN(base).resize(nb_gauss_pts, 0,
711  false);
712  }
713 
714  } else {
715  data.dataOnEntities[MBHEX][0].getN(base).resize(nb_gauss_pts, 0, false);
716  data.dataOnEntities[MBHEX][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
717  }
718 
720 }
static Index< 'n', 3 > n
#define NBFACEQUAD_DEMKOWICZ_FAMILY_HCURL(P, Q)
Number of base functions on quad for Hcurl space.
#define NBEDGE_DEMKOWICZ_HCURL(P)
#define NBVOLUMEHEX_DEMKOWICZ_HCURL(P)
#define NBFACEQUAD_DEMKOWICZ_HCURL(P)
#define NBVOLUMEHEX_DEMKOWICZ_FAMILY_HCURL(P, Q, R)
FTensor::Index< 'j', 3 > j
MoFEMErrorCode Hcurl_EdgeShapeFunctions_ONHEX(int *sense, int *p, double *N, double *N_diff, double *edgeN[12], double *diff_edgeN[12], int nb_integration_pts)
MoFEMErrorCode Hcurl_FaceShapeFunctions_ONHEX(int *face_nodes, int *face_nodes_order, int *p, double *N, double *N_diff, double *faceN[6][2], double *diff_faceN[6][2], int nb_integration_pts)
MoFEMErrorCode Hcurl_InteriorShapeFunctions_ONHEX(int *p, double *N, double *N_diff, double *volN[3], double *diff_volN[3], int nb_integration_pts)
std::array< MatrixDouble, 6 > faceFamily
std::array< MatrixDouble, 6 > diffFaceFamily

◆ getValueHdiv()

MoFEMErrorCode HexPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Get base functions for Hdiv space.

Parameters
ptsmatrix of intergation pts
Returns
MoFEMErrorCode
Note
matrix of integration points on rows has local coordinates of finite element on columns are integration pts.

Definition at line 268 of file HexPolynomialBase.cpp.

268  {
270 
271  switch (cTx->bAse) {
274  break;
275  default:
276  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
277  }
278 
280 }
MoFEMErrorCode getValueHdivDemkowiczBase(MatrixDouble &pts)

◆ getValueHdivDemkowiczBase()

MoFEMErrorCode HexPolynomialBase::getValueHdivDemkowiczBase ( MatrixDouble pts)
private

Definition at line 282 of file HexPolynomialBase.cpp.

282  {
284 
285  auto &data = cTx->dAta;
286  const auto base = cTx->bAse;
287  const auto copy_base = cTx->copyNodeBase;
288  const int nb_gauss_pts = pts.size2();
289 
290  auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
291  auto &copy_diff_base_fun =
292  data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
293 
294  if (nb_gauss_pts != copy_base_fun.size1())
295  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
296  "Wrong number of gauss pts");
297 
298  if (nb_gauss_pts != copy_diff_base_fun.size1())
299  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
300  "Wrong number of gauss pts");
301 
302  // Quad
303  if (data.spacesOnEntities[MBQUAD].test(HDIV)) {
304 
305  if (data.dataOnEntities[MBQUAD].size() != 6)
306  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
307  "Expected six data structures on Hex");
308 
309  std::array<int, 6> order;
310  std::array<double *, 6> hdiv_face_n;
311  std::array<double *, 6> diff_hdiv_face_n;
312 
313  bool sum_nb_dofs = false;
314  for (int ff = 0; ff != 6; ff++) {
315  if (data.dataOnEntities[MBQUAD][ff].getSense() == 0)
316  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
317  "Sense pn quad <%d> not set", ff);
318 
319  order[ff] = data.dataOnEntities[MBQUAD][ff].getOrder();
320  if (data.facesNodes.size1() != 6)
321  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
322  "Expected six faces");
323  if (data.facesNodes.size2() != 4)
324  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
325  "Expected four nodes on face");
326 
327  const int nb_dofs = NBFACEQUAD_DEMKOWICZ_HDIV(order[ff]);
328  auto &face_n = data.dataOnEntities[MBQUAD][ff].getN(base);
329  auto &diff_face_n = data.dataOnEntities[MBQUAD][ff].getDiffN(base);
330  face_n.resize(nb_gauss_pts, 3 * nb_dofs, false);
331  diff_face_n.resize(nb_gauss_pts, 9 * nb_dofs, false);
332 
333  hdiv_face_n[ff] = &*face_n.data().begin();
334  diff_hdiv_face_n[ff] = &*diff_face_n.data().begin();
335 
336  sum_nb_dofs |= (nb_dofs > 0);
337  }
338 
339  if (sum_nb_dofs) {
341  &*data.facesNodes.data().begin(),
342  &*data.facesNodesOrder.data().begin(), order.data(),
343  &*copy_base_fun.data().begin(), &*copy_diff_base_fun.data().begin(),
344  hdiv_face_n.data(), diff_hdiv_face_n.data(), nb_gauss_pts);
345 
346  } else {
347  for (int ff = 0; ff != 6; ff++) {
348  data.dataOnEntities[MBQUAD][ff].getN(base).resize(nb_gauss_pts, 0,
349  false);
350  data.dataOnEntities[MBQUAD][ff].getDiffN(base).resize(nb_gauss_pts, 0,
351  false);
352  }
353  }
354 
355  } else {
356  for (int ff = 0; ff != 6; ff++) {
357  data.dataOnEntities[MBQUAD][ff].getN(base).resize(nb_gauss_pts, 0, false);
358  data.dataOnEntities[MBQUAD][ff].getDiffN(base).resize(nb_gauss_pts, 0,
359  false);
360  }
361  }
362 
363  // Hex
364  if (data.spacesOnEntities[MBHEX].test(HDIV)) {
365 
366  const int order = data.dataOnEntities[MBHEX][0].getOrder();
367  const int nb_dofs_family =
369 
370  volFamily.resize(3, 3 * nb_dofs_family * nb_gauss_pts);
371  diffVolFamily.resize(3, 9 * nb_dofs_family * nb_gauss_pts);
372  if (nb_dofs_family) {
373 
374  std::array<double *, 3> family_ptr = {&volFamily(0, 0), &volFamily(1, 0),
375  &volFamily(2, 0)};
376  std::array<double *, 3> diff_family_ptr = {
377  &diffVolFamily(0, 0), &diffVolFamily(1, 0), &diffVolFamily(2, 0)};
378 
379  std::array<int, 3> p{order, order, order};
381  p.data(), &*copy_base_fun.data().begin(),
382  &*copy_diff_base_fun.data().begin(), family_ptr.data(),
383  diff_family_ptr.data(), nb_gauss_pts);
384 
385  const int nb_vol_dofs = NBVOLUMEHEX_DEMKOWICZ_HDIV(order);
386  auto &face_n = data.dataOnEntities[MBHEX][0].getN(base);
387  auto &diff_face_n = data.dataOnEntities[MBHEX][0].getDiffN(base);
388  face_n.resize(nb_gauss_pts, 3 * nb_vol_dofs, false);
389  diff_face_n.resize(nb_gauss_pts, 9 * nb_vol_dofs, false);
390 
391  auto ptr_f0 = &(volFamily(0, 0));
392  auto ptr_f1 = &(volFamily(1, 0));
393  auto ptr_f2 = &(volFamily(2, 0));
394  double *ptr = &face_n(0, 0);
395  for (int n = 0; n != volFamily.size2() / 3; ++n) {
396  for (int j = 0; j != 3; ++j) {
397  *ptr = *ptr_f0;
398  ++ptr;
399  ++ptr_f0;
400  }
401  for (int j = 0; j != 3; ++j) {
402  *ptr = *ptr_f1;
403  ++ptr;
404  ++ptr_f1;
405  }
406  for (int j = 0; j != 3; ++j) {
407  *ptr = *ptr_f2;
408  ++ptr;
409  ++ptr_f2;
410  }
411  }
412 
413  auto diff_ptr_f0 = &(diffVolFamily(0, 0));
414  auto diff_ptr_f1 = &(diffVolFamily(1, 0));
415  auto diff_ptr_f2 = &(diffVolFamily(2, 0));
416  double *diff_ptr = &diff_face_n(0, 0);
417  for (int n = 0; n != diffVolFamily.size2() / 9; ++n) {
418  for (int j = 0; j != 9; ++j) {
419  *diff_ptr = *diff_ptr_f0;
420  ++diff_ptr;
421  ++diff_ptr_f0;
422  }
423  for (int j = 0; j != 9; ++j) {
424  *diff_ptr = *diff_ptr_f1;
425  ++diff_ptr;
426  ++diff_ptr_f1;
427  }
428  for (int j = 0; j != 9; ++j) {
429  *diff_ptr = *diff_ptr_f2;
430  ++diff_ptr;
431  ++diff_ptr_f2;
432  }
433  }
434  } else {
435  data.dataOnEntities[MBHEX][0].getN(base).resize(nb_gauss_pts, 0, false);
436  data.dataOnEntities[MBHEX][0].getDiffN(base).resize(nb_gauss_pts, 0,
437  false);
438  }
439 
440  } else {
441  data.dataOnEntities[MBHEX][0].getN(base).resize(nb_gauss_pts, 0, false);
442  data.dataOnEntities[MBHEX][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
443  }
444 
446 }
#define NBFACEQUAD_DEMKOWICZ_HDIV(P)
#define NBVOLUMEHEX_DEMKOWICZ_HDIV(P)
#define NBVOLUMEHEX_DEMKOWICZ_FAMILY_HDIV(P, Q, R)
MoFEMErrorCode Hdiv_FaceShapeFunctions_ONHEX(int *face_nodes, int *face_nodes_order, int *p, double *N, double *diffN, double *faceN[6], double *div_faceN[6], int nb_integration_pts)
MoFEMErrorCode Hdiv_InteriorShapeFunctions_ONHEX(int *p, double *N, double *N_diff, double *bubleN[3], double *div_bubleN[3], int nb_integration_pts)

◆ getValueL2()

MoFEMErrorCode HexPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Get base functions for L2 space.

Parameters
ptsmatrix of intergation pts
Returns
MoFEMErrorCode
Note
matrix of integration points on rows has local coordinates of finite element on columns are integration pts.

Definition at line 214 of file HexPolynomialBase.cpp.

214  {
216 
217  switch (cTx->bAse) {
220  break;
221  default:
222  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
223  }
224 
226 }
MoFEMErrorCode getValueL2DemkowiczBase(MatrixDouble &pts)

◆ getValueL2DemkowiczBase()

MoFEMErrorCode HexPolynomialBase::getValueL2DemkowiczBase ( MatrixDouble pts)
private

Definition at line 228 of file HexPolynomialBase.cpp.

228  {
230 
231  auto &data = cTx->dAta;
232  const auto base = cTx->bAse;
233  const auto copy_base = cTx->copyNodeBase;
234 
235  int nb_gauss_pts = pts.size2();
236 
237  auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
238  auto &copy_diff_base_fun =
239  data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
240 
241  if (nb_gauss_pts != copy_base_fun.size1())
242  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
243  "Wrong number of gauss pts");
244 
245  if (nb_gauss_pts != copy_diff_base_fun.size1())
246  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
247  "Wrong number of gauss pts");
248 
249  auto &base_fun = data.dataOnEntities[MBHEX][0].getN(base);
250  auto &diff_base_fun = data.dataOnEntities[MBHEX][0].getDiffN(base);
251  const int vol_order = data.dataOnEntities[MBHEX][0].getOrder();
252 
253  const int nb_dofs = NBVOLUMEHEX_L2(vol_order);
254  base_fun.resize(nb_gauss_pts, nb_dofs, false);
255  diff_base_fun.resize(nb_gauss_pts, 3 * nb_dofs, false);
256 
257  if (nb_dofs) {
258  const std::array<int, 3> p{vol_order, vol_order, vol_order};
260  p.data(), &*copy_base_fun.data().begin(),
261  &*copy_diff_base_fun.data().begin(), &*base_fun.data().begin(),
262  &*diff_base_fun.data().begin(), nb_gauss_pts);
263  }
264 
266 }
#define NBVOLUMEHEX_L2(P)
Number of base functions on hexahedron for L2 space.
MoFEMErrorCode L2_InteriorShapeFunctions_ONHEX(const int *p, double *N, double *N_diff, double *volN, double *diff_volN, int nb_integration_pts)

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 25 of file HexPolynomialBase.cpp.

26  {
28  *iface = const_cast<HexPolynomialBase *>(this);
30 }
Calculate base functions on tetrahedral.

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* MoFEM::HexPolynomialBase::cTx
private

Definition at line 42 of file HexPolynomialBase.hpp.

◆ diffFaceFamily

std::array<MatrixDouble, 6> MoFEM::HexPolynomialBase::diffFaceFamily
private

Definition at line 95 of file HexPolynomialBase.hpp.

◆ diffVolFamily

MatrixDouble MoFEM::HexPolynomialBase::diffVolFamily
private

Definition at line 98 of file HexPolynomialBase.hpp.

◆ faceFamily

std::array<MatrixDouble, 6> MoFEM::HexPolynomialBase::faceFamily
private

Definition at line 94 of file HexPolynomialBase.hpp.

◆ volFamily

MatrixDouble MoFEM::HexPolynomialBase::volFamily
private

Definition at line 97 of file HexPolynomialBase.hpp.


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