v0.14.0
Functions | Variables
base_functions.c File Reference
#include <cblas.h>
#include <petscsys.h>
#include <phg-quadrule/quad.h>
#include <definitions.h>
#include <base_functions.h>

Go to the source code of this file.

Functions

PetscErrorCode Legendre_polynomials (int p, double s, double *diff_s, double *L, double *diffL, const int dim)
 Calculate Legendre approximation basis. More...
 
PetscErrorCode Jacobi_polynomials (int p, double alpha, double x, double t, double *diff_x, double *diff_t, double *L, double *diffL, const int dim)
 Calculate Jacobi approximation basis. More...
 
PetscErrorCode IntegratedJacobi_polynomials (int p, double alpha, double x, double t, double *diff_x, double *diff_t, double *L, double *diffL, const int dim)
 Calculate integrated Jacobi approximation basis. More...
 
PetscErrorCode Lobatto_polynomials (int p, double s, double *diff_s, double *L, double *diffL, const int dim)
 Calculate Lobatto base functions [28]. More...
 
static double f_phi0 (double x)
 
static double f_phi1 (double x)
 
static double f_phi2 (double x)
 
static double f_phi3 (double x)
 
static double f_phi4 (double x)
 
static double f_phi5 (double x)
 
static double f_phi6 (double x)
 
static double f_phi7 (double x)
 
static double f_phi8 (double x)
 
static double f_phi9 (double x)
 
static double f_phi0x (double x)
 
static double f_phi1x (double x)
 
static double f_phi2x (double x)
 
static double f_phi3x (double x)
 
static double f_phi4x (double x)
 
static double f_phi5x (double x)
 
static double f_phi6x (double x)
 
static double f_phi7x (double x)
 
static double f_phi8x (double x)
 
static double f_phi9x (double x)
 
PetscErrorCode LobattoKernel_polynomials (int p, double s, double *diff_s, double *L, double *diffL, const int dim)
 Calculate Kernel Lobatto base functions. More...
 

Variables

static PetscErrorCode ierr
 
static double(* f_phi [])(double x)
 
static double(* f_phix [])(double x)
 

Function Documentation

◆ f_phi0()

static double f_phi0 ( double  x)
static

Definition at line 233 of file base_functions.c.

233 { return LOBATTO_PHI0(x); }

◆ f_phi0x()

static double f_phi0x ( double  x)
static

Definition at line 247 of file base_functions.c.

247 { return LOBATTO_PHI0X(x); }

◆ f_phi1()

static double f_phi1 ( double  x)
static

Definition at line 234 of file base_functions.c.

234 { return LOBATTO_PHI1(x); }

◆ f_phi1x()

static double f_phi1x ( double  x)
static

Definition at line 248 of file base_functions.c.

248 { return LOBATTO_PHI1X(x); }

◆ f_phi2()

static double f_phi2 ( double  x)
static

Definition at line 235 of file base_functions.c.

235 { return LOBATTO_PHI2(x); }

◆ f_phi2x()

static double f_phi2x ( double  x)
static

Definition at line 249 of file base_functions.c.

249 { return LOBATTO_PHI2X(x); }

◆ f_phi3()

static double f_phi3 ( double  x)
static

Definition at line 236 of file base_functions.c.

236 { return LOBATTO_PHI3(x); }

◆ f_phi3x()

static double f_phi3x ( double  x)
static

Definition at line 250 of file base_functions.c.

250 { return LOBATTO_PHI3X(x); }

◆ f_phi4()

static double f_phi4 ( double  x)
static

Definition at line 237 of file base_functions.c.

237 { return LOBATTO_PHI4(x); }

◆ f_phi4x()

static double f_phi4x ( double  x)
static

Definition at line 251 of file base_functions.c.

251 { return LOBATTO_PHI4X(x); }

◆ f_phi5()

static double f_phi5 ( double  x)
static

Definition at line 238 of file base_functions.c.

238 { return LOBATTO_PHI5(x); }

◆ f_phi5x()

static double f_phi5x ( double  x)
static

Definition at line 252 of file base_functions.c.

252 { return LOBATTO_PHI5X(x); }

