#include "petscts.h" PetscErrorCode TSSetIJacobian(TS ts,Mat A,Mat B,TSIJacobian f,void *ctx)Logically Collective on TS
ts | - the TS context obtained from TSCreate() | |
A | - Jacobian matrix | |
B | - preconditioning matrix for A (may be same as A) | |
f | - the Jacobian evaluation routine | |
ctx | - user-defined context for private data for the Jacobian evaluation routine (may be PETSC_NULL) |
f(TS ts,PetscReal t,Vec U,Vec U_t,PetscReal a,Mat *A,Mat *B,MatStructure *flag,void *ctx);
t | - time at step/stage being solved | |
U | - state vector | |
U_t | - time derivative of state vector | |
a | - shift | |
A | - Jacobian of G(U) = F(t,U,W+a*U), equivalent to dF/dU + a*dF/dU_t | |
B | - preconditioning matrix for A, may be same as A | |
flag | - flag indicating information about the preconditioner matrix structure (same as flag in KSPSetOperators()) | |
ctx | - [optional] user-defined context for matrix evaluation routine |
The matrix dF/dU + a*dF/dU_t you provide turns out to be the Jacobian of G(U) = 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
Level:beginner
Location:src/ts/interface/ts.c
Index of all TS routines
Table of Contents for all manual pages
Index of all manual pages