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