Actual source code: snesj2.c
2: #include <private/snesimpl.h> /*I "petscsnes.h" I*/
6: /*@C
7: SNESDefaultComputeJacobianColor - Computes the Jacobian using
8: finite differences and coloring to exploit matrix sparsity.
9:
10: Collective on SNES
12: Input Parameters:
13: + snes - nonlinear solver object
14: . x1 - location at which to evaluate Jacobian
15: - ctx - coloring context, where ctx must have type MatFDColoring,
16: as created via MatFDColoringCreate()
18: Output Parameters:
19: + J - Jacobian matrix (not altered in this routine)
20: . B - newly computed Jacobian matrix to use with preconditioner (generally the same as J)
21: - flag - flag indicating whether the matrix sparsity structure has changed
23: Level: intermediate
25: .keywords: SNES, finite differences, Jacobian, coloring, sparse
27: .seealso: SNESSetJacobian(), SNESTestJacobian(), SNESDefaultComputeJacobian()
28: TSDefaultComputeJacobianColor(), MatFDColoringCreate(),
29: MatFDColoringSetFunction()
31: @*/
33: PetscErrorCode SNESDefaultComputeJacobianColor(SNES snes,Vec x1,Mat *J,Mat *B,MatStructure *flag,void *ctx)
34: {
35: MatFDColoring color = (MatFDColoring) ctx;
37: Vec f;
38: PetscErrorCode (*ff)(void),(*fd)(void);
42: *flag = SAME_NONZERO_PATTERN;
43: SNESGetFunction(snes,&f,(PetscErrorCode (**)(SNES,Vec,Vec,void*))&ff,0);
44: MatFDColoringGetFunction(color,&fd,PETSC_NULL);
45: if (fd == ff) { /* reuse function value computed in SNES */
46: MatFDColoringSetF(color,f);
47: }
48: MatFDColoringApply(*B,color,x1,flag,snes);
49: if (*J != *B) {
50: MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
51: MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
52: }
53: return(0);
54: }