50            "Polynomial type not set");
 
   51  PetscErrorCode (*base_polynomials)(
int p, 
double s, 
double *diff_s, 
double *
L,
 
   52                                     double *diffL, 
const int dim) =
 
   55  auto ©_base_fun = data.
dataOnEntities[MBVERTEX][0].getN(copy_base);
 
   56  auto ©_diff_base_fun =
 
   59  int nb_gauss_pts = pts.size2();
 
   66              "should be four edges on quad");
 
   68    int sense[4], 
order[4];
 
   69    double *H1edgeN[4], *diffH1edgeN[4];
 
   70    for (
int ee = 0; ee != 4; ++ee) {
 
   72      if (ent_dat.getSense() == 0)
 
   75      sense[ee] = ent_dat.getSense();
 
   76      order[ee] = ent_dat.getOrder();
 
   77      int nb_dofs = 
NBEDGE_H1(ent_dat.getOrder());
 
   78      ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, 
false);
 
   79      ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, 
false);
 
   80      H1edgeN[ee] = &*ent_dat.getN(base).data().begin();
 
   81      diffH1edgeN[ee] = &*ent_dat.getDiffN(base).data().begin();
 
   84        sense, 
order, &*copy_base_fun.data().begin(),
 
   85        &*copy_diff_base_fun.data().begin(), H1edgeN, diffH1edgeN, nb_gauss_pts,
 
   94              "should be one quad to store bubble base on quad");
 
   98    ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, 
false);
 
   99    ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, 
false);
 
  100    int face_nodes[] = {0, 1, 2, 3};
 
  102        face_nodes, ent_dat.getOrder(), &*vert_dat.getN(base).data().begin(),
 
  103        &*vert_dat.getDiffN(base).data().begin(),
 
  104        &*ent_dat.getN(base).data().begin(),
 
  105        &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts,
 
 
  119  int nb_gauss_pts = pts.size2();
 
  120  auto ©_base_fun = data.
dataOnEntities[MBVERTEX][0].getN(copy_base);
 
  121  auto ©_diff_base_fun =
 
  128              "should be four edges on quad");
 
  130    int sense[4], 
order[4];
 
  131    double *H1edgeN[4], *diffH1edgeN[4];
 
  132    for (
int ee = 0; ee != 4; ++ee) {
 
  134      if (ent_dat.getSense() == 0)
 
  137      sense[ee] = ent_dat.getSense();
 
  138      order[ee] = ent_dat.getOrder();
 
  139      int nb_dofs = 
NBEDGE_H1(ent_dat.getOrder());
 
  140      ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, 
false);
 
  141      ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, 
false);
 
  142      H1edgeN[ee] = &*ent_dat.getN(base).data().begin();
 
  143      diffH1edgeN[ee] = &*ent_dat.getDiffN(base).data().begin();
 
  146        sense, 
order, &*copy_base_fun.data().begin(),
 
  147        &*copy_diff_base_fun.data().begin(), H1edgeN, diffH1edgeN,
 
  156              "should be one quad to store bubble base on quad");
 
  160    int p = ent_dat.getOrder();
 
  161    int order[2] = {p, p};
 
  162    ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, 
false);
 
  163    ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, 
false);
 
  165    int face_nodes[] = {0, 1, 2, 3};
 
  168          face_nodes, 
order, &*copy_base_fun.data().begin(),
 
  169          &*copy_diff_base_fun.data().begin(),
 
  170          &*ent_dat.getN(base).data().begin(),
 
  171          &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts);
 
 
  187            "Ainsworth not implemented on quad for L2 base");
 
 
  207  int nb_gauss_pts = pts.size2();
 
  208  auto ©_base_fun = data.
dataOnEntities[MBVERTEX][0].getN(copy_base);
 
  209  auto ©_diff_base_fun =
 
  213  int p = ent_dat.getOrder();
 
  214  int order[2] = {p, p};
 
  216  ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, 
false);
 
  217  ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, 
