Actual source code: snesj2.c
petsc-dev 2014-02-02
2: #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/
3: #include <petscdm.h> /*I "petscdm.h" I*/
7: /*@C
8: SNESComputeJacobianDefaultColor - Computes the Jacobian using
9: finite differences and coloring to exploit matrix sparsity.
11: Collective on SNES
13: Input Parameters:
14: + snes - nonlinear solver object
15: . x1 - location at which to evaluate Jacobian
16: - ctx - MatFDColoring context or NULL
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: .notes: If the coloring is not provided through the context, this will first try to get the
26: coloring from the DM. If the DM type has no coloring routine, then it will try to
27: get the coloring from the matrix. This requires that the matrix have nonzero entries
28: precomputed. This is discouraged, as MatColoringApply() is not parallel by default.
30: .keywords: SNES, finite differences, Jacobian, coloring, sparse
32: .seealso: SNESSetJacobian(), SNESTestJacobian(), SNESComputeJacobianDefault()
33: MatFDColoringCreate(), MatFDColoringSetFunction()
35: @*/
37: PetscErrorCode SNESComputeJacobianDefaultColor(SNES snes,Vec x1,Mat *J,Mat *B,MatStructure *flag,void *ctx)
38: {
39: MatFDColoring color = (MatFDColoring)ctx;
41: DM dm;
42: PetscErrorCode (*func)(SNES,Vec,Vec,void*);
43: Vec F;
44: void *funcctx;
45: MatColoring mc;
46: ISColoring iscoloring;
47: PetscBool hascolor;
48: PetscBool solvec,matcolor = PETSC_FALSE;
52: else {PetscObjectQuery((PetscObject)*B,"SNESMatFDColoring",(PetscObject*)&color);}
53: *flag = SAME_NONZERO_PATTERN;
54: SNESGetFunction(snes,&F,&func,&funcctx);
55: if (!color) {
56: SNESGetDM(snes,&dm);
57: DMHasColoring(dm,&hascolor);
58: matcolor = PETSC_FALSE;
59: PetscOptionsGetBool(((PetscObject)snes)->prefix,"-snes_fd_color_use_mat",&matcolor,NULL);
60: if (hascolor && !matcolor) {
61: DMCreateColoring(dm,IS_COLORING_GLOBAL,&iscoloring);
62: MatFDColoringCreate(*B,iscoloring,&color);
63: MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))func,funcctx);
64: MatFDColoringSetFromOptions(color);
65: MatFDColoringSetUp(*B,iscoloring,color);
66: ISColoringDestroy(&iscoloring);
67: } else {
68: MatColoringCreate(*B,&mc);
69: MatColoringSetDistance(mc,2);
70: MatColoringSetType(mc,MATCOLORINGSL);
71: MatColoringSetFromOptions(mc);
72: MatColoringApply(mc,&iscoloring);
73: MatColoringDestroy(&mc);
74: MatFDColoringCreate(*B,iscoloring,&color);
75: MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))func,(void*)funcctx);
76: MatFDColoringSetFromOptions(color);
77: MatFDColoringSetUp(*B,iscoloring,color);
78: ISColoringDestroy(&iscoloring);
79: }
80: PetscObjectCompose((PetscObject)*B,"SNESMatFDColoring",(PetscObject)color);
81: PetscObjectDereference((PetscObject)color);
82: }
84: /* F is only usable if there is no RHS on the SNES and the full solution corresponds to x1 */
85: VecEqual(x1,snes->vec_sol,&solvec);
86: if (!snes->vec_rhs && solvec) {
87: MatFDColoringSetF(color,F);
88: }
89: MatFDColoringApply(*B,color,x1,flag,snes);
90: if (*J != *B) {
91: MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
92: MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
93: }
94: return(0);
95: }