1b79b07cfSJed Brown 2af0996ceSBarry Smith #include <petsc/private/snesimpl.h> 3b79b07cfSJed Brown 4b79b07cfSJed Brown #undef __FUNCT__ 5b79b07cfSJed Brown #define __FUNCT__ "SNESSolve_KSPONLY" 6b79b07cfSJed Brown static PetscErrorCode SNESSolve_KSPONLY(SNES snes) 7b79b07cfSJed Brown { 8b79b07cfSJed Brown PetscErrorCode ierr; 9b79b07cfSJed Brown PetscInt lits; 10b79b07cfSJed Brown Vec Y,X,F; 11b79b07cfSJed Brown 12b79b07cfSJed Brown PetscFunctionBegin; 136c4ed002SBarry Smith if (snes->xl || snes->xu || snes->ops->computevariablebounds) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name); 14c579b300SPatrick Farrell 15b79b07cfSJed Brown snes->numFailures = 0; 16b79b07cfSJed Brown snes->numLinearSolveFailures = 0; 17b79b07cfSJed Brown snes->reason = SNES_CONVERGED_ITERATING; 18b79b07cfSJed Brown snes->iter = 0; 19b79b07cfSJed Brown snes->norm = 0.0; 20b79b07cfSJed Brown 21b79b07cfSJed Brown X = snes->vec_sol; 22b79b07cfSJed Brown F = snes->vec_func; 23b79b07cfSJed Brown Y = snes->vec_sol_update; 24b79b07cfSJed Brown 25b79b07cfSJed Brown ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 26a5eaddd9SBarry Smith if (snes->numbermonitors) { 27a5eaddd9SBarry Smith PetscReal fnorm; 28a5eaddd9SBarry Smith ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 29a5eaddd9SBarry Smith ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 30a5eaddd9SBarry Smith } 31b79b07cfSJed Brown 325341784dSBarry Smith /* Call general purpose update function */ 335341784dSBarry Smith if (snes->ops->update) { 345341784dSBarry Smith ierr = (*snes->ops->update)(snes, 0);CHKERRQ(ierr); 355341784dSBarry Smith } 365341784dSBarry Smith 37b79b07cfSJed Brown /* Solve J Y = F, where J is Jacobian matrix */ 38d1e9a80fSBarry Smith ierr = SNESComputeJacobian(snes,X,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 3923ee1639SBarry Smith ierr = KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 40b79b07cfSJed Brown ierr = KSPSolve(snes->ksp,F,Y);CHKERRQ(ierr); 41422a814eSBarry Smith snes->reason = SNES_CONVERGED_ITS; 42422a814eSBarry Smith SNESCheckKSPSolve(snes); 431aa26658SKarl Rupp 44b79b07cfSJed Brown ierr = KSPGetIterationNumber(snes->ksp,&lits);CHKERRQ(ierr); 45b79b07cfSJed Brown snes->linear_its += lits; 46b79b07cfSJed Brown ierr = PetscInfo2(snes,"iter=%D, linear solve iterations=%D\n",snes->iter,lits);CHKERRQ(ierr); 47b79b07cfSJed Brown snes->iter++; 48b79b07cfSJed Brown 49b79b07cfSJed Brown /* Take the computed step. */ 50b79b07cfSJed Brown ierr = VecAXPY(X,-1.0,Y);CHKERRQ(ierr); 51a5eaddd9SBarry Smith if (snes->numbermonitors) { 52a5eaddd9SBarry Smith PetscReal fnorm; 53fef909d3SLawrence Mitchell ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 54a5eaddd9SBarry Smith ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 55a5eaddd9SBarry Smith ierr = SNESMonitor(snes,1,fnorm);CHKERRQ(ierr); 56a5eaddd9SBarry Smith } 57b79b07cfSJed Brown PetscFunctionReturn(0); 58b79b07cfSJed Brown } 59b79b07cfSJed Brown 60b79b07cfSJed Brown #undef __FUNCT__ 61b79b07cfSJed Brown #define __FUNCT__ "SNESSetUp_KSPONLY" 62b79b07cfSJed Brown static PetscErrorCode SNESSetUp_KSPONLY(SNES snes) 63b79b07cfSJed Brown { 646cab3a1bSJed Brown PetscErrorCode ierr; 65b79b07cfSJed Brown 66b79b07cfSJed Brown PetscFunctionBegin; 676cab3a1bSJed Brown ierr = SNESSetUpMatrices(snes);CHKERRQ(ierr); 68b79b07cfSJed Brown PetscFunctionReturn(0); 69b79b07cfSJed Brown } 70b79b07cfSJed Brown 71b79b07cfSJed Brown #undef __FUNCT__ 72b79b07cfSJed Brown #define __FUNCT__ "SNESDestroy_KSPONLY" 73b79b07cfSJed Brown static PetscErrorCode SNESDestroy_KSPONLY(SNES snes) 74b79b07cfSJed Brown { 75b79b07cfSJed Brown 76b79b07cfSJed Brown PetscFunctionBegin; 77b79b07cfSJed Brown PetscFunctionReturn(0); 78b79b07cfSJed Brown } 79b79b07cfSJed Brown 80b79b07cfSJed Brown /* -------------------------------------------------------------------------- */ 81b79b07cfSJed Brown /*MC 82b79b07cfSJed Brown SNESKSPONLY - Nonlinear solver that only performs one Newton step and does not compute any norms. 83b79b07cfSJed Brown The main purpose of this solver is to solve linear problems using the SNES interface, without 84b79b07cfSJed Brown any additional overhead in the form of vector operations. 85b79b07cfSJed Brown 86b79b07cfSJed Brown Level: beginner 87b79b07cfSJed Brown 8804d7464bSBarry Smith .seealso: SNESCreate(), SNES, SNESSetType(), SNESNEWTONLS, SNESNEWTONTR 89b79b07cfSJed Brown M*/ 90b79b07cfSJed Brown #undef __FUNCT__ 91b79b07cfSJed Brown #define __FUNCT__ "SNESCreate_KSPONLY" 928cc058d9SJed Brown PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes) 93b79b07cfSJed Brown { 94b79b07cfSJed Brown 95b79b07cfSJed Brown PetscFunctionBegin; 96b79b07cfSJed Brown snes->ops->setup = SNESSetUp_KSPONLY; 97b79b07cfSJed Brown snes->ops->solve = SNESSolve_KSPONLY; 98b79b07cfSJed Brown snes->ops->destroy = SNESDestroy_KSPONLY; 9937596af1SLisandro Dalcin snes->ops->setfromoptions = 0; 10037596af1SLisandro Dalcin snes->ops->view = 0; 10137596af1SLisandro Dalcin snes->ops->reset = 0; 102b79b07cfSJed Brown 10342f4f86dSBarry Smith snes->usesksp = PETSC_TRUE; 10442f4f86dSBarry Smith snes->usespc = PETSC_FALSE; 10542f4f86dSBarry Smith 106*4fc747eaSLawrence Mitchell snes->alwayscomputesfinalresidual = PETSC_FALSE; 107*4fc747eaSLawrence Mitchell 108b79b07cfSJed Brown snes->data = 0; 109b79b07cfSJed Brown PetscFunctionReturn(0); 110b79b07cfSJed Brown } 111