PetscGaussLobattoLegendreElementGradientCreate#
computes the gradient for a single 1d GLL element
Synopsis#
#include "petscdt.h"
PetscErrorCode PetscGaussLobattoLegendreElementGradientCreate(PetscInt n, PetscReal nodes[], PetscReal weights[], PetscReal ***AA, PetscReal ***AAT)
Not Collective
Input Parameters#
n - the number of GLL nodes
nodes - the GLL nodes, of length
nweights - the GLL weights, of length
n
Output Parameters#
AA - the stiffness element, of dimension
nbynAAT - the transpose of AA (pass in
NULLif you do not need this array), of dimensionnbyn
Notes#
Destroy this with PetscGaussLobattoLegendreElementGradientDestroy()
You can access entries in these arrays with AA[i][j] but in memory it is stored in contiguous memory, row-oriented
See Also#
PetscDTGaussLobattoLegendreQuadrature(), PetscGaussLobattoLegendreElementLaplacianDestroy(), PetscGaussLobattoLegendreElementGradientDestroy()
Level#
beginner
Location#
Index of all DT routines
Table of Contents for all manual pages
Index of all manual pages