:orphan: # PetscFEIntegrateResidual Produce the element residual vector for a chunk of elements by quadrature integration ## Synopsis ``` #include "petscfe.h" PetscErrorCode PetscFEIntegrateResidual(PetscDS ds, PetscFormKey key, PetscInt Ne, PetscFEGeom *cgeom, const PetscScalar coefficients[], const PetscScalar coefficients_t[], PetscDS probAux, const PetscScalar coefficientsAux[], PetscReal t, PetscScalar elemVec[]) ``` Not Collective ## Input Parameters - ***ds -*** The `PetscDS` specifying the discretizations and continuum functions - ***key -*** The (label+value, field) being integrated - ***Ne -*** The number of elements in the chunk - ***cgeom -*** The cell geometry for each cell in the chunk - ***coefficients -*** The array of FEM basis coefficients for the elements - ***coefficients_t -*** The array of FEM basis time derivative coefficients for the elements - ***probAux -*** The `PetscDS` specifying the auxiliary discretizations - ***coefficientsAux -*** The array of FEM auxiliary basis coefficients for the elements - ***t -*** The time ## Output Parameter - ***elemVec -*** the element residual vectors from each element ## Note ```none Loop over batch of elements (e): Loop over quadrature points (q): Make u_q and gradU_q (loops over fields,Nb,Ncomp) and x_q Call f_0 and f_1 Loop over element vector entries (f,fc --> i): elemVec[i] += \psi^{fc}_f(q) f0_{fc}(u, \nabla u) + \nabla\psi^{fc}_f(q) \cdot f1_{fc,df}(u, \nabla u) ``` ## See Also `PetscFEIntegrateResidual()` ## Level intermediate ## Location src/dm/dt/fe/interface/fe.c ## Implementations PetscFEIntegrateResidual_Basic in src/dm/dt/fe/impls/basic/febasic.c
PetscFEIntegrateResidual_OpenCL in src/dm/dt/fe/impls/opencl/feopencl.c
--- [Edit on GitLab](https://gitlab.com/petsc/petsc/-/edit/release/src/dm/dt/fe/interface/fe.c) [Index of all FE routines](index.md) [Table of Contents for all manual pages](/manualpages/index.md) [Index of all manual pages](/manualpages/singleindex.md)