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;
61: DM dm;
62: MatColoring mc;
63: ISColoring iscoloring;
64: PetscBool hascolor;
65: PetscBool solvec,matcolor = PETSC_FALSE;
66: DMSNES dms;
69: if (!color) PetscObjectQuery((PetscObject)B,"SNESMatFDColoring",(PetscObject*)&color);
71: if (!color) {
72: SNESGetDM(snes,&dm);
73: DMHasColoring(dm,&hascolor);
74: matcolor = PETSC_FALSE;
75: PetscOptionsGetBool(((PetscObject)snes)->options,((PetscObject)snes)->prefix,"-snes_fd_color_use_mat",&matcolor,NULL);
76: if (hascolor && !matcolor) {
77: DMCreateColoring(dm,IS_COLORING_GLOBAL,&iscoloring);
78: } else {
79: MatColoringCreate(B,&mc);
80: MatColoringSetDistance(mc,2);
81: MatColoringSetType(mc,MATCOLORINGSL);
82: MatColoringSetFromOptions(mc);
83: MatColoringApply(mc,&iscoloring);
84: MatColoringDestroy(&mc);
85: }
86: MatFDColoringCreate(B,iscoloring,&color);
87: DMGetDMSNES(dm,&dms);
88: if (dms->ops->computemffunction) {
89: MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))SNESComputeMFFunctionCtx,NULL);
90: } else {
91: MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))SNESComputeFunctionCtx,NULL);
92: }
93: MatFDColoringSetFromOptions(color);
94: MatFDColoringSetUp(B,iscoloring,color);
95: ISColoringDestroy(&iscoloring);
96: PetscObjectCompose((PetscObject)B,"SNESMatFDColoring",(PetscObject)color);
97: PetscObjectDereference((PetscObject)color);
98: }
100: /* F is only usable if there is no RHS on the SNES and the full solution corresponds to x1 */
101: VecEqual(x1,snes->vec_sol,&solvec);
102: if (!snes->vec_rhs && solvec) {
103: Vec F;
104: SNESGetFunction(snes,&F,NULL,NULL);
105: MatFDColoringSetF(color,F);
106: }
107: MatFDColoringApply(B,color,x1,snes);
108: if (J != B) {
109: MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY);
110: MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY);
111: }
112: return 0;
113: }