Actual source code: snesj2.c
2: #include src/mat/matimpl.h
3: #include src/snes/snesimpl.h
7: /*@C
8: SNESDefaultComputeJacobianColor - Computes the Jacobian using
9: finite differences and coloring to exploit matrix sparsity.
10:
11: Collective on SNES
13: Input Parameters:
14: + snes - nonlinear solver object
15: . x1 - location at which to evaluate Jacobian
16: - ctx - coloring context, where ctx must have type MatFDColoring,
17: as created via MatFDColoringCreate()
19: Output Parameters:
20: + J - Jacobian matrix (not altered in this routine)
21: . B - newly computed Jacobian matrix to use with preconditioner (generally the same as J)
22: - flag - flag indicating whether the matrix sparsity structure has changed
24: Options Database Keys:
25: . -mat_fd_coloring_freq <freq> - Activates SNESDefaultComputeJacobianColor()
27: Level: intermediate
29: .keywords: SNES, finite differences, Jacobian, coloring, sparse
31: .seealso: SNESSetJacobian(), SNESTestJacobian(), SNESDefaultComputeJacobian()
32: TSDefaultComputeJacobianColor(), MatFDColoringCreate(),
33: MatFDColoringSetFunction()
35: @*/
36: PetscErrorCode SNESDefaultComputeJacobianColor(SNES snes,Vec x1,Mat *J,Mat *B,MatStructure *flag,void *ctx)
37: {
38: MatFDColoring color = (MatFDColoring) ctx;
40: PetscInt freq,it;
41: Vec f;
44: MatFDColoringGetFrequency(color,&freq);
45: SNESGetIterationNumber(snes,&it);
47: if ((freq > 1) && ((it % freq))) {
48: PetscLogInfo(color,"SNESDefaultComputeJacobianColor:Skipping Jacobian recomputation, it %D, freq %D\n",it,freq);
49: *flag = SAME_PRECONDITIONER;
50: } else {
51: PetscLogInfo(color,"SNESDefaultComputeJacobianColor:Computing Jacobian, it %D, freq %D\n",it,freq);
52: *flag = SAME_NONZERO_PATTERN;
53: SNESGetFunction(snes,&f,0,0);
54: MatFDColoringSetF(color,f);
55: MatFDColoringApply(*B,color,x1,flag,snes);
56: }
57: if (*J != *B) {
58: MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
59: MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
60: }
61: return(0);
62: }