1 2 #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/ 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "SNESDefaultComputeJacobianColor" 6 /*@C 7 SNESDefaultComputeJacobianColor - Computes the Jacobian using 8 finite differences and coloring to exploit matrix sparsity. 9 10 Collective on SNES 11 12 Input Parameters: 13 + snes - nonlinear solver object 14 . x1 - location at which to evaluate Jacobian 15 - ctx - coloring context, where ctx must have type MatFDColoring, 16 as created via MatFDColoringCreate() 17 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 22 23 Level: intermediate 24 25 .notes: If ctx is not provided, this will try to get the coloring from the DM. If the DM type 26 has no coloring routine, then it will try to get the coloring from the matrix. This 27 requires that the matrix have nonzero entries precomputed, such as in 28 snes/examples/tutorials/ex45.c. 29 30 .keywords: SNES, finite differences, Jacobian, coloring, sparse 31 32 .seealso: SNESSetJacobian(), SNESTestJacobian(), SNESDefaultComputeJacobian() 33 MatFDColoringCreate(), MatFDColoringSetFunction() 34 35 @*/ 36 37 PetscErrorCode SNESDefaultComputeJacobianColor(SNES snes,Vec x1,Mat *J,Mat *B,MatStructure *flag,void *ctx) 38 { 39 MatFDColoring color = (MatFDColoring)ctx; 40 PetscErrorCode ierr; 41 DM dm; 42 PetscErrorCode (*func)(SNES,Vec,Vec,void*); 43 Vec F; 44 void *funcctx; 45 ISColoring iscoloring; 46 PetscBool hascolor; 47 48 PetscFunctionBegin; 49 if (color) { 50 PetscValidHeaderSpecific(color,MAT_FDCOLORING_CLASSID,6); 51 } else { 52 ierr = PetscObjectQuery((PetscObject)*B,"SNESMatFDColoring",(PetscObject *)&color);CHKERRQ(ierr); 53 } 54 *flag = SAME_NONZERO_PATTERN; 55 ierr = SNESGetFunction(snes,&F,&func,&funcctx); 56 if (!color) { 57 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 58 ierr = DMHasColoring(dm,&hascolor);CHKERRQ(ierr); 59 if (hascolor) { 60 ierr = DMCreateColoring(dm,IS_COLORING_GLOBAL,MATAIJ,&iscoloring);CHKERRQ(ierr); 61 ierr = MatFDColoringCreate(*B,iscoloring,&color);CHKERRQ(ierr); 62 ierr = ISColoringDestroy(&iscoloring);CHKERRQ(ierr); 63 ierr = MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))func,funcctx);CHKERRQ(ierr); 64 ierr = MatFDColoringSetFromOptions(color);CHKERRQ(ierr); 65 } else { 66 ierr = MatGetColoring(*B,MATCOLORINGSL,&iscoloring);CHKERRQ(ierr); 67 ierr = MatFDColoringCreate(*B,iscoloring,&color);CHKERRQ(ierr); 68 ierr = ISColoringDestroy(&iscoloring);CHKERRQ(ierr); 69 ierr = MatFDColoringSetFunction(color,(PetscErrorCode (*)(void))func,(void*)funcctx);CHKERRQ(ierr); 70 ierr = MatFDColoringSetFromOptions(color);CHKERRQ(ierr); 71 } 72 ierr = PetscObjectCompose((PetscObject)*B,"SNESMatFDColoring",(PetscObject)color);CHKERRQ(ierr); 73 ierr = PetscObjectDereference((PetscObject)color);CHKERRQ(ierr); 74 } 75 ierr = MatFDColoringSetF(color,PETSC_NULL);CHKERRQ(ierr); 76 ierr = MatFDColoringApply(*B,color,x1,flag,snes);CHKERRQ(ierr); 77 if (*J != *B) { 78 ierr = MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 79 ierr = MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 80 } 81 PetscFunctionReturn(0); 82 } 83