false);
 
  220      order, &*copy_base_fun.data().begin(),
 
  221      &*copy_diff_base_fun.data().begin(), &*ent_dat.getN(base).data().begin(),
 
  222      &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts);
 
 
  234            "Ainsworth not implemented on quad for Hcurl base");
 
 
  255  int nb_gauss_pts = pts.size2();
 
  256  auto ©_base_fun = data.
dataOnEntities[MBVERTEX][0].getN(copy_base);
 
  257  auto ©_diff_base_fun =
 
  265              "wrong number of edges on quad, should be 4 but is %ld",
 
  268    int sense[4], 
order[4];
 
  269    double *hcurl_edge_n[4];
 
  270    double *diff_hcurl_edge_n[4];
 
  272    for (
int ee = 0; ee != 4; ++ee) {
 
  276                "orientation (sense) of edge is not set");
 
  286          nb_gauss_pts, 3 * 2 * nb_dofs, 
false);
 
  290      diff_hcurl_edge_n[ee] =
 
  294        sense, 
order, &*copy_base_fun.data().begin(),
 
  295        &*copy_diff_base_fun.data().begin(), hcurl_edge_n, diff_hcurl_edge_n,
 
  304              "No data struture to keep base functions on face");
 
  308    if (nb_dofs_family) {
 
  309      faceFamily.resize(2, 3 * nb_dofs_family * nb_gauss_pts, 
false);
 
  310      diffFaceFamily.resize(2, 6 * nb_dofs_family * nb_gauss_pts, 
false);
 
  312      int order[2] = {p, p};
 
  316      int face_nodes[] = {0, 1, 2, 3};
 
  318          face_nodes, 
order, &*copy_base_fun.data().begin(),
 
  319          &*copy_diff_base_fun.data().begin(), face_family_ptr,
 
  320          diff_face_family_ptr, nb_gauss_pts);
 
  327    auto &diff_face_n = data.
dataOnEntities[MBQUAD][0].getDiffN(base);
 
  328    face_n.resize(nb_gauss_pts, 3 * nb_dofs, 
false);
 
  329    diff_face_n.resize(nb_gauss_pts, 3 * 2 * nb_dofs, 
false);
 
  334      double *ptr = &face_n(0, 0);
 
  336        for (
int j = 0; 
j != 3; ++
j) {
 
  341        for (
int j = 0; 
j != 3; ++
j) {
 
  350      double *diff_ptr = &diff_face_n(0, 0);
 
  352        for (
int j = 0; 
j != 6; ++
j) {
 
  353          *diff_ptr = *diff_ptr_f0;
 
  357        for (
int j = 0; 
j != 6; ++
j) {
 
  358          *diff_ptr = *diff_ptr_f1;
 
 
  375            "Ainsworth not implemented on quad for Hdiv base");
 
 
  395  int nb_gauss_pts = pts.size2();
 
  397  auto ©_base_fun = data.
dataOnEntities[MBVERTEX][0].getN(copy_base);
 
  398  auto ©_diff_base_fun =
 
  406              "No data struture to keep base functions on face");
 
  411    auto &diff_face_n = data.
dataOnEntities[MBQUAD][0].getDiffN(base);
 
  412    face_n.resize(nb_gauss_pts, 3 * nb_dofs, 
false);
 
  413    diff_face_n.resize(nb_gauss_pts, 6 * nb_dofs, 
false);
 
  417      std::array<int, 2> 
order = {p, p};
 
  418      std::array<int, 6> face_nodes = {0, 1, 2, 3};
 
  420          face_nodes.data(), 
order.data(), &*copy_base_fun.data().begin(),
 
  421          &*copy_diff_base_fun.data().begin(), &*face_n.data().begin(),
 
  422          &*diff_face_n.data().begin(), nb_gauss_pts);
 
 
  431                             boost::shared_ptr<BaseFunctionCtx> ctx_ptr) {
 
  436  int nb_gauss_pts = pts.size2();
 
  443        "Wrong dimension of pts, should be at least 3 rows with coordinates");
 
  449            "Shape base has to be on NOBASE <%s>",
 
  456  if (base_shape.size1() != pts.size2())
 
  458            "Number of base functions integration points not equal number of " 
  459            "set integration point");
 
  460  if (base_shape.size2() != 4)
 
  462            "Number of shape functions should be four");
 
  463  if (diff_base.size1() != pts.size2())
 
  466        "Number of diff base functions integration points not equal number of " 
  467        "set integration point %ld != %ld",
 
  468        diff_base.size1(), pts.size2());
 
  469  if (diff_base.size2() != 8)
 
  471            "Number of shape functions should be four");
 
 
