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 object

  • nrepl - 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#

src/dm/dt/fe/interface/fe.c

Implementations#

PetscFECreateTabulation_Basic in src/dm/dt/fe/impls/basic/febasic.c
PetscFECreateTabulation_Composite in src/dm/dt/fe/impls/composite/fecomposite.c


Edit on GitLab

Index of all FE routines
Table of Contents for all manual pages
Index of all manual pages