◆ f_phi6()

static double f_phi6 ( double  x)
static

Definition at line 239 of file base_functions.c.

239 { return LOBATTO_PHI6(x); }

◆ f_phi6x()

static double f_phi6x ( double  x)
static

Definition at line 253 of file base_functions.c.

253 { return LOBATTO_PHI6X(x); }

◆ f_phi7()

static double f_phi7 ( double  x)
static

Definition at line 240 of file base_functions.c.

240 { return LOBATTO_PHI7(x); }

◆ f_phi7x()

static double f_phi7x ( double  x)
static

Definition at line 254 of file base_functions.c.

254 { return LOBATTO_PHI7X(x); }

◆ f_phi8()

static double f_phi8 ( double  x)
static

Definition at line 241 of file base_functions.c.

241 { return LOBATTO_PHI8(x); }

◆ f_phi8x()

static double f_phi8x ( double  x)
static

Definition at line 255 of file base_functions.c.

255 { return LOBATTO_PHI8X(x); }

◆ f_phi9()

static double f_phi9 ( double  x)
static

Definition at line 242 of file base_functions.c.

242 { return LOBATTO_PHI9(x); }

◆ f_phi9x()

static double f_phi9x ( double  x)
static

Definition at line 256 of file base_functions.c.

256 { return LOBATTO_PHI9X(x); }

◆ IntegratedJacobi_polynomials()

PetscErrorCode IntegratedJacobi_polynomials ( int  p,
double  alpha,
double  x,
double  t,
double diff_x,
double diff_t,
double L,
double diffL,
const int  dim 
)

Calculate integrated Jacobi approximation basis.

For more details see [29]

Parameters
pis approximation order
alphapolynomial parameter
xis position \(s\in[0,t]\)
trange of polynomial
diff_xderivatives of shape functions, i.e. \(\frac{\partial x}{\partial \xi_i}\)
diff_tderivatives of shape functions, i.e. \(\frac{\partial t}{\partial \xi_i}\)
Return values
Lapproximation functions
diffLderivatives, i.e. \(\frac{\partial L}{\partial \xi_i}\)
Parameters
dimdimension
Returns
error code

Definition at line 134 of file base_functions.c.

137  {
139 #ifndef NDEBUG
140  if (dim < 1)
141  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim < 1");
142  if (dim > 3)
143  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim > 3");
144  if (p < 1)
145  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "p < 1");
146 #endif // NDEBUG
147  L[0] = x;
148  if (diffL != NULL) {
149  int d = 0;
150  for (; d != dim; ++d) {
151  diffL[d * p + 0] = diff_x[d];
152  }
153  }
154  if (p == 0)
156  double jacobi[(p + 1)];
157  double diff_jacobi[(p + 1) * dim];
158  ierr = Jacobi_polynomials(p, alpha, x, t, diff_x, diff_t, jacobi, diff_jacobi,
159  dim);
160  CHKERRQ(ierr);
161  int l = 1;
162  for (; l < p; l++) {
163  int i = l + 1;
164  const double a = (i + alpha) / ((2 * i + alpha - 1) * (2 * i + alpha));
165  const double b = alpha / ((2 * i + alpha - 2) * (2 * i + alpha));
166  const double c = (i - 1) / ((2 * i + alpha - 2) * (2 * i + alpha - 1));
167  L[l] = a * jacobi[i] + b * t * jacobi[i - 1] - c * t * t * jacobi[i - 2];
168  if (diffL != NULL) {
169  int d = 0;
170  for (; d != dim; ++d) {
171  diffL[d * p + l] = a * diff_jacobi[d * (p + 1) + i] +
172  b * (t * diff_jacobi[d * (p + 1) + i - 1] +
173  diff_t[d] * jacobi[i - 1]) -
174  c * (t * t * diff_jacobi[d * (p + 1) + i - 2] +
175  2 * t * diff_t[d] * jacobi[i - 2]);
176  }
177  }
178  }
180 }

◆ Jacobi_polynomials()

PetscErrorCode Jacobi_polynomials ( int  p,
double  alpha,
double  x,
double  t,
double diff_x,
double diff_t,
double L,
double diffL,
const int  dim 
)