FieldApproximationBase
approximation base
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ L2
field with C-1 continuity
@ HCURL
field with continuous tangents
@ HDIV
field with continuous normal traction
@ CONTINUOUS
Regular field.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
static const char *const ApproximationBaseNames[]
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define NBFACEQUAD_DEMKOWICZ_HDIV(P)
#define NBFACEQUAD_DEMKOWICZ_FAMILY_HCURL(P, Q)
Number of base functions on quad for Hcurl space.
#define NBEDGE_DEMKOWICZ_HCURL(P)
#define NBFACEQUAD_H1(P)
Number of base functions on quad for H1 space.
#define NBFACEQUAD_L2(P)
Number of base functions on quad for L2 space.
#define NBFACEQUAD_DEMKOWICZ_HCURL(P)
PetscErrorCode H1_QuadShapeFunctions_MBQUAD(int *faces_nodes, int p, double *N, double *diffN, double *faceN, double *diff_faceN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
#define NBEDGE_H1(P)
Number of base function on edge for H1 space.
PetscErrorCode H1_EdgeShapeFunctions_MBQUAD(int *sense, int *p, double *N, double *diffN, double *edgeN[4], double *diff_edgeN[4], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
const double n
refractive index of diffusive medium
FTensor::Index< 'j', 3 > j
MoFEMErrorCode L2_FaceShapeFunctions_ONQUAD(int *p, double *N, double *diffN, double *face_buble, double *diff_face_bubble, int nb_integration_pts)
L2 Face base functions on Quad.
MoFEMErrorCode H1_FaceShapeFunctions_ONQUAD(int *face_nodes, int *p, double *N, double *diffN, double *faceN, double *diff_faceN, int nb_integration_pts)
H1 Face bubble functions on Quad.
MoFEMErrorCode Hcurl_EdgeShapeFunctions_ONQUAD(int *sense, int *p, double *N, double *diffN, double *edgeN[4], double *curl_edgeN[4], int nb_integration_pts)
MoFEMErrorCode Hcurl_FaceShapeFunctions_ONQUAD(int *face_nodes, int *p, double *N, double *diffN, double *faceN[], double *diff_faceN[], int nb_integration_pts)
MoFEMErrorCode H1_EdgeShapeFunctions_ONQUAD(int *sense, int *p, double *N, double *diffN, double *edgeN[4], double *edgeDiffN[4], int nb_integration_pts)
H1 Edge base functions on Quad.
MoFEMErrorCode Hdiv_FaceShapeFunctions_ONQUAD(int *face_nodes, int *p, double *N, double *diffN, double *faceN, double *diff_faceN, int nb_integration_pts)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
Class used to pass element data to calculate base functions on tet,triangle,edge.
const FieldApproximationBase copyNodeBase
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
const FieldContinuity spaceContinuity
const FieldApproximationBase bAse
data structure for finite element entity
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
Calculate base functions on triangle.
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
MoFEMErrorCode getValueL2DemkowiczBase(MatrixDouble &pts)
MoFEMErrorCode getValueH1AinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueH1(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode getValueL2(MatrixDouble &pts)
EntPolynomialBaseCtx * cTx
MatrixDouble diffFaceFamily
MoFEMErrorCode getValueH1DemkowiczBase(MatrixDouble &pts)
MoFEMErrorCode getValue(MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
MoFEMErrorCode getValueHdivDemkowiczBase(MatrixDouble &pts)
base class for all interface classes
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.