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; 13c579b300SPatrick Farrell 14c579b300SPatrick Farrell if (snes->xl || snes->xu || snes->ops->computevariablebounds) { 15c579b300SPatrick Farrell SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name); 16c579b300SPatrick Farrell } 17c579b300SPatrick Farrell 18b79b07cfSJed Brown snes->numFailures = 0; 19b79b07cfSJed Brown snes->numLinearSolveFailures = 0; 20b79b07cfSJed Brown snes->reason = SNES_CONVERGED_ITERATING; 21b79b07cfSJed Brown snes->iter = 0; 22b79b07cfSJed Brown snes->norm = 0.0; 23b79b07cfSJed Brown 24b79b07cfSJed Brown X = snes->vec_sol; 25b79b07cfSJed Brown F = snes->vec_func; 26b79b07cfSJed Brown Y = snes->vec_sol_update; 27b79b07cfSJed Brown 28b79b07cfSJed Brown ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 29a5eaddd9SBarry Smith if (snes->numbermonitors) { 30a5eaddd9SBarry Smith PetscReal fnorm; 31a5eaddd9SBarry Smith ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 32a5eaddd9SBarry Smith ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 33a5eaddd9SBarry Smith } 34b79b07cfSJed Brown 355341784dSBarry Smith /* Call general purpose update function */ 365341784dSBarry Smith if (snes->ops->update) { 375341784dSBarry Smith ierr = (*snes->ops->update)(snes, 0);CHKERRQ(ierr); 385341784dSBarry Smith } 395341784dSBarry Smith 40b79b07cfSJed Brown /* Solve J Y = F, where J is Jacobian matrix */ 41d1e9a80fSBarry Smith ierr = SNESComputeJacobian(snes,X,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 4223ee1639SBarry Smith ierr = KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 43b79b07cfSJed Brown ierr = KSPSolve(snes->ksp,F,Y);CHKERRQ(ierr); 44422a814eSBarry Smith snes->reason = SNES_CONVERGED_ITS; 45422a814eSBarry Smith SNESCheckKSPSolve(snes); 461aa26658SKarl Rupp 47b79b07cfSJed Brown ierr = KSPGetIterationNumber(snes->ksp,&lits);CHKERRQ(ierr); 48b79b07cfSJed Brown snes->linear_its += lits; 49b79b07cfSJed Brown ierr = PetscInfo2(snes,"iter=%D, linear solve iterations=%D\n",snes->iter,lits);CHKERRQ(ierr); 50b79b07cfSJed Brown snes->iter++; 51b79b07cfSJed Brown 52b79b07cfSJed Brown /* Take the computed step. */ 53b79b07cfSJed Brown ierr = VecAXPY(X,-1.0,Y);CHKERRQ(ierr); 54*1989ec39SBarry Smith ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 55a5eaddd9SBarry Smith if (snes->numbermonitors) { 56a5eaddd9SBarry Smith PetscReal fnorm; 57a5eaddd9SBarry Smith ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 58a5eaddd9SBarry Smith ierr = SNESMonitor(snes,1,fnorm);CHKERRQ(ierr); 59a5eaddd9SBarry Smith } 60b79b07cfSJed Brown PetscFunctionReturn(0); 61b79b07cfSJed Brown } 62b79b07cfSJed Brown 63b79b07cfSJed Brown #undef __FUNCT__ 64b79b07cfSJed Brown #define __FUNCT__ "SNESSetUp_KSPONLY" 65b79b07cfSJed Brown static PetscErrorCode SNESSetUp_KSPONLY(SNES snes) 66b79b07cfSJed Brown { 676cab3a1bSJed Brown PetscErrorCode ierr; 68b79b07cfSJed Brown 69b79b07cfSJed Brown PetscFunctionBegin; 706cab3a1bSJed Brown ierr = SNESSetUpMatrices(snes);CHKERRQ(ierr); 71b79b07cfSJed Brown PetscFunctionReturn(0); 72b79b07cfSJed Brown } 73b79b07cfSJed Brown 74b79b07cfSJed Brown #undef __FUNCT__ 75b79b07cfSJed Brown #define __FUNCT__ "SNESDestroy_KSPONLY" 76b79b07cfSJed Brown static PetscErrorCode SNESDestroy_KSPONLY(SNES snes) 77b79b07cfSJed Brown { 78b79b07cfSJed Brown 79b79b07cfSJed Brown PetscFunctionBegin; 80b79b07cfSJed Brown PetscFunctionReturn(0); 81b79b07cfSJed Brown } 82b79b07cfSJed Brown 83b79b07cfSJed Brown /* -------------------------------------------------------------------------- */ 84b79b07cfSJed Brown /*MC 85b79b07cfSJed Brown SNESKSPONLY - Nonlinear solver that only performs one Newton step and does not compute any norms. 86b79b07cfSJed Brown The main purpose of this solver is to solve linear problems using the SNES interface, without 87b79b07cfSJed Brown any additional overhead in the form of vector operations. 88b79b07cfSJed Brown 89b79b07cfSJed Brown Level: beginner 90b79b07cfSJed Brown 9104d7464bSBarry Smith .seealso: SNESCreate(), SNES, SNESSetType(), SNESNEWTONLS, SNESNEWTONTR 92b79b07cfSJed Brown M*/ 93b79b07cfSJed Brown #undef __FUNCT__ 94b79b07cfSJed Brown #define __FUNCT__ "SNESCreate_KSPONLY" 958cc058d9SJed Brown PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes) 96b79b07cfSJed Brown { 97b79b07cfSJed Brown 98b79b07cfSJed Brown PetscFunctionBegin; 99b79b07cfSJed Brown snes->ops->setup = SNESSetUp_KSPONLY; 100b79b07cfSJed Brown snes->ops->solve = SNESSolve_KSPONLY; 101b79b07cfSJed Brown snes->ops->destroy = SNESDestroy_KSPONLY; 10237596af1SLisandro Dalcin snes->ops->setfromoptions = 0; 10337596af1SLisandro Dalcin snes->ops->view = 0; 10437596af1SLisandro Dalcin snes->ops->reset = 0; 105b79b07cfSJed Brown 10642f4f86dSBarry Smith snes->usesksp = PETSC_TRUE; 10742f4f86dSBarry Smith snes->usespc = PETSC_FALSE; 10842f4f86dSBarry Smith 109b79b07cfSJed Brown snes->data = 0; 110b79b07cfSJed Brown PetscFunctionReturn(0); 111b79b07cfSJed Brown } 112