CeedBasis
-
static int CeedChebyshevPolynomialsAtPoint(CeedScalar x, CeedInt n, CeedScalar *chebyshev_x)
Compute Chebyshev polynomial values at a point.
Library Developer Functions
- Parameters:
x – [in] Coordinate to evaluate Chebyshev polynomials at
n – [in] Number of Chebyshev polynomials to evaluate,
n >= 2
chebyshev_x – [out] Array of Chebyshev polynomial values
- Returns:
An error code: 0 - success, otherwise - failure
-
static int CeedChebyshevDerivativeAtPoint(CeedScalar x, CeedInt n, CeedScalar *chebyshev_dx)
Compute values of the derivative of Chebyshev polynomials at a point.
Library Developer Functions
- Parameters:
x – [in] Coordinate to evaluate derivative of Chebyshev polynomials at
n – [in] Number of Chebyshev polynomials to evaluate,
n >= 2
chebyshev_dx – [out] Array of Chebyshev polynomial derivative values
- Returns:
An error code: 0 - success, otherwise - failure
-
static int CeedHouseholderReflect(CeedScalar *A, const CeedScalar *v, CeedScalar b, CeedInt m, CeedInt n, CeedInt row, CeedInt col)
Compute Householder reflection.
Computes \(A = (I - b v v^T) A\), where \(A\) is an \(m \times n\) matrix indexed as
A[i*row + j*col]
.Library Developer Functions
- Parameters:
A – [inout] Matrix to apply Householder reflection to, in place
v – [in] Householder vector
b – [in] Scaling factor
m – [in] Number of rows in
A
n – [in] Number of columns in
A
row – [in] Row stride
col – [in] Col stride
- Returns:
An error code: 0 - success, otherwise - failure
-
static int CeedGivensRotation(CeedScalar *A, CeedScalar c, CeedScalar s, CeedTransposeMode t_mode, CeedInt i, CeedInt k, CeedInt m, CeedInt n)
Compute Givens rotation.
Computes \(A = G A\) (or \(G^T A\) in transpose mode), where \(A\) is an \(m \times n\) matrix indexed as
A[i*n + j*m]
.Library Developer Functions
- Parameters:
A – [inout] Row major matrix to apply Givens rotation to, in place
c – [in] Cosine factor
s – [in] Sine factor
t_mode – [in] CEED_NOTRANSPOSE to rotate the basis counter-clockwise, which has the effect of rotating columns of
A
clockwise; CEED_TRANSPOSE for the opposite rotationi – [in] First row/column to apply rotation
k – [in] Second row/column to apply rotation
m – [in] Number of rows in
A
n – [in] Number of columns in
A
- Returns:
An error code: 0 - success, otherwise - failure
-
static int CeedScalarView(const char *name, const char *fp_fmt, CeedInt m, CeedInt n, const CeedScalar *a, FILE *stream)
View an array stored in a
CeedBasis
Library Developer Functions
- Parameters:
name – [in] Name of array
fp_fmt – [in] Printing format
m – [in] Number of rows in array
n – [in] Number of columns in array
a – [in] Array to be viewed
stream – [in] Stream to view to, e.g.,
stdout
- Returns:
An error code: 0 - success, otherwise - failure
-
static int CeedBasisCreateProjectionMatrices(CeedBasis basis_from, CeedBasis basis_to, CeedScalar **interp_project, CeedScalar **grad_project)
Create the interpolation and gradient matrices for projection from the nodes of
basis_from
to the nodes ofbasis_to
.The interpolation is given by
interp_project = interp_to^+ * interp_from
, where the pseudoinverseinterp_to^+
is given by QR factorization. The gradient is given bygrad_project = interp_to^+ * grad_from
, and is only computed for \(H^1\) spaces otherwise it should not be used.Note:
basis_from
andbasis_to
must have compatible quadrature spaces.Library Developer Functions
- Parameters:
basis_from – [in]
CeedBasis
to project frombasis_to – [in]
CeedBasis
to project tointerp_project – [out] Address of the variable where the newly created interpolation matrix will be stored
grad_project – [out] Address of the variable where the newly created gradient matrix will be stored
- Returns:
An error code: 0 - success, otherwise - failure