DMProjectField#
This projects the given function of the input fields into the function space provided, putting the coefficients in a global vector.
Synopsis#
#include "petscdm.h"
#include "petscdmplex.h"
#include "petscksp.h"
PetscErrorCode DMProjectField(DM dm, PetscReal time, Vec U, void (**funcs)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), InsertMode mode, Vec X)
Collective
Input Parameters#
dm - The
DM
time - The time
U - The input field vector
funcs - The functions to evaluate, one per field
mode - The insertion mode for values
Output Parameter#
X - The output vector
Calling sequence of func
#
void funcs(PetscInt dim, PetscInt Nf, PetscInt NfAux,
const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f[]);
dim - The spatial dimension
Nf - The number of input fields
NfAux - The number of input auxiliary fields
uOff - The offset of each field in u[]
uOff_x - The offset of each field in u_x[]
u - The field values at this point in space
u_t - The field time derivative at this point in space (or
NULL
)u_x - The field derivatives at this point in space
aOff - The offset of each auxiliary field in u[]
aOff_x - The offset of each auxiliary field in u_x[]
a - The auxiliary field values at this point in space
a_t - The auxiliary field time derivative at this point in space (or
NULL
)a_x - The auxiliary field derivatives at this point in space
t - The current time
x - The coordinates of this point
numConstants - The number of constants
constants - The value of each constant
f - The value of the function at this point in space
Note#
There are three different DM
s that potentially interact in this function. The output DM
, dm, specifies the layout of the values calculates by funcs.
The input DM
, attached to U, may be different. For example, you can input the solution over the full domain, but output over a piece of the boundary, or
a subdomain. You can also output a different number of fields than the input, with different discretizations. Last the auxiliary DM
, attached to the
auxiliary field vector, which is attached to dm, can also be different. It can have a different topology, number of fields, and discretizations.
See Also#
KSP: Linear System Solvers, DM
, DMProjectFieldLocal()
, DMProjectFieldLabelLocal()
, DMProjectFunction()
, DMComputeL2Diff()
Level#
advanced
Location#
Examples#
src/snes/tutorials/ex12.c
src/snes/tutorials/ex13.c
src/ts/tutorials/ex76.c
Index of all KSP routines
Table of Contents for all manual pages
Index of all manual pages