:orphan:
# TSSetIJacobian
Set the function to compute the matrix dF/dU + a*dF/dU_t where F(t,U,U_t) is the function provided with `TSSetIFunction()`.
## Synopsis
```
#include "petscts.h"
PetscErrorCode TSSetIJacobian(TS ts, Mat Amat, Mat Pmat, TSIJacobian f, void *ctx)
```
Logically Collective
## Input Parameters
- ***ts -*** the `TS` context obtained from `TSCreate()`
- ***Amat -*** (approximate) matrix to store Jacobian entries computed by `f`
- ***Pmat -*** matrix used to compute preconditioner (usually the same as `Amat`)
- ***f -*** the Jacobian evaluation routine
- ***ctx -*** user-defined context for private data for the Jacobian evaluation routine (may be `NULL`)
## Calling sequence of `f`
```none
PetscErrorCode f(TS ts, PetscReal t, Vec U, Vec U_t, PetscReal a, Mat Amat, Mat Pmat, void *ctx)
```
- ***ts -*** the `TS` context obtained from `TSCreate()`
- ***t -*** time at step/stage being solved
- ***U -*** state vector
- ***U_t -*** time derivative of state vector
- ***a -*** shift
- ***Amat -*** (approximate) Jacobian of F(t,U,W+a*U), equivalent to dF/dU + a*dF/dU_t
- ***Pmat -*** matrix used for constructing preconditioner, usually the same as `Amat`
- ***ctx -*** [optional] user-defined context for matrix evaluation routine
## Notes
The matrices `Amat` and `Pmat` are exactly the matrices that are used by `SNES` for the nonlinear solve.
If you know the operator Amat has a null space you can use `MatSetNullSpace()` and `MatSetTransposeNullSpace()` to supply the null
space to `Amat` and the `KSP` solvers will automatically use that null space as needed during the solution process.
The matrix dF/dU + a*dF/dU_t you provide turns out to be
the Jacobian of F(t,U,W+a*U) where F(t,U,U_t) = 0 is the DAE to be solved.
The time integrator internally approximates U_t by W+a*U where the positive "shift"
a and vector W depend on the integration method, step size, and past states. For example with
the backward Euler method a = 1/dt and W = -a*U(previous timestep) so
W + a*U = a*(U - U(previous timestep)) = (U - U(previous timestep))/dt
You must set all the diagonal entries of the matrices, if they are zero you must still set them with a zero value
The TS solver may modify the nonzero structure and the entries of the matrices `Amat` and `Pmat` between the calls to `f`
You should not assume the values are the same in the next call to `f` as you set them in the previous call.
## See Also
[](ch_ts), `TS`, `TSSetIFunction()`, `TSSetRHSJacobian()`, `SNESComputeJacobianDefaultColor()`, `SNESComputeJacobianDefault()`, `TSSetRHSFunction()`
## Level
beginner
## Location
src/ts/interface/ts.c
## Examples
src/ts/tutorials/ex10.c
src/ts/tutorials/ex14.c
src/ts/tutorials/ex15.c
src/ts/tutorials/ex16.c
src/ts/tutorials/ex17.c
src/ts/tutorials/ex19.c
src/ts/tutorials/ex20.c
src/ts/tutorials/ex20adj.c
src/ts/tutorials/ex20fwd.c
src/ts/tutorials/ex20opt_ic.c
src/ts/tutorials/ex20opt_p.c
---
[Edit on GitLab](https://gitlab.com/petsc/petsc/-/edit/release/src/ts/interface/ts.c)
[Index of all TS routines](index.md)
[Table of Contents for all manual pages](/manualpages/index.md)
[Index of all manual pages](/manualpages/singleindex.md)