v0.8.17
Macros | Functions
h1_hdiv_hcurl_l2.h File Reference

Functions to approximate hierarchical spaces. More...

Go to the source code of this file.

Macros

#define DEPRECATED
 
#define NBVOLUMETET_L2(P)   ((P + 1) * (P + 2) * (P + 3) / 6)
 Number of base functions on tetrahedron for L2 space. More...
 
#define NBFACETRI_L2(P)   ((P + 1) * (P + 2) / 2)
 Number of base functions on triangle for L2 space. More...
 
#define NBEDGE_L2(P)   (P + 1)
 Number of base functions on edge fro L2 space. More...
 
#define NBEDGE_H1(P)   (((P) > 0) ? (P - 1) : 0)
 Numer of base function on edge for H1 space. More...
 
#define NBFACETRI_H1(P)   (((P) > 1) ? ((P - 2) * (P - 1) / 2) : 0)
 Number of base function on triangle for H1 space. More...
 
#define NBFACEQUAD_H1(P)   (((P) > 2) ? ((P - 3) * (P - 2) / 2) : 0)
 Number of base functions on quad for H1 space. More...
 
#define NBVOLUMETET_H1(P)   (((P) > 2) ? ((P - 3) * (P - 2) * (P - 1) / 6) : 0)
 Number of base functions on tetrahedron fro H1 space. More...
 
#define NBVOLUMEPRISM_H1(P)   (((P) > 4) ? ((P - 5) * (P - 4) * (P - 3) / 6) : 0)
 Number of base functions on prism for H1 space. More...
 
#define NBEDGE_AINSWORTH_HCURL(P)   (((P) > 0) ? (P + 1) : 0)
 
#define NBFACETRI_AINSWORTH_EDGE_HCURL(P)   (((P) > 1) ? P - 1 : 0)
 
#define NBFACETRI_AINSWORTH_FACE_HCURL(P)   (((P) > 2) ? (P - 1) * (P - 2) : 0)
 
#define NBFACETRI_AINSWORTH_HCURL(P)   ((P > 1) ? ((P)-1) * (P + 1) : 0)
 
#define NBVOLUMETET_AINSWORTH_FACE_HCURL(P)   (((P) > 2) ? (2 * (P - 1) * (P - 2)) : 0)
 
#define NBVOLUMETET_AINSWORTH_TET_HCURL(P)   (((P) > 3) ? ((P - 3) * (P - 2) * (P - 1) / 2) : 0)
 
#define NBVOLUMETET_AINSWORTH_HCURL(P)   (((P) > 2) ? (P - 2) * (P - 1) * (P + 1) / 2 : 0)
 
#define NBEDGE_DEMKOWICZ_HCURL(P)   (((P) > 0) ? (P) : 0)
 
#define NBFACETRI_DEMKOWICZ_HCURL(P)   (((P) > 1) ? (P) * ((P)-1) : 0)
 
#define NBVOLUMETET_DEMKOWICZ_HCURL(P)   (((P) > 2) ? ((P) * ((P)-1) * ((P)-2) / 2) : 0)
 
#define NBEDGE_HDIV(P)   (0)
 
#define NBFACETRI_AINSWORTH_EDGE_HDIV(P)   (((P) > 0) ? (P) : 0)
 
#define NBFACETRI_AINSWORTH_FACE_HDIV(P)   (((P) > 2) ? (P - 1) * (P - 2) / 2 : 0)
 
#define NBFACETRI_AINSWORTH_HDIV(P)   (((P) > 0) ? (P + 1) * (P + 2) / 2 : 0)
 
#define NBVOLUMETET_AINSWORTH_EDGE_HDIV(P)   (((P) > 1) ? (P - 1) : 0)
 
#define NBVOLUMETET_AINSWORTH_FACE_HDIV(P)   (((P) > 2) ? (P - 1) * (P - 2) : 0)
 
#define NBVOLUMETET_AINSWORTH_VOLUME_HDIV(P)   (((P) > 3) ? (P - 3) * (P - 2) * (P - 1) / 2 : 0)
 
#define NBVOLUMETET_AINSWORTH_HDIV(P)   (((P) > 1) ? (P - 1) * (P + 1) * (P + 2) / 2 : 0)
 
#define NBFACETRI_DEMKOWICZ_HDIV(P)   ((P > 0) ? (P) * (P + 1) / 2 : 0)
 
#define NBVOLUMETET_DEMKOWICZ_HDIV(P)   (((P) > 1) ? (P) * (P - 1) * (P + 1) / 2 : 0)
 

Functions

PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTRI (int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Get base functions on triangle for L2 space. More...
 
PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTET (int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 Get base functions on tetrahedron for L2 space. More...
 
DEPRECATED PetscErrorCode L2_ShapeFunctions_MBTRI (int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 
DEPRECATED PetscErrorCode L2_ShapeFunctions_MBTET (int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 
PetscErrorCode L2_VolumeShapeDiffMBTETinvJ (int base_p, int p, double *volume_diffN, double *invJac, double *volume_diffNinvJac, int GDIM)
 
PetscErrorCode H1_EdgeShapeFunctions_MBTRI (int *sense, int *p, double *N, double *diffN, double *edgeN[3], double *diff_edgeN[3], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 H1_EdgeShapeFunctions_MBTRI. More...
 
PetscErrorCode H1_FaceShapeFunctions_MBTRI (const int *face_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))
 
PetscErrorCode H1_EdgeShapeFunctions_MBTET (int *sense, int *p, double *N, double *diffN, double *edgeN[], double *diff_edgeN[], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 
PetscErrorCode H1_FaceShapeFunctions_MBTET (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))
 
PetscErrorCode H1_VolumeShapeFunctions_MBTET (int p, double *N, double *diffN, double *volumeN, double *diff_volumeN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 
PetscErrorCode H1_EdgeShapeDiffMBTETinvJ (int *base_p, int *p, double *edge_diffN[], double *invJac, double *edge_diffNinvJac[], int GDIM)
 
PetscErrorCode H1_FaceShapeDiffMBTETinvJ (int *base_p, int *p, double *face_diffN[], double *invJac, double *face_diffNinvJac[], int GDIM)
 
PetscErrorCode H1_VolumeShapeDiffMBTETinvJ (int base_p, int p, double *volume_diffN, double *invJac, double *volume_diffNinvJac, int GDIM)
 
PetscErrorCode H1_EdgeGradientOfDeformation_hierarchical (int p, double *diffN, double *dofs, double *F)
 
PetscErrorCode H1_FaceGradientOfDeformation_hierarchical (int p, double *diffN, double *dofs, double *F)
 
PetscErrorCode H1_VolumeGradientOfDeformation_hierarchical (int p, double *diffN, double *dofs, double *F)
 
DEPRECATED PetscErrorCode H1_EdgeGradientOfDeformation_hierachical (int p, double *diffN, double *dofs, double *F)
 
DEPRECATED PetscErrorCode H1_FaceGradientOfDeformation_hierachical (int p, double *diffN, double *dofs, double *F)
 
DEPRECATED PetscErrorCode H1_VolumeGradientOfDeformation_hierachical (int p, double *diffN, double *dofs, double *F)
 
PetscErrorCode H1_QuadShapeFunctions_MBPRISM (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))
 
PetscErrorCode H1_VolumeShapeFunctions_MBPRISM (int p, double *N, double *diffN, double *volumeN, double *diff_volumeN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
 

Detailed Description

Functions to approximate hierarchical spaces.

: Name Shape Functions is used, in that context is more appropriate to use base functions. Need to be changed.

Definition in file h1_hdiv_hcurl_l2.h.

Macro Definition Documentation

◆ DEPRECATED

#define DEPRECATED

Definition at line 23 of file h1_hdiv_hcurl_l2.h.

◆ NBEDGE_AINSWORTH_HCURL

#define NBEDGE_AINSWORTH_HCURL (   P)    (((P) > 0) ? (P + 1) : 0)

Definition at line 79 of file h1_hdiv_hcurl_l2.h.

◆ NBEDGE_DEMKOWICZ_HCURL

#define NBEDGE_DEMKOWICZ_HCURL (   P)    (((P) > 0) ? (P) : 0)

Definition at line 90 of file h1_hdiv_hcurl_l2.h.

◆ NBEDGE_H1

#define NBEDGE_H1 (   P)    (((P) > 0) ? (P - 1) : 0)

Numer of base function on edge for H1 space.

Definition at line 55 of file h1_hdiv_hcurl_l2.h.

◆ NBEDGE_HDIV

#define NBEDGE_HDIV (   P)    (0)

Definition at line 97 of file h1_hdiv_hcurl_l2.h.

◆ NBEDGE_L2

#define NBEDGE_L2 (   P)    (P + 1)

Number of base functions on edge fro L2 space.

Definition at line 48 of file h1_hdiv_hcurl_l2.h.

◆ NBFACEQUAD_H1

#define NBFACEQUAD_H1 (   P)    (((P) > 2) ? ((P - 3) * (P - 2) / 2) : 0)

Number of base functions on quad for H1 space.

Definition at line 65 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_AINSWORTH_EDGE_HCURL

#define NBFACETRI_AINSWORTH_EDGE_HCURL (   P)    (((P) > 1) ? P - 1 : 0)

Definition at line 80 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_AINSWORTH_EDGE_HDIV

#define NBFACETRI_AINSWORTH_EDGE_HDIV (   P)    (((P) > 0) ? (P) : 0)

Definition at line 98 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_AINSWORTH_FACE_HCURL

#define NBFACETRI_AINSWORTH_FACE_HCURL (   P)    (((P) > 2) ? (P - 1) * (P - 2) : 0)

Definition at line 81 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_AINSWORTH_FACE_HDIV

#define NBFACETRI_AINSWORTH_FACE_HDIV (   P)    (((P) > 2) ? (P - 1) * (P - 2) / 2 : 0)

Definition at line 99 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_AINSWORTH_HCURL

#define NBFACETRI_AINSWORTH_HCURL (   P)    ((P > 1) ? ((P)-1) * (P + 1) : 0)

Definition at line 82 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_AINSWORTH_HDIV

#define NBFACETRI_AINSWORTH_HDIV (   P)    (((P) > 0) ? (P + 1) * (P + 2) / 2 : 0)

Definition at line 100 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_DEMKOWICZ_HCURL

#define NBFACETRI_DEMKOWICZ_HCURL (   P)    (((P) > 1) ? (P) * ((P)-1) : 0)

Definition at line 91 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_DEMKOWICZ_HDIV

#define NBFACETRI_DEMKOWICZ_HDIV (   P)    ((P > 0) ? (P) * (P + 1) / 2 : 0)
Examples:
forces_and_sources_testing_users_base.cpp.

Definition at line 107 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_H1

#define NBFACETRI_H1 (   P)    (((P) > 1) ? ((P - 2) * (P - 1) / 2) : 0)

Number of base function on triangle for H1 space.

Definition at line 60 of file h1_hdiv_hcurl_l2.h.

◆ NBFACETRI_L2

#define NBFACETRI_L2 (   P)    ((P + 1) * (P + 2) / 2)

Number of base functions on triangle for L2 space.

Definition at line 42 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMEPRISM_H1

#define NBVOLUMEPRISM_H1 (   P)    (((P) > 4) ? ((P - 5) * (P - 4) * (P - 3) / 6) : 0)

Number of base functions on prism for H1 space.

Definition at line 75 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_AINSWORTH_EDGE_HDIV

#define NBVOLUMETET_AINSWORTH_EDGE_HDIV (   P)    (((P) > 1) ? (P - 1) : 0)

Definition at line 101 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_AINSWORTH_FACE_HCURL

#define NBVOLUMETET_AINSWORTH_FACE_HCURL (   P)    (((P) > 2) ? (2 * (P - 1) * (P - 2)) : 0)

Definition at line 83 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_AINSWORTH_FACE_HDIV

#define NBVOLUMETET_AINSWORTH_FACE_HDIV (   P)    (((P) > 2) ? (P - 1) * (P - 2) : 0)

Definition at line 102 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_AINSWORTH_HCURL

#define NBVOLUMETET_AINSWORTH_HCURL (   P)    (((P) > 2) ? (P - 2) * (P - 1) * (P + 1) / 2 : 0)

Definition at line 87 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_AINSWORTH_HDIV

#define NBVOLUMETET_AINSWORTH_HDIV (   P)    (((P) > 1) ? (P - 1) * (P + 1) * (P + 2) / 2 : 0)

Definition at line 105 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_AINSWORTH_TET_HCURL

#define NBVOLUMETET_AINSWORTH_TET_HCURL (   P)    (((P) > 3) ? ((P - 3) * (P - 2) * (P - 1) / 2) : 0)

Definition at line 85 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_AINSWORTH_VOLUME_HDIV

#define NBVOLUMETET_AINSWORTH_VOLUME_HDIV (   P)    (((P) > 3) ? (P - 3) * (P - 2) * (P - 1) / 2 : 0)

Definition at line 103 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_DEMKOWICZ_HCURL

#define NBVOLUMETET_DEMKOWICZ_HCURL (   P)    (((P) > 2) ? ((P) * ((P)-1) * ((P)-2) / 2) : 0)

Definition at line 92 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_DEMKOWICZ_HDIV

#define NBVOLUMETET_DEMKOWICZ_HDIV (   P)    (((P) > 1) ? (P) * (P - 1) * (P + 1) / 2 : 0)
Examples:
forces_and_sources_testing_users_base.cpp.

Definition at line 108 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_H1

#define NBVOLUMETET_H1 (   P)    (((P) > 2) ? ((P - 3) * (P - 2) * (P - 1) / 6) : 0)

Number of base functions on tetrahedron fro H1 space.

Definition at line 70 of file h1_hdiv_hcurl_l2.h.

◆ NBVOLUMETET_L2

#define NBVOLUMETET_L2 (   P)    ((P + 1) * (P + 2) * (P + 3) / 6)

Number of base functions on tetrahedron for L2 space.

Definition at line 37 of file h1_hdiv_hcurl_l2.h.

Function Documentation

◆ H1_EdgeGradientOfDeformation_hierachical()

DEPRECATED PetscErrorCode H1_EdgeGradientOfDeformation_hierachical ( int  p,
double diffN,
double dofs,
double F 
)
Deprecated:
use H1_EdgeGradientOfDeformation_hierarchical

Definition at line 779 of file h1.c.

781  {
782  return H1_EdgeGradientOfDeformation_hierarchical(p, diffN, dofs, F);
783 }
PetscErrorCode H1_EdgeGradientOfDeformation_hierarchical(int p, double *diffN, double *dofs, double *F)
Definition: h1.c:553

◆ H1_EdgeGradientOfDeformation_hierarchical()

PetscErrorCode H1_EdgeGradientOfDeformation_hierarchical ( int  p,
double diffN,
double dofs,
double F 
)

Definition at line 553 of file h1.c.

555  {
557  int col, row = 0;
558  for (; row < 3; row++)
559  for (col = 0; col < 3; col++)
560  F[3 * row + col] =
561  cblas_ddot(NBEDGE_H1(p), &diffN[col], 3, &dofs[row], 3);
563 }
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
double cblas_ddot(const int N, const double *X, const int incX, const double *Y, const int incY)
Definition: cblas_ddot.c:12

◆ H1_EdgeShapeDiffMBTETinvJ()

PetscErrorCode H1_EdgeShapeDiffMBTETinvJ ( int base_p,
int p,
double edge_diffN[],
double invJac,
double edge_diffNinvJac[],
int  GDIM 
)

Definition at line 500 of file h1.c.

502  {
504  int ee = 0, ii, gg;
505  for (; ee < 6; ee++) {
506  for (ii = 0; ii < NBEDGE_H1(p[ee]); ii++) {
507  for (gg = 0; gg < GDIM; gg++) {
508  int shift1 = NBEDGE_H1(base_p[ee]) * gg;
509  int shift2 = NBEDGE_H1(p[ee]) * gg;
510  cblas_dgemv(CblasRowMajor, CblasTrans, 3, 3, 1., invJac, 3,
511  &(edge_diffN[ee])[3 * shift1 + 3 * ii], 1, 0.,
512  &(edge_diffNinvJac[ee])[3 * shift2 + 3 * ii], 1);
513  }
514  }
515  }
517 }
void cblas_dgemv(const enum CBLAS_ORDER order, const enum CBLAS_TRANSPOSE TransA, const int M, const int N, const double alpha, const double *A, const int lda, const double *X, const int incX, const double beta, double *Y, const int incY)
Definition: cblas_dgemv.c:11
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498

◆ H1_EdgeShapeFunctions_MBTET()

PetscErrorCode H1_EdgeShapeFunctions_MBTET ( int sense,
int p,
double N,
double diffN,
double edgeN[],
double diff_edgeN[],
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Definition at line 245 of file h1.c.

250  {
252 
253  int P[6];
254  int ee = 0;
255  for (; ee < 6; ee++)
256  P[ee] = NBEDGE_H1(p[ee]);
257  int edges_nodes[2 * 6] = {0, 1, 1, 2, 2, 0, 0, 3, 1, 3, 2, 3};
258  double diff_ksi01[3], diff_ksi12[3], diff_ksi20[3];
259  double diff_ksi03[3], diff_ksi13[3], diff_ksi23[3];
260  double *edges_diff_ksi[6] = {diff_ksi01, diff_ksi12, diff_ksi20,
261  diff_ksi03, diff_ksi13, diff_ksi23};
262  for (ee = 0; ee < 6; ee++) {
263  int dd = 0;
264  for (; dd < 3; dd++) {
265  edges_diff_ksi[ee][dd] = (diffN[edges_nodes[2 * ee + 1] * 3 + dd] -
266  diffN[edges_nodes[2 * ee + 0] * 3 + dd]) *
267  sense[ee];
268  }
269  }
270  int ii = 0;
271  for (; ii < GDIM; ii++) {
272  int node_shift = ii * 4;
273  double edges_ksi[6];
274  for (ee = 0; ee < 6; ee++) {
275  edges_ksi[ee] = (N[node_shift + edges_nodes[2 * ee + 1]] -
276  N[node_shift + edges_nodes[2 * ee + 0]]) *
277  sense[ee];
278  }
279  for (ee = 0; ee < 6; ee++) {
280  if (P[ee] == 0)
281  continue;
282  double L[p[ee] + 1], diffL[3 * (p[ee] + 1)];
283  ierr = base_polynomials(p[ee], edges_ksi[ee], edges_diff_ksi[ee], L,
284  diffL, 3);
285  CHKERRQ(ierr);
286  double v = N[node_shift + edges_nodes[2 * ee + 0]] *
287  N[node_shift + edges_nodes[2 * ee + 1]];
288  if (edgeN != NULL)
289  if (edgeN[ee] != NULL) {
290  int shift = ii * P[ee];
291  cblas_dcopy(P[ee], L, 1, &edgeN[ee][shift], 1);
292  cblas_dscal(P[ee],
293  N[node_shift + edges_nodes[2 * ee + 0]] *
294  N[node_shift + edges_nodes[2 * ee + 1]],
295  &edgeN[ee][shift], 1);
296  }
297  if (diff_edgeN != NULL)
298  if (diff_edgeN[ee] != NULL) {
299  int shift = ii * P[ee];
300  bzero(&diff_edgeN[ee][3 * shift], sizeof(double) * 3 * P[ee]);
301  int dd = 0;
302  for (; dd < 3; dd++) {
303  cblas_daxpy(P[ee], v, &diffL[dd * (p[ee] + 1)], 1,
304  &diff_edgeN[ee][3 * shift + dd], 3);
305  cblas_daxpy(P[ee],
306  diffN[3 * edges_nodes[2 * ee + 0] + dd] *
307  N[node_shift + edges_nodes[2 * ee + 1]] +
308  N[node_shift + edges_nodes[2 * ee + 0]] *
309  diffN[3 * edges_nodes[2 * ee + 1] + dd],
310  L, 1, &diff_edgeN[ee][3 * shift + dd], 3);
311  }
312  }
313  }
314  }
316 }
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
void cblas_daxpy(const int N, const double alpha, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_daxpy.c:11
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
void cblas_dcopy(const int N, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_dcopy.c:11
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
void cblas_dscal(const int N, const double alpha, double *X, const int incX)
Definition: cblas_dscal.c:11
CHKERRQ(ierr)
static PetscErrorCode ierr
Definition: h1.c:25
const int N
Definition: speed_test.cpp:3

◆ H1_EdgeShapeFunctions_MBTRI()

PetscErrorCode H1_EdgeShapeFunctions_MBTRI ( int sense,
int p,
double N,
double diffN,
double edgeN[3],
double diff_edgeN[3],
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

H1_EdgeShapeFunctions_MBTRI.

Parameters
senseof edges, it is array of integers dim 3 (3-edges of triangle)
pof edges

Definition at line 27 of file h1.c.

32  {
34 
35  double *edgeN01 = NULL, *edgeN12 = NULL, *edgeN20 = NULL;
36  if (edgeN != NULL) {
37  edgeN01 = edgeN[0];
38  edgeN12 = edgeN[1];
39  edgeN20 = edgeN[2];
40  }
41  double *diff_edgeN01 = NULL, *diff_edgeN12 = NULL, *diff_edgeN20 = NULL;
42  if (diff_edgeN != NULL) {
43  diff_edgeN01 = diff_edgeN[0];
44  diff_edgeN12 = diff_edgeN[1];
45  diff_edgeN20 = diff_edgeN[2];
46  }
47  int P[3];
48  int ee = 0;
49  for (; ee < 3; ee++) {
50  P[ee] = NBEDGE_H1(p[ee]);
51  }
52  int dd = 0;
53  double diff_ksi01[2], diff_ksi12[2], diff_ksi20[2];
54  if (diff_edgeN != NULL) {
55  for (; dd < 2; dd++) {
56  diff_ksi01[dd] = (diffN[1 * 2 + dd] - diffN[0 * 2 + dd]) * sense[0];
57  diff_ksi12[dd] = (diffN[2 * 2 + dd] - diffN[1 * 2 + dd]) * sense[1];
58  diff_ksi20[dd] = (diffN[0 * 2 + dd] - diffN[2 * 2 + dd]) * sense[2];
59  }
60  }
61  int ii = 0;
62  for (; ii < GDIM; ii++) {
63  int node_shift = ii * 3;
64  double ksi01 = (N[node_shift + 1] - N[node_shift + 0]) * sense[0];
65  double ksi12 = (N[node_shift + 2] - N[node_shift + 1]) * sense[1];
66  double ksi20 = (N[node_shift + 0] - N[node_shift + 2]) * sense[2];
67  double L01[p[0] + 1], L12[p[1] + 1], L20[p[2] + 1];
68  double diffL01[2 * (p[0] + 1)], diffL12[2 * (p[1] + 1)],
69  diffL20[2 * (p[2] + 1)];
70  ierr = base_polynomials(p[0], ksi01, diff_ksi01, L01, diffL01, 2);
71  CHKERRQ(ierr);
72  ierr = base_polynomials(p[1], ksi12, diff_ksi12, L12, diffL12, 2);
73  CHKERRQ(ierr);
74  ierr = base_polynomials(p[2], ksi20, diff_ksi20, L20, diffL20, 2);
75  CHKERRQ(ierr);
76  int shift;
77  if (edgeN != NULL) {
78  // edge01
79  shift = ii * (P[0]);
80  cblas_dcopy(P[0], L01, 1, &edgeN01[shift], 1);
81  cblas_dscal(P[0], N[node_shift + 0] * N[node_shift + 1], &edgeN01[shift],
82  1);
83  // edge12
84  shift = ii * (P[1]);
85  cblas_dcopy(P[1], L12, 1, &edgeN12[shift], 1);
86  cblas_dscal(P[1], N[node_shift + 1] * N[node_shift + 2], &edgeN12[shift],
87  1);
88  // edge20
89  shift = ii * (P[2]);
90  cblas_dcopy(P[2], L20, 1, &edgeN20[shift], 1);
91  cblas_dscal(P[2], N[node_shift + 0] * N[node_shift + 2], &edgeN20[shift],
92  1);
93  }
94  if (diff_edgeN != NULL) {
95  if (P[0] > 0) {
96  // edge01
97  shift = ii * (P[0]);
98  bzero(&diff_edgeN01[2 * shift], sizeof(double) * 2 * (P[0]));
99  // diffX
100  cblas_daxpy(P[0], N[node_shift + 0] * N[node_shift + 1],
101  &diffL01[0 * (p[0] + 1)], 1, &diff_edgeN01[2 * shift + 0],
102  2);
103  cblas_daxpy(P[0],
104  diffN[2 * 0 + 0] * N[node_shift + 1] +
105  N[node_shift + 0] * diffN[2 * 1 + 0],
106  L01, 1, &diff_edgeN01[2 * shift + 0], 2);
107  // diff Y
108  cblas_daxpy(P[0], N[node_shift + 0] * N[node_shift + 1],
109  &diffL01[1 * (p[0] + 1)], 1, &diff_edgeN01[2 * shift + 1],
110  2);
111  cblas_daxpy(P[0],
112  diffN[2 * 0 + 1] * N[node_shift + 1] +
113  N[node_shift + 0] * diffN[2 * 1 + 1],
114  L01, 1, &diff_edgeN01[2 * shift + 1], 2);
115  }
116  if (P[1] > 0) {
117  // edge12
118  shift = ii * (P[1]);
119  bzero(&diff_edgeN12[2 * shift], sizeof(double) * 2 * (P[1]));
120  // diffX
121  cblas_daxpy(P[1], N[node_shift + 1] * N[node_shift + 2],
122  &diffL12[0 * (p[1] + 1)], 1, &diff_edgeN12[2 * shift + 0],
123  2);
124  cblas_daxpy(P[1],
125  diffN[2 * 1 + 0] * N[node_shift + 2] +
126  N[node_shift + 1] * diffN[2 * 2 + 0],
127  L12, 1, &diff_edgeN12[2 * shift + 0], 2);
128  // diffY
129  cblas_daxpy(P[1], N[node_shift + 1] * N[node_shift + 2],
130  &diffL12[1 * (p[1] + 1)], 1, &diff_edgeN12[2 * shift + 1],
131  2);
132  cblas_daxpy(P[1],
133  diffN[2 * 1 + 1] * N[node_shift + 2] +
134  N[node_shift + 1] * diffN[2 * 2 + 1],
135  L12, 1, &diff_edgeN12[2 * shift + 1], 2);
136  }
137  if (P[2] > 0) {
138  // edge20
139  shift = ii * (P[2]);
140  bzero(&diff_edgeN20[2 * shift], sizeof(double) * 2 * (P[2]));
141  // diffX
142  cblas_daxpy(P[2], N[node_shift + 0] * N[node_shift + 2],
143  &diffL20[0 * (p[2] + 1)], 1, &diff_edgeN20[2 * shift + 0],
144  2);
145  cblas_daxpy(P[2],
146  diffN[2 * 0 + 0] * N[node_shift + 2] +
147  N[node_shift + 0] * diffN[2 * 2 + 0],
148  L20, 1, &diff_edgeN20[2 * shift + 0], 2);
149  // diffY
150  cblas_daxpy(P[2], N[node_shift + 0] * N[node_shift + 2],
151  &diffL20[1 * (p[2] + 1)], 1, &diff_edgeN20[2 * shift + 1],
152  2);
153  cblas_daxpy(P[2],
154  diffN[2 * 0 + 1] * N[node_shift + 2] +
155  N[node_shift + 0] * diffN[2 * 2 + 1],
156  L20, 1, &diff_edgeN20[2 * shift + 1], 2);
157  }
158  }
159  }
161 }
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
void cblas_daxpy(const int N, const double alpha, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_daxpy.c:11
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
void cblas_dcopy(const int N, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_dcopy.c:11
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
void cblas_dscal(const int N, const double alpha, double *X, const int incX)
Definition: cblas_dscal.c:11
CHKERRQ(ierr)
static PetscErrorCode ierr
Definition: h1.c:25
const int N
Definition: speed_test.cpp:3

◆ H1_FaceGradientOfDeformation_hierachical()

DEPRECATED PetscErrorCode H1_FaceGradientOfDeformation_hierachical ( int  p,
double diffN,
double dofs,
double F 
)
Deprecated:
use H1_FaceGradientOfDeformation_hierarchical

Definition at line 784 of file h1.c.

786  {
787  return H1_FaceGradientOfDeformation_hierarchical(p, diffN, dofs, F);
788 }
PetscErrorCode H1_FaceGradientOfDeformation_hierarchical(int p, double *diffN, double *dofs, double *F)
Definition: h1.c:564

◆ H1_FaceGradientOfDeformation_hierarchical()

PetscErrorCode H1_FaceGradientOfDeformation_hierarchical ( int  p,
double diffN,
double dofs,
double F 
)

Definition at line 564 of file h1.c.

566  {
568  int col, row = 0;
569  for (; row < 3; row++)
570  for (col = 0; col < 3; col++)
571  F[3 * row + col] =
572  cblas_ddot(NBFACETRI_H1(p), &diffN[col], 3, &dofs[row], 3);
574 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
double cblas_ddot(const int N, const double *X, const int incX, const double *Y, const int incY)
Definition: cblas_ddot.c:12
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.

◆ H1_FaceShapeDiffMBTETinvJ()

PetscErrorCode H1_FaceShapeDiffMBTETinvJ ( int base_p,
int p,
double face_diffN[],
double invJac,
double face_diffNinvJac[],
int  GDIM 
)

Definition at line 518 of file h1.c.

520  {
522  int ff = 0, ii, gg;
523  for (; ff < 4; ff++) {
524  for (ii = 0; ii < NBFACETRI_H1(p[ff]); ii++) {
525  for (gg = 0; gg < GDIM; gg++) {
526  int shift1 = NBFACETRI_H1(base_p[ff]) * gg;
527  int shift2 = NBFACETRI_H1(p[ff]) * gg;
528  cblas_dgemv(CblasRowMajor, CblasTrans, 3, 3, 1., invJac, 3,
529  &(face_diffN[ff])[3 * shift1 + 3 * ii], 1, 0.,
530  &(face_diffNinvJac[ff])[3 * shift2 + 3 * ii], 1);
531  }
532  }
533  }
535 }
void cblas_dgemv(const enum CBLAS_ORDER order, const enum CBLAS_TRANSPOSE TransA, const int M, const int N, const double alpha, const double *A, const int lda, const double *X, const int incX, const double beta, double *Y, const int incY)
Definition: cblas_dgemv.c:11
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.

◆ H1_FaceShapeFunctions_MBTET()

PetscErrorCode H1_FaceShapeFunctions_MBTET ( int faces_nodes,
int p,
double N,
double diffN,
double faceN[],
double diff_faceN[],
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Definition at line 317 of file h1.c.

322  {
324 
325  int P[4];
326  int ff = 0;
327  for (; ff < 4; ff++)
328  P[ff] = NBFACETRI_H1(p[ff]);
329  double diff_ksiL0F0[3], diff_ksiL1F0[3];
330  double diff_ksiL0F1[3], diff_ksiL1F1[3];
331  double diff_ksiL0F2[3], diff_ksiL1F2[3];
332  double diff_ksiL0F3[3], diff_ksiL1F3[3];
333  double *diff_ksi_faces[] = {diff_ksiL0F0, diff_ksiL1F0, diff_ksiL0F1,
334  diff_ksiL1F1, diff_ksiL0F2, diff_ksiL1F2,
335  diff_ksiL0F3, diff_ksiL1F3};
336  int dd = 0;
337  for (; dd < 3; dd++) {
338  for (ff = 0; ff < 4; ff++) {
339  diff_ksi_faces[ff * 2 + 0][dd] =
340  (diffN[faces_nodes[3 * ff + 1] * 3 + dd] -
341  diffN[faces_nodes[3 * ff + 0] * 3 + dd]);
342  diff_ksi_faces[ff * 2 + 1][dd] =
343  (diffN[faces_nodes[3 * ff + 2] * 3 + dd] -
344  diffN[faces_nodes[3 * ff + 0] * 3 + dd]);
345  }
346  }
347  int ii = 0;
348  for (; ii < GDIM; ii++) {
349  int node_shift = ii * 4;
350  double ksi_faces[8];
351  for (ff = 0; ff < 4; ff++) {
352  ksi_faces[2 * ff + 0] = N[node_shift + faces_nodes[3 * ff + 1]] -
353  N[node_shift + faces_nodes[3 * ff + 0]];
354  ksi_faces[2 * ff + 1] = N[node_shift + faces_nodes[3 * ff + 2]] -
355  N[node_shift + faces_nodes[3 * ff + 0]];
356  }
357  int shift;
358  for (ff = 0; ff < 4; ff++) {
359  if (P[ff] == 0)
360  continue;
361  double L0[p[ff] + 1], L1[p[ff] + 1];
362  double diffL0[3 * (p[ff] + 1)], diffL1[3 * (p[ff] + 1)];
363  ierr = base_polynomials(p[ff], ksi_faces[ff * 2 + 0],
364  diff_ksi_faces[ff * 2 + 0], L0, diffL0, 3);
365  CHKERRQ(ierr);
366  ierr = base_polynomials(p[ff], ksi_faces[ff * 2 + 1],
367  diff_ksi_faces[ff * 2 + 1], L1, diffL1, 3);
368  CHKERRQ(ierr);
369  double v = N[node_shift + faces_nodes[3 * ff + 0]] *
370  N[node_shift + faces_nodes[3 * ff + 1]] *
371  N[node_shift + faces_nodes[3 * ff + 2]];
372  double v2[3] = {0, 0, 0};
373  dd = 0;
374  double n1n2 = N[node_shift + faces_nodes[3 * ff + 1]] *
375  N[node_shift + faces_nodes[3 * ff + 2]];
376  double n0n2 = N[node_shift + faces_nodes[3 * ff + 0]] *
377  N[node_shift + faces_nodes[3 * ff + 2]];
378  double n0n1 = N[node_shift + faces_nodes[3 * ff + 0]] *
379  N[node_shift + faces_nodes[3 * ff + 1]];
380  for (; dd < 3; dd++) {
381  v2[dd] = diffN[3 * faces_nodes[3 * ff + 0] + dd] * n1n2 +
382  diffN[3 * faces_nodes[3 * ff + 1] + dd] * n0n2 +
383  diffN[3 * faces_nodes[3 * ff + 2] + dd] * n0n1;
384  }
385  shift = ii * P[ff];
386  int jj = 0;
387  int oo = 0;
388  for (; oo <= (p[ff] - 3); oo++) {
389  int pp0 = 0;
390  for (; pp0 <= oo; pp0++) {
391  int pp1 = oo - pp0;
392  if (pp1 >= 0) {
393  if (faceN != NULL)
394  if (faceN[ff] != NULL) {
395  faceN[ff][shift + jj] = v * L0[pp0] * L1[pp1];
396  }
397  if (diff_faceN != NULL)
398  if (diff_faceN[ff] != NULL) {
399  dd = 0;
400  double L0L1 = L0[pp0] * L1[pp1];
401  for (; dd < 3; dd++) {
402  diff_faceN[ff][3 * shift + 3 * jj + dd] =
403  (L0[pp0] * diffL1[dd * (p[ff] + 1) + pp1] +
404  diffL0[dd * (p[ff] + 1) + pp0] * L1[pp1]) *
405  v;
406  diff_faceN[ff][3 * shift + 3 * jj + dd] += L0L1 * v2[dd];
407  }
408  }
409  jj++;
410  }
411  }
412  }
413  if (jj != P[ff])
414  SETERRQ2(PETSC_COMM_SELF, 1, "wrong order %d != %d", jj, P[ff]);
415  }
416  }
418 }
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
CHKERRQ(ierr)
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.
static PetscErrorCode ierr
Definition: h1.c:25
const int N
Definition: speed_test.cpp:3

◆ H1_FaceShapeFunctions_MBTRI()

PetscErrorCode H1_FaceShapeFunctions_MBTRI ( const int face_nodes,
int  p,
double N,
double diffN,
double faceN,
double diff_faceN,
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Definition at line 162 of file h1.c.

167  {
169 
170  int P = NBFACETRI_H1(p);
171  if (P == 0)
173  double diff_ksiL0[2], diff_ksiL1[2];
174  double *diff_ksi_faces[] = {diff_ksiL0, diff_ksiL1};
175  int dd = 0;
176  if (diff_faceN != NULL) {
177  for (; dd < 2; dd++) {
178  diff_ksi_faces[0][dd] =
179  diffN[face_nodes[1] * 2 + dd] - diffN[face_nodes[0] * 2 + dd];
180  diff_ksi_faces[1][dd] =
181  diffN[face_nodes[2] * 2 + dd] - diffN[face_nodes[0] * 2 + dd];
182  }
183  }
184  int ii = 0;
185  for (; ii < GDIM; ii++) {
186  int node_shift = ii * 3;
187  double ksi_faces[2];
188  ksi_faces[0] =
189  N[node_shift + face_nodes[1]] - N[node_shift + face_nodes[0]];
190  ksi_faces[1] =
191  N[node_shift + face_nodes[2]] - N[node_shift + face_nodes[0]];
192  double L0[p + 1], L1[p + 1];
193  double diffL0[2 * (p + 1)], diffL1[2 * (p + 1)];
194  ierr = base_polynomials(p, ksi_faces[0], diff_ksi_faces[0], L0, diffL0, 2);
195  CHKERRQ(ierr);
196  ierr = base_polynomials(p, ksi_faces[1], diff_ksi_faces[1], L1, diffL1, 2);
197  CHKERRQ(ierr);
198  double v = N[node_shift + face_nodes[0]] * N[node_shift + face_nodes[1]] *
199  N[node_shift + face_nodes[2]];
200  double v2[2] = {0, 0};
201  if (diff_faceN != NULL) {
202  double n1n2 =
203  N[node_shift + face_nodes[1]] * N[node_shift + face_nodes[2]];
204  double n0n2 =
205  N[node_shift + face_nodes[0]] * N[node_shift + face_nodes[2]];
206  double n0n1 =
207  N[node_shift + face_nodes[0]] * N[node_shift + face_nodes[1]];
208  dd = 0;
209  for (; dd < 2; dd++) {
210  v2[dd] = diffN[face_nodes[0] * 2 + dd] * n1n2 +
211  diffN[face_nodes[1] * 2 + dd] * n0n2 +
212  diffN[face_nodes[2] * 2 + dd] * n0n1;
213  }
214  }
215  int shift = ii * P;
216  int jj = 0;
217  int oo = 0;
218  for (; oo <= (p - 3); oo++) {
219  int pp0 = 0;
220  for (; pp0 <= oo; pp0++) {
221  int pp1 = oo - pp0;
222  if (pp1 >= 0) {
223  if (faceN != NULL) {
224  faceN[shift + jj] = v * L0[pp0] * L1[pp1];
225  }
226  if (diff_faceN != NULL) {
227  dd = 0;
228  for (; dd < 2; dd++) {
229  double *val = &diff_faceN[2 * shift + 2 * jj + dd];
230  *val = (L0[pp0] * diffL1[dd * (p + 1) + pp1] +
231  diffL0[dd * (p + 1) + pp0] * L1[pp1]) *
232  v;
233  *val += L0[pp0] * L1[pp1] * v2[dd];
234  }
235  }
236  jj++;
237  }
238  }
239  }
240  if (jj != P)
241  SETERRQ2(PETSC_COMM_SELF, 1, "wrong order %d != %d", jj, P);
242  }
244 }
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
CHKERRQ(ierr)
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.
static PetscErrorCode ierr
Definition: h1.c:25
const int N
Definition: speed_test.cpp:3

◆ H1_QuadShapeFunctions_MBPRISM()

PetscErrorCode H1_QuadShapeFunctions_MBPRISM ( int faces_nodes,
int p,
double N,
double diffN,
double faceN[],
double diff_faceN[],
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Definition at line 586 of file h1.c.

591  {
593 
594  int P[3];
595  int ff = 0;
596  for (; ff < 3; ff++) {
597  P[ff] = NBFACEQUAD_H1(p[ff]);
598  }
599  double ksi_faces[6];
600  double diff_ksiL0F0[3], diff_ksiL3F0[3];
601  double diff_ksiL0F1[3], diff_ksiL3F1[3];
602  double diff_ksiL0F2[3], diff_ksiL3F2[3];
603  double *diff_ksi_faces[] = {diff_ksiL0F0, diff_ksiL3F0, diff_ksiL0F1,
604  diff_ksiL3F1, diff_ksiL0F2, diff_ksiL3F2};
605  int ii = 0;
606  for (; ii < GDIM; ii++) {
607  int node_shift = ii * 6;
608  int node_diff_shift = 3 * node_shift;
609  int ff = 0;
610  for (; ff < 3; ff++) {
611  if (P[ff] == 0)
612  continue;
613  int n0 = faces_nodes[4 * ff + 0];
614  int n1 = faces_nodes[4 * ff + 1];
615  int n2 = faces_nodes[4 * ff + 2];
616  int n3 = faces_nodes[4 * ff + 3];
617  int e0 = 2 * ff + 0;
618  int e1 = 2 * ff + 1;
619  {
620  ksi_faces[e0] = 2. * (N[node_shift + n1] + N[node_shift + n2] -
621  N[node_shift + n0] - N[node_shift + n3]) -
622  1.;
623  ksi_faces[e1] = 2. * (N[node_shift + n2] + N[node_shift + n3] -
624  N[node_shift + n0] - N[node_shift + n1]) -
625  1.;
626  }
627  // if(ksi_faces[e0]<-1 || ksi_faces[e0]>1) {
628  // SETERRQ(PETSC_COMM_SELF,1,"Data inconsistency");
629  // }
630  // if(ksi_faces[e1]<-1 || ksi_faces[e1]>1) {
631  // SETERRQ(PETSC_COMM_SELF,1,"Data inconsistency");
632  // }
633  int dd = 0;
634  for (; dd < 3; dd++) {
635  diff_ksi_faces[e0][dd] = 2. * (diffN[node_diff_shift + 3 * n1 + dd] +
636  diffN[node_diff_shift + 3 * n2 + dd] -
637  diffN[node_diff_shift + 3 * n0 + dd] -
638  diffN[node_diff_shift + 3 * n3 + dd]);
639  diff_ksi_faces[e1][dd] = 2. * (diffN[node_diff_shift + 3 * n2 + dd] +
640  diffN[node_diff_shift + 3 * n3 + dd] -
641  diffN[node_diff_shift + 3 * n0 + dd] -
642  diffN[node_diff_shift + 3 * n1 + dd]);
643  }
644  double L0[p[ff] + 1], L1[p[ff] + 1];
645  double diffL0[3 * (p[ff] + 1)], diffL1[3 * (p[ff] + 1)];
646  ierr = base_polynomials(p[ff], ksi_faces[e0], diff_ksi_faces[e0], L0,
647  diffL0, 3);
648  CHKERRQ(ierr);
649  ierr = base_polynomials(p[ff], ksi_faces[e1], diff_ksi_faces[e1], L1,
650  diffL1, 3);
651  CHKERRQ(ierr);
652  double v = v = N[node_shift + n0] * N[node_shift + n2];
653  double v2[3] = {0, 0, 0};
654  dd = 0;
655  for (; dd < 3; dd++) {
656  v2[dd] = diffN[node_diff_shift + 3 * n0 + dd] * N[node_shift + n2] +
657  N[node_shift + n0] * diffN[node_diff_shift + 3 * n2 + dd];
658  }
659  int shift;
660  shift = ii * P[ff];
661  int jj = 0;
662  int oo = 0;
663  for (; oo <= (p[ff] - 4); oo++) {
664  int pp0 = 0;
665  for (; pp0 <= oo; pp0++) {
666  int pp1 = oo - pp0;
667  if (pp1 >= 0) {
668  if (faceN != NULL) {
669  if (faceN[ff] != NULL) {
670  faceN[ff][shift + jj] = v * L0[pp0] * L1[pp1];
671  }
672  }
673  if (diff_faceN != NULL) {
674  if (diff_faceN[ff] != NULL) {
675  dd = 0;
676  for (; dd < 3; dd++) {
677  diff_faceN[ff][3 * shift + 3 * jj + dd] =
678  (L0[pp0] * diffL1[dd * (p[ff] + 1) + pp1] +
679  diffL0[dd * (p[ff] + 1) + pp0] * L1[pp1]) *
680  v;
681  diff_faceN[ff][3 * shift + 3 * jj + dd] +=
682  L0[pp0] * L1[pp1] * v2[dd];
683  }
684  }
685  }
686  jj++;
687  }
688  }
689  }
690  if (jj != P[ff])
691  SETERRQ2(PETSC_COMM_SELF, 1, "wrong order %d != %d", jj, P[ff]);
692  }
693  }
695 }
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
CHKERRQ(ierr)
#define NBFACEQUAD_H1(P)
Number of base functions on quad for H1 space.
static PetscErrorCode ierr
Definition: h1.c:25
const int N
Definition: speed_test.cpp:3

◆ H1_VolumeGradientOfDeformation_hierachical()

DEPRECATED PetscErrorCode H1_VolumeGradientOfDeformation_hierachical ( int  p,
double diffN,
double dofs,
double F 
)
Deprecated:
use H1_VolumeGradientOfDeformation_hierarchical

Definition at line 789 of file h1.c.

791  {
792  return H1_VolumeGradientOfDeformation_hierarchical(p, diffN, dofs, F);
793 }
PetscErrorCode H1_VolumeGradientOfDeformation_hierarchical(int p, double *diffN, double *dofs, double *F)
Definition: h1.c:575

◆ H1_VolumeGradientOfDeformation_hierarchical()

PetscErrorCode H1_VolumeGradientOfDeformation_hierarchical ( int  p,
double diffN,
double dofs,
double F 
)

Definition at line 575 of file h1.c.

577  {
579  int col, row = 0;
580  for (; row < 3; row++)
581  for (col = 0; col < 3; col++)
582  F[3 * row + col] =
583  cblas_ddot(NBVOLUMETET_H1(p), &diffN[col], 3, &dofs[row], 3);
585 }
#define NBVOLUMETET_H1(P)
Number of base functions on tetrahedron fro H1 space.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
double cblas_ddot(const int N, const double *X, const int incX, const double *Y, const int incY)
Definition: cblas_ddot.c:12

◆ H1_VolumeShapeDiffMBTETinvJ()

PetscErrorCode H1_VolumeShapeDiffMBTETinvJ ( int  base_p,
int  p,
double volume_diffN,
double invJac,
double volume_diffNinvJac,
int  GDIM 
)

Definition at line 536 of file h1.c.

539  {
541  int ii, gg;
542  for (ii = 0; ii < NBVOLUMETET_H1(p); ii++) {
543  for (gg = 0; gg < GDIM; gg++) {
544  int shift1 = NBVOLUMETET_H1(base_p) * gg;
545  int shift2 = NBVOLUMETET_H1(p) * gg;
546  cblas_dgemv(CblasRowMajor, CblasTrans, 3, 3, 1., invJac, 3,
547  &(volume_diffN)[3 * shift1 + 3 * ii], 1, 0.,
548  &(volume_diffNinvJac)[3 * shift2 + 3 * ii], 1);
549  }
550  }
552 }
void cblas_dgemv(const enum CBLAS_ORDER order, const enum CBLAS_TRANSPOSE TransA, const int M, const int N, const double alpha, const double *A, const int lda, const double *X, const int incX, const double beta, double *Y, const int incY)
Definition: cblas_dgemv.c:11
#define NBVOLUMETET_H1(P)
Number of base functions on tetrahedron fro H1 space.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498

◆ H1_VolumeShapeFunctions_MBPRISM()

PetscErrorCode H1_VolumeShapeFunctions_MBPRISM ( int  p,
double N,
double diffN,
double volumeN,
double diff_volumeN,
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Definition at line 696 of file h1.c.

701  {
703 
704  int P = NBVOLUMEPRISM_H1(p);
705  if (P == 0)
707  double diff_ksiL0[3], diff_ksiL1[3], diff_ksiL2[3];
708  int ii = 0;
709  for (; ii < GDIM; ii++) {
710  int node_shift = ii * 6;
711  int node_diff_shift = ii * 18;
712  double ksiL0 = N[node_shift + 1] - N[node_shift + 0];
713  double ksiL1 = N[node_shift + 2] - N[node_shift + 0];
714  double ksiL2 =
715  N[node_shift + 3] - N[node_shift + 0]; // 2*gauss_pts[2*GDIM+ii]-1.;
716  int dd = 0;
717  for (; dd < 3; dd++) {
718  diff_ksiL0[dd] = (diffN[node_diff_shift + 1 * 3 + dd] -
719  diffN[node_diff_shift + 0 * 3 + dd]);
720  diff_ksiL1[dd] = (diffN[node_diff_shift + 2 * 3 + dd] -
721  diffN[node_diff_shift + 0 * 3 + dd]);
722  diff_ksiL2[dd] = 2;
723  }
724  double L0[p + 1], L1[p + 1], L2[p + 1];
725  double diffL0[3 * (p + 1)], diffL1[3 * (p + 1)], diffL2[3 * (p + 1)];
726  ierr = base_polynomials(p, ksiL0, diff_ksiL0, L0, diffL0, 3);
727  CHKERRQ(ierr);
728  ierr = base_polynomials(p, ksiL1, diff_ksiL1, L1, diffL1, 3);
729  CHKERRQ(ierr);
730  ierr = base_polynomials(p, ksiL2, diff_ksiL2, L2, diffL2, 3);
731  CHKERRQ(ierr);
732  double t0 = (N[node_shift + 2] + N[node_shift + 5]);
733  double v = N[node_shift + 0] * N[node_shift + 4] * t0;
734  double v2[3] = {0, 0, 0};
735  dd = 0;
736  for (; dd < 3; dd++) {
737  v2[dd] = diffN[node_diff_shift + 3 * 0 + dd] * N[node_shift + 4] * t0 +
738  N[node_shift + 0] * diffN[node_diff_shift + 3 * 4 + dd] * t0 +
739  N[node_shift + 0] * N[node_shift + 4] *
740  (diffN[node_diff_shift + 3 * 2 + dd] +
741  diffN[node_diff_shift + 3 * 5 + dd]);
742  }
743  int shift = ii * P;
744  int jj = 0;
745  int oo = 0;
746  for (; oo <= (p - 6); oo++) {
747  int pp0 = 0;
748  for (; pp0 <= oo; pp0++) {
749  int pp1 = 0;
750  for (; (pp0 + pp1) <= oo; pp1++) {
751  int pp2 = oo - pp0 - pp1;
752  if (pp2 >= 0) {
753  if (volumeN != NULL) {
754  volumeN[shift + jj] = L0[pp0] * L1[pp1] * L2[pp2] * v;
755  }
756  if (diff_volumeN != NULL) {
757  dd = 0;
758  for (; dd < 3; dd++) {
759  diff_volumeN[3 * shift + 3 * jj + dd] =
760  (diffL0[dd * (p + 1) + pp0] * L1[pp1] * L2[pp2] +
761  L0[pp0] * diffL1[dd * (p + 1) + pp1] * L2[pp2] +
762  L0[pp0] * L1[pp1] * diffL2[dd * (p + 1) + pp2]) *
763  v;
764  diff_volumeN[3 * shift + 3 * jj + dd] +=
765  L0[pp0] * L1[pp1] * L2[pp2] * v2[dd];
766  }
767  }
768  jj++;
769  }
770  }
771  }
772  }
773  if (jj != P)
774  SETERRQ1(PETSC_COMM_SELF, 1, "wrong order %d", jj);
775  }
777 }
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
CHKERRQ(ierr)
#define NBVOLUMEPRISM_H1(P)
Number of base functions on prism for H1 space.
static PetscErrorCode ierr
Definition: h1.c:25
const int N
Definition: speed_test.cpp:3
field with C-1 continuity
Definition: definitions.h:171

◆ H1_VolumeShapeFunctions_MBTET()

PetscErrorCode H1_VolumeShapeFunctions_MBTET ( int  p,
double N,
double diffN,
double volumeN,
double diff_volumeN,
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Definition at line 419 of file h1.c.

424  {
426 
427  int P = NBVOLUMETET_H1(p);
428  if (P == 0)
430  double diff_ksiL0[3], diff_ksiL1[3], diff_ksiL2[3];
431  int dd = 0;
432  for (; dd < 3; dd++) {
433  diff_ksiL0[dd] = (diffN[1 * 3 + dd] - diffN[0 * 3 + dd]);
434  diff_ksiL1[dd] = (diffN[2 * 3 + dd] - diffN[0 * 3 + dd]);
435  diff_ksiL2[dd] = (diffN[3 * 3 + dd] - diffN[0 * 3 + dd]);
436  }
437  int ii = 0;
438  for (; ii < GDIM; ii++) {
439  int node_shift = ii * 4;
440  double ksiL0 = N[node_shift + 1] - N[node_shift + 0];
441  double ksiL1 = N[node_shift + 2] - N[node_shift + 0];
442  double ksiL2 = N[node_shift + 3] - N[node_shift + 0];
443  double L0[p + 1], L1[p + 1], L2[p + 1];
444  double diffL0[3 * (p + 1)], diffL1[3 * (p + 1)], diffL2[3 * (p + 1)];
445  ierr = base_polynomials(p, ksiL0, diff_ksiL0, L0, diffL0, 3);
446  CHKERRQ(ierr);
447  ierr = base_polynomials(p, ksiL1, diff_ksiL1, L1, diffL1, 3);
448  CHKERRQ(ierr);
449  ierr = base_polynomials(p, ksiL2, diff_ksiL2, L2, diffL2, 3);
450  CHKERRQ(ierr);
451  double v = N[node_shift + 0] * N[node_shift + 1] * N[node_shift + 2] *
452  N[node_shift + 3];
453  double v2[3] = {0, 0, 0};
454  dd = 0;
455  for (; dd < 3; dd++) {
456  v2[dd] = diffN[3 * 0 + dd] * N[node_shift + 1] * N[node_shift + 2] *
457  N[node_shift + 3] +
458  N[node_shift + 0] * diffN[3 * 1 + dd] * N[node_shift + 2] *
459  N[node_shift + 3] +
460  N[node_shift + 0] * N[node_shift + 1] * diffN[3 * 2 + dd] *
461  N[node_shift + 3] +
462  N[node_shift + 0] * N[node_shift + 1] * N[node_shift + 2] *
463  diffN[3 * 3 + dd];
464  }
465  int shift = ii * P;
466  int jj = 0;
467  int oo = 0;
468  for (; oo <= (p - 4); oo++) {
469  int pp0 = 0;
470  for (; pp0 <= oo; pp0++) {
471  int pp1 = 0;
472  for (; (pp0 + pp1) <= oo; pp1++) {
473  int pp2 = oo - pp0 - pp1;
474  if (pp2 >= 0) {
475  if (volumeN != NULL) {
476  volumeN[shift + jj] = L0[pp0] * L1[pp1] * L2[pp2] * v;
477  }
478  if (diff_volumeN != NULL) {
479  dd = 0;
480  for (; dd < 3; dd++) {
481  diff_volumeN[3 * shift + 3 * jj + dd] =
482  (diffL0[dd * (p + 1) + pp0] * L1[pp1] * L2[pp2] +
483  L0[pp0] * diffL1[dd * (p + 1) + pp1] * L2[pp2] +
484  L0[pp0] * L1[pp1] * diffL2[dd * (p + 1) + pp2]) *
485  v;
486  diff_volumeN[3 * shift + 3 * jj + dd] +=
487  L0[pp0] * L1[pp1] * L2[pp2] * v2[dd];
488  }
489  }
490  jj++;
491  }
492  }
493  }
494  }
495  if (jj != P)
496  SETERRQ1(PETSC_COMM_SELF, 1, "wrong order %d", jj);
497  }
499 }
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define NBVOLUMETET_H1(P)
Number of base functions on tetrahedron fro H1 space.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
CHKERRQ(ierr)
static PetscErrorCode ierr
Definition: h1.c:25
const int N
Definition: speed_test.cpp:3
field with C-1 continuity
Definition: definitions.h:171

◆ L2_Ainsworth_ShapeFunctions_MBTET()

PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTET ( int  p,
double N,
double diffN,
double L2N,
double diff_L2N,
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Get base functions on tetrahedron for L2 space.

Parameters
ppolynomial order
Nbarycentric coordinates (shape functions) at integration points
diffNdirevatives of barycentric coordinates, i.e. direvatives of shape functions
L2Nvalues of L2 base at integration points
diff_L2Ndirvatives of base functions at integration points
GDIMnumber of integration points
base_polynomialspolynomial base used to construct L2 base on element
Returns
PetscErrorCode

Definition at line 84 of file l2.c.

88  {
90 
91  int P = NBVOLUMETET_L2(p);
92  if (P == 0)
94  double diff_ksiL0[3], diff_ksiL1[3], diff_ksiL2[3];
95  int dd = 0;
96  if (diffN != NULL) {
97  for (; dd < 3; dd++) {
98  diff_ksiL0[dd] = (diffN[1 * 3 + dd] - diffN[0 * 3 + dd]);
99  diff_ksiL1[dd] = (diffN[2 * 3 + dd] - diffN[0 * 3 + dd]);
100  diff_ksiL2[dd] = (diffN[3 * 3 + dd] - diffN[0 * 3 + dd]);
101  }
102  }
103  int ii = 0;
104  for (; ii != GDIM; ++ii) {
105  int node_shift = ii * 4;
106  double ksiL0 = N[node_shift + 1] - N[node_shift + 0];
107  double ksiL1 = N[node_shift + 2] - N[node_shift + 0];
108  double ksiL2 = N[node_shift + 3] - N[node_shift + 0];
109  double L0[p + 1], L1[p + 1], L2[p + 1];
110  double diffL0[3 * (p + 1)], diffL1[3 * (p + 1)], diffL2[3 * (p + 1)];
111  if (diffN != NULL) {
112  ierr = base_polynomials(p, ksiL0, diff_ksiL0, L0, diffL0, 3);
113  CHKERRQ(ierr);
114  ierr = base_polynomials(p, ksiL1, diff_ksiL1, L1, diffL1, 3);
115  CHKERRQ(ierr);
116  ierr = base_polynomials(p, ksiL2, diff_ksiL2, L2, diffL2, 3);
117  CHKERRQ(ierr);
118  } else {
119  ierr = base_polynomials(p, ksiL0, NULL, L0, NULL, 3);
120  CHKERRQ(ierr);
121  ierr = base_polynomials(p, ksiL1, NULL, L1, NULL, 3);
122  CHKERRQ(ierr);
123  ierr = base_polynomials(p, ksiL2, NULL, L2, NULL, 3);
124  CHKERRQ(ierr);
125  }
126  int shift = ii * P;
127  int jj = 0;
128  int oo = 0;
129  for (; oo <= p; oo++) {
130  int pp0 = 0;
131  for (; pp0 <= oo; pp0++) {
132  int pp1 = 0;
133  for (; (pp0 + pp1) <= oo; pp1++) {
134  int pp2 = oo - pp0 - pp1;
135  if (pp2 >= 0) {
136  if (L2N != NULL) {
137  L2N[shift + jj] = L0[pp0] * L1[pp1] * L2[pp2];
138  }
139  if (diff_L2N != NULL) {
140  int dd = 0;
141  for (; dd < 3; dd++) {
142  diff_L2N[3 * shift + 3 * jj + dd] =
143  diffL0[dd * (p + 1) + pp0] * L1[pp1] * L2[pp2] +
144  L0[pp0] * diffL1[dd * (p + 1) + pp1] * L2[pp2] +
145  L0[pp0] * L1[pp1] * diffL2[dd * (p + 1) + pp2];
146  }
147  }
148  jj++;
149  }
150  }
151  }
152  }
153  if (jj != P)
154  SETERRQ2(PETSC_COMM_SELF, 1, "wrong order %d != %d", jj, P);
155  }
157 }
#define NBVOLUMETET_L2(P)
Number of base functions on tetrahedron for L2 space.
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
static PetscErrorCode ierr
Definition: l2.c:27
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
CHKERRQ(ierr)
const int N
Definition: speed_test.cpp:3
field with C-1 continuity
Definition: definitions.h:171

◆ L2_Ainsworth_ShapeFunctions_MBTRI()

PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTRI ( int  p,
double N,
double diffN,
double L2N,
double diff_L2N,
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)

Get base functions on triangle for L2 space.

Parameters
ppolynomial order
Nbarycentric coordinates (shape functions) at integration points
diffNdirevatives of barycentric coordinates, i.e. direvatives of shape functions
L2Nvalues of L2 base at integration points
diff_L2Ndirvatives of base functions at integration points
GDIMnumber of integration points
base_polynomialspolynomial base used to construct L2 base on element
Returns
PetscErrorCode

Definition at line 29 of file l2.c.

33  {
35 
36  int P = NBFACETRI_L2(p);
37  if (P == 0)
39  double diff_ksiL01[2], diff_ksiL20[2];
40  int dd = 0;
41  for (; dd < 2; dd++) {
42  diff_ksiL01[dd] = (diffN[1 * 2 + dd] - diffN[0 * 2 + dd]);
43  diff_ksiL20[dd] = (diffN[2 * 2 + dd] - diffN[0 * 2 + dd]);
44  }
45  int ii = 0;
46  for (; ii != GDIM; ++ii) {
47  int node_shift = ii * 3;
48  double ksiL01 = N[node_shift + 1] - N[node_shift + 0];
49  double ksiL20 = N[node_shift + 2] - N[node_shift + 0];
50  double L01[p + 1], L20[p + 1];
51  double diffL01[2 * (p + 1)], diffL20[2 * (p + 1)];
52  ierr = base_polynomials(p, ksiL01, diff_ksiL01, L01, diffL01, 2);
53  CHKERRQ(ierr);
54  ierr = base_polynomials(p, ksiL20, diff_ksiL20, L20, diffL20, 2);
55  CHKERRQ(ierr);
56  int shift = ii * P;
57  int jj = 0;
58  int oo = 0;
59  for (; oo <= p; oo++) {
60  int pp0 = 0;
61  for (; pp0 <= oo; pp0++) {
62  int pp1 = oo - pp0;
63  if (pp1 >= 0) {
64  if (L2N != NULL) {
65  L2N[shift + jj] = L01[pp0] * L20[pp1];
66  }
67  if (diff_L2N != NULL) {
68  int dd = 0;
69  for (; dd < 2; dd++) {
70  diff_L2N[2 * shift + 2 * jj + dd] =
71  diffL01[dd * (p + 1) + pp0] * L20[pp1] +
72  L01[pp0] * diffL20[dd * (p + 1) + pp1];
73  }
74  }
75  jj++;
76  }
77  }
78  }
79  if (jj != P)
80  SETERRQ1(PETSC_COMM_SELF, 1, "wrong order %d", jj);
81  }
83 }
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
static PetscErrorCode ierr
Definition: l2.c:27
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
CHKERRQ(ierr)
#define NBFACETRI_L2(P)
Number of base functions on triangle for L2 space.
const int N
Definition: speed_test.cpp:3

◆ L2_ShapeFunctions_MBTET()

DEPRECATED PetscErrorCode L2_ShapeFunctions_MBTET ( int  p,
double N,
double diffN,
double L2N,
double diff_L2N,
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)
Deprecated:
Use L2_Ainsworth_ShapeFunctions_MBTET

Definition at line 185 of file l2.c.

189  {
190  return L2_ShapeFunctions_MBTET(p, N, diffN, L2N, diff_L2N, GDIM,
191  base_polynomials);
192 }
const int N
Definition: speed_test.cpp:3
PetscErrorCode L2_ShapeFunctions_MBTET(int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: l2.c:185

◆ L2_ShapeFunctions_MBTRI()

DEPRECATED PetscErrorCode L2_ShapeFunctions_MBTRI ( int  p,
double N,
double diffN,
double L2N,
double diff_L2N,
int  GDIM,
PetscErrorCode(*)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)  base_polynomials 
)
Deprecated:
Use L2_Ainsworth_ShapeFunctions_MBTRI

Definition at line 176 of file l2.c.

180  {
181  return L2_Ainsworth_ShapeFunctions_MBTRI(p, N, diffN, L2N, diff_L2N, GDIM,
182  base_polynomials);
183 }
PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTRI(int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Get base functions on triangle for L2 space.
Definition: l2.c:29
const int N
Definition: speed_test.cpp:3

◆ L2_VolumeShapeDiffMBTETinvJ()

PetscErrorCode L2_VolumeShapeDiffMBTETinvJ ( int  base_p,
int  p,
double volume_diffN,
double invJac,
double volume_diffNinvJac,
int  GDIM 
)

Definition at line 158 of file l2.c.

161  {
163  int ii, gg;
164  for (ii = 0; ii != NBVOLUMETET_L2(p); ++ii) {
165  for (gg = 0; gg < GDIM; gg++) {
166  int shift1 = NBVOLUMETET_L2(base_p) * gg;
167  int shift2 = NBVOLUMETET_L2(p) * gg;
168  cblas_dgemv(CblasRowMajor, CblasTrans, 3, 3, 1., invJac, 3,
169  &(volume_diffN)[3 * shift1 + 3 * ii], 1, 0.,
170  &(volume_diffNinvJac)[3 * shift2 + 3 * ii], 1);
171  }
172  }
174 }
void cblas_dgemv(const enum CBLAS_ORDER order, const enum CBLAS_TRANSPOSE TransA, const int M, const int N, const double alpha, const double *A, const int lda, const double *X, const int incX, const double beta, double *Y, const int incY)
Definition: cblas_dgemv.c:11
#define NBVOLUMETET_L2(P)
Number of base functions on tetrahedron for L2 space.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:491
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498