1 2 #include <petsc-private/snesimpl.h> 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "SNESSolve_KSPONLY" 6 static PetscErrorCode SNESSolve_KSPONLY(SNES snes) 7 { 8 PetscErrorCode ierr; 9 PetscInt lits; 10 MatStructure flg = DIFFERENT_NONZERO_PATTERN; 11 Vec Y,X,F; 12 KSPConvergedReason kspreason; 13 14 PetscFunctionBegin; 15 snes->numFailures = 0; 16 snes->numLinearSolveFailures = 0; 17 snes->reason = SNES_CONVERGED_ITERATING; 18 snes->iter = 0; 19 snes->norm = 0.0; 20 21 X = snes->vec_sol; 22 F = snes->vec_func; 23 Y = snes->vec_sol_update; 24 25 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 26 if (snes->domainerror) { 27 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 28 PetscFunctionReturn(0); 29 } 30 if (snes->numbermonitors) { 31 PetscReal fnorm; 32 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 33 ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 34 } 35 36 /* Solve J Y = F, where J is Jacobian matrix */ 37 ierr = SNESComputeJacobian(snes,X,&snes->jacobian,&snes->jacobian_pre,&flg);CHKERRQ(ierr); 38 ierr = KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre,flg);CHKERRQ(ierr); 39 ierr = KSPSolve(snes->ksp,F,Y);CHKERRQ(ierr); 40 ierr = KSPGetConvergedReason(snes->ksp,&kspreason);CHKERRQ(ierr); 41 if (kspreason < 0 && ++snes->numLinearSolveFailures >= snes->maxLinearSolveFailures) { 42 ierr = PetscInfo2(snes,"iter=%D, number linear solve failures %D greater than current SNES allowed, stopping solve\n",snes->iter,snes->numLinearSolveFailures);CHKERRQ(ierr); 43 snes->reason = SNES_DIVERGED_LINEAR_SOLVE; 44 } else { 45 snes->reason = SNES_CONVERGED_ITS; 46 } 47 ierr = KSPGetIterationNumber(snes->ksp,&lits);CHKERRQ(ierr); 48 snes->linear_its += lits; 49 ierr = PetscInfo2(snes,"iter=%D, linear solve iterations=%D\n",snes->iter,lits);CHKERRQ(ierr); 50 snes->iter++; 51 52 /* Take the computed step. */ 53 ierr = VecAXPY(X,-1.0,Y);CHKERRQ(ierr); 54 if (snes->numbermonitors) { 55 PetscReal fnorm; 56 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 57 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 58 ierr = SNESMonitor(snes,1,fnorm);CHKERRQ(ierr); 59 } 60 PetscFunctionReturn(0); 61 } 62 63 #undef __FUNCT__ 64 #define __FUNCT__ "SNESSetUp_KSPONLY" 65 static PetscErrorCode SNESSetUp_KSPONLY(SNES snes) 66 { 67 PetscErrorCode ierr; 68 69 PetscFunctionBegin; 70 ierr = SNESSetUpMatrices(snes);CHKERRQ(ierr); 71 PetscFunctionReturn(0); 72 } 73 74 #undef __FUNCT__ 75 #define __FUNCT__ "SNESDestroy_KSPONLY" 76 static PetscErrorCode SNESDestroy_KSPONLY(SNES snes) 77 { 78 79 PetscFunctionBegin; 80 PetscFunctionReturn(0); 81 } 82 83 /* -------------------------------------------------------------------------- */ 84 /*MC 85 SNESKSPONLY - Nonlinear solver that only performs one Newton step and does not compute any norms. 86 The main purpose of this solver is to solve linear problems using the SNES interface, without 87 any additional overhead in the form of vector operations. 88 89 Level: beginner 90 91 .seealso: SNESCreate(), SNES, SNESSetType(), SNESNEWTONLS, SNESNEWTONTR 92 M*/ 93 EXTERN_C_BEGIN 94 #undef __FUNCT__ 95 #define __FUNCT__ "SNESCreate_KSPONLY" 96 PetscErrorCode SNESCreate_KSPONLY(SNES snes) 97 { 98 99 PetscFunctionBegin; 100 snes->ops->setup = SNESSetUp_KSPONLY; 101 snes->ops->solve = SNESSolve_KSPONLY; 102 snes->ops->destroy = SNESDestroy_KSPONLY; 103 snes->ops->setfromoptions = 0; 104 snes->ops->view = 0; 105 snes->ops->reset = 0; 106 107 snes->usesksp = PETSC_TRUE; 108 snes->usespc = PETSC_FALSE; 109 110 snes->data = 0; 111 PetscFunctionReturn(0); 112 } 113 EXTERN_C_END 114