Actual source code: snesj2.c
petsc-3.10.5 2019-03-28
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: }
14: /*@C
15: SNESComputeJacobianDefaultColor - Computes the Jacobian using
16: finite differences and coloring to exploit matrix sparsity.
18: Collective on SNES
20: Input Parameters:
21: + snes - nonlinear solver object
22: . x1 - location at which to evaluate Jacobian
23: - ctx - MatFDColoring context or NULL
25: Output Parameters:
26: + J - Jacobian matrix (not altered in this routine)
27: - B - newly computed Jacobian matrix to use with preconditioner (generally the same as J)
29: Level: intermediate
31: Options Database Key:
32: + -snes_fd_color_use_mat - use a matrix coloring from the explicit matrix nonzero pattern instead of from the DM providing the matrix
33: . -snes_fd_color - Activates SNESComputeJacobianDefaultColor() in SNESSetFromOptions()
34: . -mat_fd_coloring_err <err> - Sets <err> (square root of relative error in the function)
35: . -mat_fd_coloring_umin <umin> - Sets umin, the minimum allowable u-value magnitude
36: - -mat_fd_type - Either wp or ds (see MATMFFD_WP or MATMFFD_DS)
38: Notes:
39: If the coloring is not provided through the context, this will first try to get the
40: coloring from the DM. If the DM type has no coloring routine, then it will try to
41: get the coloring from the matrix. This requires that the matrix have nonzero entries
42: precomputed.
44: .keywords: SNES, finite differences, Jacobian, coloring, sparse
46: .seealso: SNESSetJacobian(), SNESTestJacobian(), SNESComputeJacobianDefault()
47: MatFDColoringCreate(), MatFDColoringSetFunction()
49: @*/
51: PetscErrorCode SNESComputeJacobianDefaultColor(SNES snes,Vec x1,Mat J,Mat B,void *ctx)
52: {
53: MatFDColoring color = (MatFDColoring)ctx;
55: DM dm;
56: MatColoring mc;
57: ISColoring iscoloring;
58: PetscBool hascolor;
59: PetscBool solvec,matcolor = PETSC_FALSE;
63: if (!color) {PetscObjectQuery((PetscObject)B,"SNESMatFDColoring",(PetscObject*)&color);}
65: if (!color) {
66: SNESGetDM(snes,&dm);
67: DMHasColoring(dm,&hascolor);
68: matcolor = PETSC_FALSE;
69: PetscOptionsGetBool(((PetscObject)snes)->options,((PetscObject)snes)->prefix,"-snes_fd_color_use_mat",&matcolor,NULL);
70: if (hascolor && !matcolor) {
71: DMCreateColoring(dm,IS_COLORING_GLOBAL,&iscoloring);
72: MatFDColoringCreate(B,iscoloring,&color);
73: MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))SNESComputeFunctionCtx,NULL);
74: MatFDColoringSetFromOptions(color);
75: MatFDColoringSetUp(B,iscoloring,color);
76: ISColoringDestroy(&iscoloring);
77: } else {
78: MatColoringCreate(B,&mc);
79: MatColoringSetDistance(mc,2);
80: MatColoringSetType(mc,MATCOLORINGSL);
81: MatColoringSetFromOptions(mc);
82: MatColoringApply(mc,&iscoloring);
83: MatColoringDestroy(&mc);
84: MatFDColoringCreate(B,iscoloring,&color);
85: MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))SNESComputeFunctionCtx,NULL);
86: MatFDColoringSetFromOptions(color);
87: MatFDColoringSetUp(B,iscoloring,color);
88: ISColoringDestroy(&iscoloring);
89: }
90: PetscObjectCompose((PetscObject)B,"SNESMatFDColoring",(PetscObject)color);
91: PetscObjectDereference((PetscObject)color);
92: }
94: /* F is only usable if there is no RHS on the SNES and the full solution corresponds to x1 */
95: VecEqual(x1,snes->vec_sol,&solvec);
96: if (!snes->vec_rhs && solvec) {
97: Vec F;
98: SNESGetFunction(snes,&F,NULL,NULL);
99: MatFDColoringSetF(color,F);
100: }
101: MatFDColoringApply(B,color,x1,snes);
102: if (J != B) {
103: MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY);
104: MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY);
105: }
106: return(0);
107: }