Calculate Jacobi approximation basis.

For more details see [29]

Parameters
pis approximation order
alphapolynomial parameter
xis position \(s\in[0,t]\)
trange of polynomial
diff_xderivatives of shape functions, i.e. \(\frac{\partial x}{\partial \xi_i}\)
diff_tderivatives of shape functions, i.e. \(\frac{\partial t}{\partial \xi_i}\)
Return values
Lapproximation functions
diffLderivatives, i.e. \(\frac{\partial L}{\partial \xi_i}\)
Parameters
dimdimension
Returns
error code

Definition at line 67 of file base_functions.c.

69  {
71 #ifndef NDEBUG
72  if (dim < 1)
73  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim < 1");
74  if (dim > 3)
75  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "dim > 3");
76  if (p < 0)
77  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "p < 0");
78 
79  if (diffL != NULL) {
80  if (diff_x == NULL) {
81  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "diff_s == NULL");
82  }
83  }
84 
85 #endif // NDEBUG
86  L[0] = 1;
87  if (diffL != NULL) {
88  diffL[0 * (p + 1) + 0] = 0;
89  if (dim >= 2) {
90  diffL[1 * (p + 1) + 0] = 0;
91  if (dim == 3) {
92  diffL[2 * (p + 1) + 0] = 0;
93  }
94  }
95  }
96  if (p == 0)
98  L[1] = 2 * x - t + alpha * x;
99  if (diffL != NULL) {
100  int d = 0;
101  for (; d < dim; ++d) {
102  double d_t = (diff_t) ? diff_t[d] : 0;
103  diffL[d * (p + 1) + 1] = (2 + alpha) * diff_x[d] - d_t;
104  }
105  }
106  if (p == 1)
108  int l = 1;
109  for (; l < p; l++) {
110  int lp1 = l + 1;
111  double a = 2 * lp1 * (lp1 + alpha) * (2 * lp1 + alpha - 2);
112  double b = 2 * lp1 + alpha - 1;
113  double c = (2 * lp1 + alpha) * (2 * lp1 + alpha - 2);
114  double d = 2 * (lp1 + alpha - 1) * (lp1 - 1) * (2 * lp1 + alpha);
115  double A = b * (c * (2 * x - t) + alpha * alpha * t) / a;
116  double B = d * t * t / a;
117  L[lp1] = A * L[l] - B * L[l - 1];
118  if (diffL != NULL) {
119  int z = 0;
120  for (; z < dim; ++z) {
121  double d_t = (diff_t) ? diff_t[z] : 0;
122  double diffA =
123  b * (c * (2 * diff_x[z] - d_t) + alpha * alpha * d_t) / a;
124  double diffB = d * 2 * t * d_t / a;
125  diffL[z * (p + 1) + lp1] = A * diffL[z * (p + 1) + l] -
126  B * diffL[z * (p + 1) + l - 1] +
127  diffA * L[l] - diffB * L[l - 1];
128  }
129  }
130  }
132 }

Variable Documentation

◆ f_phi

double(* f_phi[])(double x)
static
Initial value:

Definition at line 244 of file base_functions.c.

◆ f_phix

double(* f_phix[])(double x)
static
Initial value:

Definition at line 258 of file base_functions.c.

◆ ierr

PetscErrorCode ierr
static

Definition at line 13 of file base_functions.c.

MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
LOBATTO_PHI4
#define LOBATTO_PHI4(x)
Definition: base_functions.h:132
f_phi5x
static double f_phi5x(double x)
Definition: base_functions.c:252
LOBATTO_PHI1X
#define LOBATTO_PHI1X(x)
Definition: base_functions.h:158
f_phi7
static double f_phi7(double x)
Definition: base_functions.c:240
LOBATTO_PHI7
#define LOBATTO_PHI7(x)
Definition: base_functions.h:141
LOBATTO_PHI9X
#define LOBATTO_PHI9X(x)
Definition: base_functions.h:179
f_phi7x
static double f_phi7x(double x)
Definition: base_functions.c:254
ierr
static PetscErrorCode ierr
Definition: base_functions.c:13
LOBATTO_PHI8
#define LOBATTO_PHI8(x)
Definition: base_functions.h:144
f_phi1x
static double f_phi1x(double x)
Definition: base_functions.c:248
f_phi4x
static double f_phi4x(double x)
Definition: base_functions.c:251
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
LOBATTO_PHI6X
#define LOBATTO_PHI6X(x)
Definition: base_functions.h:169
LOBATTO_PHI1
#define LOBATTO_PHI1(x)
Definition: base_functions.h:127
LOBATTO_PHI2X
#define LOBATTO_PHI2X(x)
Definition: base_functions.h:159
f_phi8x
static double f_phi8x(double x)
Definition: base_functions.c:255
c
const double c
speed of light (cm/ns)
Definition: initial_diffusion.cpp:39
f_phi3
static double f_phi3(double x)
Definition: base_functions.c:236
Jacobi_polynomials
PetscErrorCode Jacobi_polynomials(int p, double alpha, double x, double t, double *diff_x, double *diff_t, double *L, double *diffL, const int dim)
Calculate Jacobi approximation basis.
Definition: base_functions.c:67
f_phi4
static double f_phi4(double x)
Definition: base_functions.c:237
LOBATTO_PHI9
#define LOBATTO_PHI9(x)
Definition: base_functions.h:149
a
constexpr double a
Definition: approx_sphere.cpp:30
f_phi8
static double f_phi8(double x)
Definition: base_functions.c:241
f_phi2x
static double f_phi2x(double x)
Definition: base_functions.c:249
LOBATTO_PHI5X
#define LOBATTO_PHI5X(x)
Definition: base_functions.h:166
f_phi9
static double f_phi9(double x)
Definition: base_functions.c:242
LOBATTO_PHI8X
#define LOBATTO_PHI8X(x)
Definition: base_functions.h:175
MoFEM::L
VectorDouble L
Definition: Projection10NodeCoordsOnField.cpp:124
f_phi9x
static double f_phi9x(double x)
Definition: base_functions.c:256
f_phi0
static double f_phi0(double x)
Definition: base_functions.c:233
t
constexpr double t
plate stiffness
Definition: plate.cpp:58
f_phi6x
static double f_phi6x(double x)
Definition: base_functions.c:253
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
f_phi3x
static double f_phi3x(double x)
Definition: base_functions.c:250
f_phi6
static double f_phi6(double x)
Definition: base_functions.c:239
LOBATTO_PHI0
#define LOBATTO_PHI0(x)
Definitions taken from Hermes2d code.
Definition: base_functions.h:126
LOBATTO_PHI3
#define LOBATTO_PHI3(x)
Definition: base_functions.h:130
f_phi1
static double f_phi1(double x)
Definition: base_functions.c:234
LOBATTO_PHI4X
#define LOBATTO_PHI4X(x)
Definition: base_functions.h:163
sdf_hertz_2d_axisymm_plane.d
float d
Definition: sdf_hertz_2d_axisymm_plane.py:4
f_phi0x
static double f_phi0x(double x)
Definition: base_functions.c:247
LOBATTO_PHI2
#define LOBATTO_PHI2(x)
Definition: base_functions.h:128
LOBATTO_PHI7X
#define LOBATTO_PHI7X(x)
Definition: base_functions.h:172
f_phi5
static double f_phi5(double x)
Definition: base_functions.c:238
f_phi2
static double f_phi2(double x)
Definition: base_functions.c:235
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
LOBATTO_PHI5
#define LOBATTO_PHI5(x)
Definition: base_functions.h:135
LOBATTO_PHI6
#define LOBATTO_PHI6(x)
Definition: base_functions.h:138
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
MOFEM_INVALID_DATA
@ MOFEM_INVALID_DATA
Definition: definitions.h:36
LOBATTO_PHI3X
#define LOBATTO_PHI3X(x)
Definition: base_functions.h:161
LOBATTO_PHI0X
#define LOBATTO_PHI0X(x)
Derivatives of kernel functions for Lobbatto base.
Definition: base_functions.h:157