16                             boost::shared_ptr<BaseFunctionCtx> ctx_ptr) {
 
   21  int nb_gauss_pts = pts.size2();
 
   28        "Wrong dimension of pts, should be at least 3 rows with coordinates");
 
   35      data.
dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 2,
 
   39          &pts(0, 0), nb_gauss_pts);
 
   45        (
unsigned int)nb_gauss_pts)
 
   47               "Base functions or nodes has wrong number of integration points " 
   51    data.
dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 2 * 1,
 
   53    for (
auto gg = 0; gg != nb_gauss_pts; ++gg)
 
 
  102  const int nb_gauss_pts = pts.size2();
 
  105    auto &ptr = data.getBBAlphaIndicesSharedPtr(
field_name);
 
  119    auto &ptr = data.getBBDiffNSharedPtr(
field_name);
 
  126  auto &get_diff_n = get_diff_base(data.
dataOnEntities[MBVERTEX][0]);
 
  127  get_n.resize(nb_gauss_pts, 2, 
false);
 
  128  get_diff_n.resize(nb_gauss_pts, 2, 
false);
 
  133      (
unsigned int)nb_gauss_pts)
 
  135             "Base functions or nodes has wrong number of integration points " 
  141  vertex_alpha.resize(2, 2, 
false);
 
  142  vertex_alpha(0, 0) = data.
dataOnEntities[MBVERTEX][0].getBBNodeOrder()[0];
 
  143  vertex_alpha(0, 1) = 0;
 
  144  vertex_alpha(1, 0) = 0;
 
  145  vertex_alpha(1, 1) = data.
dataOnEntities[MBVERTEX][0].getBBNodeOrder()[1];
 
  148      1, nb_gauss_pts, vertex_alpha.size1(), &vertex_alpha(0, 0), &
lambda(0, 0),
 
  150  std::array<double, 2> f = {
 
  151      boost::math::factorial<double>(
 
  153      boost::math::factorial<double>(
 
  156  for (
int g = 0; 
g != nb_gauss_pts; ++
g)
 
  157    for (
int n = 0; 
n != 2; ++
n) {
 
  159      get_diff_n(
g, 
n) *= f[
n];
 
  165              "Wrong size ent of ent data");
 
  172    get_n.resize(nb_gauss_pts, nb_dofs, 
false);
 
  173    get_diff_n.resize(nb_gauss_pts, nb_dofs, 
false);
 
  177      edge_alpha.resize(nb_dofs, 2);
 
  180          order, nb_gauss_pts, edge_alpha.size1(), &edge_alpha(0, 0),
 
  185    get_n.resize(nb_gauss_pts, 0, 
false);
 
  186    get_diff_n.resize(nb_gauss_pts, 0, 
false);
 
 
  197  PetscErrorCode (*base_polynomials)(
int p, 
double s, 
double *diff_s, 
double *
L,
 
  198                                     double *
diffL, 
const int dim) =
 
  201  int nb_gauss_pts = pts.size2();
 
  203  const int side_number = 0;
 
  220    for (
int gg = 0; gg != nb_gauss_pts; gg++) {
 
  222      double s = 2 * pts(0, gg) - 1; 
 
  230                              &*
L.data().begin(), &*
diffL.data().begin(), 1);
 
  232      for (
unsigned int pp = 0;
 
 
  262  int nb_gauss_pts = pts.size2();
 
  264  const int side_number = 0;
 
  275  double diff_shape[] = {-1, 1};
 
  279    for (
int gg = 0; gg != nb_gauss_pts; gg++) {
 
  280      const double mu0 = 1.0 - pts(0, gg); 
 
  281      const double mu1 = pts(0, gg);
 
  287        order, &*shape.data().begin(), diff_shape,
 
  288        &*data.
dataOnEntities[MBEDGE][side_number].getN(base).data().begin(),
 
 
  321  PetscErrorCode (*base_polynomials)(
int p, 
double s, 
double *diff_s, 
double *
L,
 
  322                                     double *
diffL, 
const int dim) =
 
  328  int nb_gauss_pts = pts.size2();
 
  330  constexpr int side_number = 0;
 
  339      &*data.
dataOnEntities[MBEDGE][side_number].getN(base).data().begin();
 
  341      &*data.
dataOnEntities[MBEDGE][side_number].getDiffN(base).data().begin();
 
  348    for (
int gg = 0; gg != nb_gauss_pts; gg++) {
 
  349      double mu = 2 * pts(0, gg) - 1;
 
  353        fun_n[qd_shift + 
n] = 
l[
n];
 
  354        diff_fun_n[qd_shift + 
n] = diff_l[
n];
 
 
  369          "Make no sense, unless problem is 2d (2d not implemented yet)");
 
 
  379  PetscErrorCode (*base_polynomials)(
int p, 
double s, 
double *diff_s, 
double *
L,
 
  380                                     double *
diffL, 
const int dim) =
 
  383  int nb_gauss_pts = pts.size2();
 
  393    data.
dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
 
  395    data.
dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
 
  401        &*data.
dataOnEntities[MBEDGE][0].getN(base).data().begin(), 
nullptr,
 
  402        nb_gauss_pts, base_polynomials);
 
  405    data.
dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, 
false);
 
  406    data.
dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
 
 
  420  int nb_gauss_pts = pts.size2();
 
  425              "No data structure to store base functions");
 
  431    data.
dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
 
  433    data.
dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
 
  444    data.
dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0, 
false);
 
  445    data.
dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
 
 
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
#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 MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
PetscErrorCode Lobatto_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Lobatto basis functions .
#define NBEDGE_DEMKOWICZ_HCURL(P)
#define NBEDGE_H1(P)
Number of base function on edge for H1 space.
#define NBEDGE_L2(P)
Number of base functions on edge from L2 space.
#define NBEDGE_AINSWORTH_HCURL(P)
const double v
phase velocity of light in medium (cm/ns)
const double n
refractive index of diffusive medium
FTensor::Index< 'l', 3 > l
MoFEMErrorCode H1_BubbleShapeFunctions_ONSEGMENT(int p, double *N, double *diffN, double *bubbleN, double *diff_bubbleN, int nb_integration_pts)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasMatrix< double > MatrixDouble
UBlasMatrix< int > MatrixInt
implementation of Data Operators for Forces and Sources
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.
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.
constexpr auto field_name
static MoFEMErrorCode baseFunctionsEdge(const int N, const int gdim, const int n_alpha, const int *alpha, const double *lambda, const double *grad_lambda, double *base, double *grad_base)
static MoFEMErrorCode generateIndicesEdgeEdge(const int N, int *alpha)
Calculate base functions on tetrahedral.
MoFEMErrorCode getValueH1AinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValue(MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
MoFEMErrorCode getValueH1DemkowiczBase(MatrixDouble &pts)
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
MoFEMErrorCode getValueL2AinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode getValueH1BernsteinBezierBase(MatrixDouble &pts)
MoFEMErrorCode getValueL2(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueH1(MatrixDouble &pts)
MoFEMErrorCode getValueL2DemkowiczBase(MatrixDouble &pts)
Class used to pass element data to calculate base functions on tet,triangle,edge.
const FieldApproximationBase copyNodeBase
const std::string fieldName
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
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
base class for all interface classes
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.