v0.9.0
Classes | Functions
Base functions

Calculation of base functions at integration points. More...

Collaboration diagram for Base functions:

Classes

struct  MoFEM::BaseFunctionCtx
 Base class used to exchange data between element data structures and class calculating base functions. More...
 
struct  MoFEM::BaseFunction
 Base class if inherited used to calculate base functions. More...
 
struct  MoFEM::EdgePolynomialBase
 Calculate base functions on tetrahedral. More...
 
struct  MoFEM::EntPolynomialBaseCtx
 Class used to pass element data to calculate base functions on tet,triangle,edge. More...
 
struct  MoFEM::FatPrismPolynomialBaseCtx
 Class used to pass element data to calculate base functions on fat prism. More...
 
struct  MoFEM::FatPrismPolynomialBase
 Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (that not perfect) More...
 
struct  MoFEM::FlatPrismPolynomialBaseCtx
 Class used to pass element data to calculate base functions on flat prism. More...
 
struct  MoFEM::FlatPrismPolynomialBase
 Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (that not perfect) More...
 
struct  MoFEM::JacobiPolynomialCtx
 Class used to give arguments to Legendre base functions. More...
 
struct  MoFEM::JacobiPolynomial
 Calculating Legendre base functions. More...
 
struct  MoFEM::LegendrePolynomialCtx
 Class used to give arguments to Legendre base functions. More...
 
struct  MoFEM::LegendrePolynomial
 Calculating Legendre base functions. More...
 
struct  MoFEM::LobattoPolynomialCtx
 Class used to give arguments to Lobatto base functions. More...
 
struct  MoFEM::LobattoPolynomial
 Calculating Lobatto base functions. More...
 
struct  MoFEM::KernelLobattoPolynomialCtx
 Class used to give arguments to Kernel Lobatto base functions. More...
 
struct  MoFEM::KernelLobattoPolynomial
 Calculating Lobatto base functions. More...
 
struct  MoFEM::QuadPolynomialBase
 Calculate base functions on triangle. More...
 
struct  MoFEM::TetPolynomialBase
 Calculate base functions on tetrahedral. More...
 
struct  MoFEM::TriPolynomialBase
 Calculate base functions on triangle. More...
 

Functions

PetscErrorCode Legendre_polynomials (int p, double s, double *diff_s, double *L, double *diffL, const int dim)
 Calculate Legendre approximation basis. More...
 
PetscErrorCode Lobatto_polynomials (int p, double s, double *diff_s, double *L, double *diffL, const int dim)
 Calculate Lobatto base functions. More...
 
PetscErrorCode LobattoKernel_polynomials (int p, double s, double *diff_s, double *L, double *diffL, const int dim)
 Calculate Kernel Lobatto base functions. More...
 

Detailed Description

Calculation of base functions at integration points.

Function Documentation

◆ Legendre_polynomials()

PetscErrorCode Legendre_polynomials ( int  p,
double  s,
double diff_s,
double L,
double diffL,
const int  dim 
)

Calculate Legendre approximation basis.

MoFEM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with MoFEM. If not, see http://www.gnu.org/licenses/

Lagrange polynomial is given by

\[ L_0(s)=1;\quad L_1(s) = s \]

and following terms are generated inductively

\[ L_{l+1}=\frac{2l+1}{l+1}sL_l(s)-\frac{l}{l+1}L_{l-1}(s) \]

Note that:

\[ s\in[-1,1] \quad \textrm{and}\; s=s(\xi_0,\xi_1,\xi_2) \]

where \(\xi_i\) are barycentric coordinates of element.

Parameters
pis approximation order
sis position \(s\in[-1,1]\)
diff_sderivatives of shape functions, i.e. \(\frac{\partial s}{\partial \xi_i}\)
Return values
Lapproximation functions
diffLderivatives, i.e. \(\frac{\partial L}{\partial \xi_i}\)
Parameters
dimdimension
Returns
error code
Examples
edge_and_bubble_shape_functions_on_quad.cpp.

Definition at line 25 of file base_functions.c.

26  {
28  if (dim < 1)
29  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim < 1");
30  if (dim > 3)
31  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim > 3");
32  if (p < 0)
33  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "p < 0");
34  L[0] = 1;
35  if (diffL != NULL) {
36  diffL[0 * (p + 1) + 0] = 0;
37  if (dim >= 2) {
38  diffL[1 * (p + 1) + 0] = 0;
39  if (dim == 3) {
40  diffL[2 * (p + 1) + 0] = 0;
41  }
42  }
43  }
44  if (p == 0)
46  L[1] = s;
47  if (diffL != NULL) {
48  if (diff_s == NULL) {
49  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "diff_s == NULL");
50  }
51  diffL[0 * (p + 1) + 1] = diff_s[0];
52  if (dim >= 2) {
53  diffL[1 * (p + 1) + 1] = diff_s[1];
54  if (dim == 3) {
55  diffL[2 * (p + 1) + 1] = diff_s[2];
56  }
57  }
58  }
59  if (p == 1)
61  int l = 1;
62  for (; l < p; l++) {
63  double A = ((2 * (double)l + 1) / ((double)l + 1));
64  double B = ((double)l / ((double)l + 1));
65  L[l + 1] = A * s * L[l] - B * L[l - 1];
66  if (diffL != NULL) {
67  diffL[0 * (p + 1) + l + 1] =
68  A * (s * diffL[0 * (p + 1) + l] + diff_s[0] * L[l]) -
69  B * diffL[0 * (p + 1) + l - 1];
70  if (dim >= 2) {
71  diffL[1 * (p + 1) + l + 1] =
72  A * (s * diffL[1 * (p + 1) + l] + diff_s[1] * L[l]) -
73  B * diffL[1 * (p + 1) + l - 1];
74  if (dim == 3) {
75  diffL[2 * (p + 1) + l + 1] =
76  A * (s * diffL[2 * (p + 1) + l] + diff_s[2] * L[l]) -
77  B * diffL[2 * (p + 1) + l - 1];
78  }
79  }
80  }
81  }
83 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ Lobatto_polynomials()

