PetscFECreateTabulation#
Tabulates the basis functions, and perhaps derivatives, at the points provided.
Synopsis#
#include "petscfe.h"
PetscErrorCode PetscFECreateTabulation(PetscFE fem, PetscInt nrepl, PetscInt npoints, const PetscReal points[], PetscInt K, PetscTabulation *T)
Not Collective
Input Parameters#
fem - The
PetscFE
objectnrepl - The number of replicas
npoints - The number of tabulation points in a replica
points - The tabulation point coordinates
K - The number of derivatives calculated
Output Parameter#
T - The basis function values and derivatives at tabulation points
Note#
T->T[0] = B[(p*pdim + i)*Nc + c] is the value at point p for basis function i and component c
T->T[1] = D[((p*pdim + i)*Nc + c)*dim + d] is the derivative value at point p for basis function i, component c, in direction d
T->T[2] = H[(((p*pdim + i)*Nc + c)*dim + d)*dim + e] is the value at point p for basis function i, component c, in directions d and e
.seealso: `PetscTabulation`, `PetscFEGetCellTabulation()`, `PetscTabulationDestroy()`
@*/
PetscErrorCode PetscFECreateTabulation(PetscFE fem, PetscInt nrepl, PetscInt npoints, const PetscReal points[], PetscInt K, PetscTabulation *T)
{
DM dm;
PetscDualSpace Q;
PetscInt Nb; /* Dimension of FE space P */
PetscInt Nc; /* Field components */
PetscInt cdim; /* Reference coordinate dimension */
PetscInt k;
PetscFunctionBegin;
if (!npoints || !fem->dualSpace || K < 0) {
*T = NULL;
PetscFunctionReturn(PETSC_SUCCESS);
}
PetscValidHeaderSpecific(fem, PETSCFE_CLASSID, 1);
PetscValidRealPointer(points, 4);
PetscValidPointer(T, 6);
PetscCall(PetscFEGetDualSpace(fem, &Q));
PetscCall(PetscDualSpaceGetDM(Q, &dm));
PetscCall(DMGetDimension(dm, &cdim));
PetscCall(PetscDualSpaceGetDimension(Q, &Nb));
PetscCall(PetscFEGetNumComponents(fem, &Nc));
PetscCall(PetscMalloc1(1, T));
(*T)->K = !cdim ? 0 : K;
(*T)->Nr = nrepl;
(*T)->Np = npoints;
(*T)->Nb = Nb;
(*T)->Nc = Nc;
(*T)->cdim = cdim;
PetscCall(PetscMalloc1((*T)->K + 1, &(*T)->T));
for (k = 0; k <= (*T)->K; ++k) PetscCall(PetscMalloc1(nrepl * npoints * Nb * Nc * PetscPowInt(cdim, k), &(*T)->T[k]));
PetscUseTypeMethod(fem, createtabulation, nrepl * npoints, points, K, *T);
PetscFunctionReturn(PETSC_SUCCESS);
}
/*@C
PetscFEComputeTabulation - Tabulates the basis functions, and perhaps derivatives, at the points provided.
Not Collective
Input Parameters:
+ fem - The `PetscFE` object
. npoints - The number of tabulation points
. points - The tabulation point coordinates
. K - The number of derivatives calculated
- T - An existing tabulation object with enough allocated space
Output Parameter:
. T - The basis function values and derivatives at tabulation points
Note:
.vb
T->T[0] = B[(p*pdim + i)*Nc + c] is the value at point p for basis function i and component c
T->T[1] = D[((p*pdim + i)*Nc + c)*dim + d] is the derivative value at point p for basis function i, component c, in direction d
T->T[2] = H[(((p*pdim + i)*Nc + c)*dim + d)*dim + e] is the value at point p for basis function i, component c, in directions d and e
See Also#
PetscTabulation
, PetscFEGetCellTabulation()
, PetscTabulationDestroy()
Level#
intermediate
Location#
Implementations#
PetscFECreateTabulation_Basic in src/dm/dt/fe/impls/basic/febasic.c
PetscFECreateTabulation_Composite in src/dm/dt/fe/impls/composite/fecomposite.c
Index of all FE routines
Table of Contents for all manual pages
Index of all manual pages