Actual source code: snesj2.c
2: #include <petsc/private/snesimpl.h>
3: #include <petscdm.h>
5: /*
6: MatFDColoringSetFunction() takes a function with four arguments, we want to use SNESComputeFunction()
7: since it logs function computation information.
8: */
9: static PetscErrorCode SNESComputeFunctionCtx(SNES snes,Vec x,Vec f,void *ctx)
10: {
11: return SNESComputeFunction(snes,x,f);
12: }
13: static PetscErrorCode SNESComputeMFFunctionCtx(SNES snes,Vec x,Vec f,void *ctx)
14: {
15: return SNESComputeMFFunction(snes,x,f);
16: }
18: /*@C
19: SNESComputeJacobianDefaultColor - Computes the Jacobian using
20: finite differences and coloring to exploit matrix sparsity.
22: Collective on SNES
24: Input Parameters:
25: + snes - nonlinear solver object
26: . x1 - location at which to evaluate Jacobian
27: - ctx - MatFDColoring context or NULL
29: Output Parameters:
30: + J - Jacobian matrix (not altered in this routine)
31: - B - newly computed Jacobian matrix to use with preconditioner (generally the same as J)
33: Level: intermediate
35: Options Database Key:
36: + -snes_fd_color_use_mat - use a matrix coloring from the explicit matrix nonzero pattern instead of from the DM providing the matrix
37: . -snes_fd_color - Activates SNESComputeJacobianDefaultColor() in SNESSetFromOptions()
38: . -mat_fd_coloring_err <err> - Sets <err> (square root of relative error in the function)
39: . -mat_fd_coloring_umin <umin> - Sets umin, the minimum allowable u-value magnitude
40: . -mat_fd_type - Either wp or ds (see MATMFFD_WP or MATMFFD_DS)
41: . -snes_mf_operator - Use matrix free application of Jacobian
42: - -snes_mf - Use matrix free Jacobian with no explicit Jacobian representation
44: Notes:
45: If the coloring is not provided through the context, this will first try to get the
46: coloring from the DM. If the DM type has no coloring routine, then it will try to
47: get the coloring from the matrix. This requires that the matrix have nonzero entries
48: precomputed.
50: SNES supports three approaches for computing (approximate) Jacobians: user provided via SNESSetJacobian(), matrix free via SNESSetUseMatrixFree(),
51: and computing explicitly with finite differences and coloring using MatFDColoring. It is also possible to use automatic differentiation
52: and the MatFDColoring object, see src/ts/tutorials/autodiff/ex16adj_tl.cxx
54: .seealso: SNESSetJacobian(), SNESTestJacobian(), SNESComputeJacobianDefault(), SNESSetUseMatrixFree(),
55: MatFDColoringCreate(), MatFDColoringSetFunction()
57: @*/
58: PetscErrorCode SNESComputeJacobianDefaultColor(SNES snes,Vec x1,Mat J,Mat B,void *ctx)
59: {
60: MatFDColoring color = (MatFDColoring)ctx;
62: DM dm;
63: MatColoring mc;
64: ISColoring iscoloring;
65: PetscBool hascolor;
66: PetscBool solvec,matcolor = PETSC_FALSE;
67: DMSNES dms;
71: if (!color) {PetscObjectQuery((PetscObject)B,"SNESMatFDColoring",(PetscObject*)&color);}
73: if (!color) {
74: SNESGetDM(snes,&dm);
75: DMHasColoring(dm,&hascolor);
76: matcolor = PETSC_FALSE;
77: PetscOptionsGetBool(((PetscObject)snes)->options,((PetscObject)snes)->prefix,"-snes_fd_color_use_mat",&matcolor,NULL);
78: if (hascolor && !matcolor) {
79: DMCreateColoring(dm,IS_COLORING_GLOBAL,&iscoloring);
80: } else {
81: MatColoringCreate(B,&mc);
82: MatColoringSetDistance(mc,2);
83: MatColoringSetType(mc,MATCOLORINGSL);
84: MatColoringSetFromOptions(mc);
85: MatColoringApply(mc,&iscoloring);
86: MatColoringDestroy(&mc);
87: }
88: MatFDColoringCreate(B,iscoloring,&color);
89: DMGetDMSNES(dm,&dms);
90: if (dms->ops->computemffunction) {
91: MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))SNESComputeMFFunctionCtx,NULL);
92: } else {
93: MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))SNESComputeFunctionCtx,NULL);
94: }
95: MatFDColoringSetFromOptions(color);
96: MatFDColoringSetUp(B,iscoloring,color);
97: ISColoringDestroy(&iscoloring);
98: PetscObjectCompose((PetscObject)B,"SNESMatFDColoring",(PetscObject)color);
99: PetscObjectDereference((PetscObject)color);
100: }
102: /* F is only usable if there is no RHS on the SNES and the full solution corresponds to x1 */
103: VecEqual(x1,snes->vec_sol,&solvec);
104: if (!snes->vec_rhs && solvec) {
105: Vec F;
106: SNESGetFunction(snes,&F,NULL,NULL);
107: MatFDColoringSetF(color,F);
108: }
109: MatFDColoringApply(B,color,x1,snes);
110: if (J != B) {
111: MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY);
112: MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY);
113: }
114: return(0);
115: }