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 /* Call general purpose update function */ 37 if (snes->ops->update) { 38 ierr = (*snes->ops->update)(snes, 0);CHKERRQ(ierr); 39 } 40 41 /* Solve J Y = F, where J is Jacobian matrix */ 42 ierr = SNESComputeJacobian(snes,X,snes->jacobian,snes->jacobian_pre,&flg);CHKERRQ(ierr); 43 if (flg == SAME_PRECONDITIONER) {ierr = KSPSetReusePreconditioner(snes->ksp,PETSC_TRUE);CHKERRQ(ierr);} 44 else {ierr = KSPSetReusePreconditioner(snes->ksp,PETSC_FALSE);CHKERRQ(ierr);} 45 ierr = KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 46 ierr = KSPSolve(snes->ksp,F,Y);CHKERRQ(ierr); 47 ierr = KSPGetConvergedReason(snes->ksp,&kspreason);CHKERRQ(ierr); 48 if (kspreason < 0 && ++snes->numLinearSolveFailures >= snes->maxLinearSolveFailures) { 49 ierr = PetscInfo2(snes,"iter=%D, number linear solve failures %D greater than current SNES allowed, stopping solve\n",snes->iter,snes->numLinearSolveFailures);CHKERRQ(ierr); 50 snes->reason = SNES_DIVERGED_LINEAR_SOLVE; 51 } else snes->reason = SNES_CONVERGED_ITS; 52 53 ierr = KSPGetIterationNumber(snes->ksp,&lits);CHKERRQ(ierr); 54 snes->linear_its += lits; 55 ierr = PetscInfo2(snes,"iter=%D, linear solve iterations=%D\n",snes->iter,lits);CHKERRQ(ierr); 56 snes->iter++; 57 58 /* Take the computed step. */ 59 ierr = VecAXPY(X,-1.0,Y);CHKERRQ(ierr); 60 if (snes->numbermonitors) { 61 PetscReal fnorm; 62 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 63 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 64 ierr = SNESMonitor(snes,1,fnorm);CHKERRQ(ierr); 65 } 66 PetscFunctionReturn(0); 67 } 68 69 #undef __FUNCT__ 70 #define __FUNCT__ "SNESSetUp_KSPONLY" 71 static PetscErrorCode SNESSetUp_KSPONLY(SNES snes) 72 { 73 PetscErrorCode ierr; 74 75 PetscFunctionBegin; 76 ierr = SNESSetUpMatrices(snes);CHKERRQ(ierr); 77 PetscFunctionReturn(0); 78 } 79 80 #undef __FUNCT__ 81 #define __FUNCT__ "SNESDestroy_KSPONLY" 82 static PetscErrorCode SNESDestroy_KSPONLY(SNES snes) 83 { 84 85 PetscFunctionBegin; 86 PetscFunctionReturn(0); 87 } 88 89 /* -------------------------------------------------------------------------- */ 90 /*MC 91 SNESKSPONLY - Nonlinear solver that only performs one Newton step and does not compute any norms. 92 The main purpose of this solver is to solve linear problems using the SNES interface, without 93 any additional overhead in the form of vector operations. 94 95 Level: beginner 96 97 .seealso: SNESCreate(), SNES, SNESSetType(), SNESNEWTONLS, SNESNEWTONTR 98 M*/ 99 #undef __FUNCT__ 100 #define __FUNCT__ "SNESCreate_KSPONLY" 101 PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes) 102 { 103 104 PetscFunctionBegin; 105 snes->ops->setup = SNESSetUp_KSPONLY; 106 snes->ops->solve = SNESSolve_KSPONLY; 107 snes->ops->destroy = SNESDestroy_KSPONLY; 108 snes->ops->setfromoptions = 0; 109 snes->ops->view = 0; 110 snes->ops->reset = 0; 111 112 snes->usesksp = PETSC_TRUE; 113 snes->usespc = PETSC_FALSE; 114 115 snes->data = 0; 116 PetscFunctionReturn(0); 117 } 118