TSSetIJacobian

Set the function to compute the Jacobian of G(U) = F(t,U,U0+a*U) where F(t,U,U_t) = 0 is the DAE to be solved.

Synopsis

#include "petscts.h"  
PetscErrorCode PETSCTS_DLLEXPORT TSSetIJacobian(TS ts,Mat A,Mat B,TSIJacobian f,void *ctx)
Collective on TS

Input Parameters

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)

Calling sequence of f

 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,U0+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

Notes

The matrices A and B are exactly the matrices that are used by SNES for the nonlinear solve.

Keywords

TS, timestep, DAE, Jacobian

See Also

TSSetIFunction(), TSSetRHSJacobian()

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

Examples

src/ts/examples/tutorials/ex8.c.html