PetscErrorCode Lobatto_polynomials ( int  p,
double  s,
double diff_s,
double L,
double diffL,
const int  dim 
)

Calculate Lobatto base functions.

Order of first function is 2.

Parameters
pis approximation order
sis position \(s\in[-1,1]\)
diff_sderivatives of shape functions, i.e. \(\frac{\partial s}{\partial \xi_i}\)
Return values
Lapproximation functions
diffLderivatives, i.e. \(\frac{\partial L}{\partial \xi_i}\)
Parameters
dimdimension
Returns
error code

Definition at line 211 of file base_functions.c.

212  {
213 
215  if (dim < 1)
216  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim < 1");
217  if (dim > 3)
218  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim > 3");
219  if (p < 2)
220  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "p < 2");
221  p -= 2; // polynomial order starts from 2
222  double l[p + 2];
223  ierr = Legendre_polynomials(p + 1, s, NULL, l, NULL, 1);
224  CHKERRQ(ierr);
225  {
226  // Derivatives
227  int k = 0;
228  for (; k <= p; k++) {
229  if (diffL != NULL) {
230  if (diff_s == NULL) {
231  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "diff_s == NULL");
232  }
233  double a = l[k + 1];
234  diffL[0 * (p + 1) + k] = a * diff_s[0];
235  if (dim >= 2) {
236  diffL[1 * (p + 1) + k] = a * diff_s[1];
237  if (dim == 3) {
238  diffL[2 * (p + 1) + k] = a * diff_s[2];
239  }
240  }
241  }
242  }
243  }
244  {
245  // Functions
246  bzero(L, (p + 1) * sizeof(double));
247  int nb_gauss_pts = QUAD_1D_TABLE[p + 2]->npoints;
248  double *points = QUAD_1D_TABLE[p + 2]->points;
249  double *weights = QUAD_1D_TABLE[p + 2]->weights;
250  s = s + 1;
251  int gg = 0;
252  for (; gg != nb_gauss_pts; gg++) {
253  double ksi = points[2 * gg + 1];
254  double zeta = s * ksi - 1;
255  ierr = Legendre_polynomials(p + 1, zeta, NULL, l, NULL, 1);
256  CHKERRQ(ierr);
257  double w = s * weights[gg];
258  cblas_daxpy(p + 1, w, &l[1], 1, &L[0], 1);
259  }
260  }
261  {
262  int k = 0;
263  for (; k <= p; k++) {
264  double a = 4 * sqrt(k + 2 - 0.5);
265  if (L != NULL)
266  L[k] *= a;
267  if (diffL != NULL)
268  diffL[k] *= a;
269  }
270  }
272 }
PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Legendre approximation basis.
double * points
Definition: quad.h:30
static QUAD *const QUAD_1D_TABLE[]
Definition: quad.h:164
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:501
int npoints
Definition: quad.h:29
static PetscErrorCode ierr
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
CHKERRQ(ierr)
double * weights
Definition: quad.h:31

◆ LobattoKernel_polynomials()

PetscErrorCode LobattoKernel_polynomials ( int  p,
double  s,
double diff_s,
double L,
double diffL,
const int  dim 
)

Calculate Kernel Lobatto base functions.

This is implemented using definitions from Hermes2d https://github.com/hpfem/hermes following book by Pavel Solin et al [solin2003higher].

Parameters
pis approximation order
sis position \(s\in[-1,1]\)
diff_sderivatives of shape functions, i.e. \(\frac{\partial s}{\partial \xi_i}\)
Return values
Lapproximation functions
diffLderivatives, i.e. \(\frac{\partial L}{\partial \xi_i}\)
Parameters
dimdimension
Returns
error code

Definition at line 369 of file base_functions.c.

371  {
373  if (dim < 1)
374  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim < 1");
375  if (dim > 3)
376  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim > 3");
377  if (p < 0)
378  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "p < 0");
379  if (p > 9)
380  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
381  "Polynomial beyond order 9 is not implemented");
382  if (L) {
383  int l = 0;
384  for (; l != p + 1; l++) {
385  L[l] = f_phi[l](s);
386  }
387  }
388  if (diffL != NULL) {
389  if (diff_s == NULL) {
390  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "diff_s == NULL");
391  }
392  int l = 0;
393  for (; l != p + 1; l++) {
394  double a = f_phix[l](s);
395  diffL[0 * (p + 1) + l] = diff_s[0] * a;
396  if (dim >= 2) {
397  diffL[1 * (p + 1) + l] = diff_s[1] * a;
398  if (dim == 3) {
399  diffL[2 * (p + 1) + l] = diff_s[2] * a;
400  }
401  }
402  }
403  }
405 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
static double(* f_phix[])(double x)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
static double(* f_phi[])(double x)