1 #include <petsc/private/snesimpl.h> /*I "petscsnes.h" I*/ 2 3 typedef struct { 4 PetscBool transpose_solve; 5 } SNES_KSPONLY; 6 7 static PetscErrorCode SNESSolve_KSPONLY(SNES snes) 8 { 9 SNES_KSPONLY *ksponly = (SNES_KSPONLY*)snes->data; 10 PetscErrorCode ierr; 11 PetscInt lits; 12 Vec Y,X,F; 13 14 PetscFunctionBegin; 15 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); 16 17 snes->numFailures = 0; 18 snes->numLinearSolveFailures = 0; 19 snes->reason = SNES_CONVERGED_ITERATING; 20 snes->iter = 0; 21 snes->norm = 0.0; 22 23 X = snes->vec_sol; 24 F = snes->vec_func; 25 Y = snes->vec_sol_update; 26 27 if (!snes->vec_func_init_set) { 28 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 29 } else snes->vec_func_init_set = PETSC_FALSE; 30 31 if (snes->numbermonitors) { 32 PetscReal fnorm; 33 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 34 ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 35 } 36 37 /* Call general purpose update function */ 38 if (snes->ops->update) { 39 ierr = (*snes->ops->update)(snes, 0);CHKERRQ(ierr); 40 } 41 42 /* Solve J Y = F, where J is Jacobian matrix */ 43 ierr = SNESComputeJacobian(snes,X,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 44 45 SNESCheckJacobianDomainerror(snes); 46 47 ierr = KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 48 if (ksponly->transpose_solve) { 49 ierr = KSPSolveTranspose(snes->ksp,F,Y);CHKERRQ(ierr); 50 } else { 51 ierr = KSPSolve(snes->ksp,F,Y);CHKERRQ(ierr); 52 } 53 snes->reason = SNES_CONVERGED_ITS; 54 SNESCheckKSPSolve(snes); 55 56 ierr = KSPGetIterationNumber(snes->ksp,&lits);CHKERRQ(ierr); 57 ierr = PetscInfo2(snes,"iter=%D, linear solve iterations=%D\n",snes->iter,lits);CHKERRQ(ierr); 58 snes->iter++; 59 60 /* Take the computed step. */ 61 ierr = VecAXPY(X,-1.0,Y);CHKERRQ(ierr); 62 if (snes->numbermonitors) { 63 PetscReal fnorm; 64 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 65 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 66 ierr = SNESMonitor(snes,1,fnorm);CHKERRQ(ierr); 67 } 68 PetscFunctionReturn(0); 69 } 70 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 static PetscErrorCode SNESDestroy_KSPONLY(SNES snes) 81 { 82 PetscErrorCode ierr; 83 84 PetscFunctionBegin; 85 ierr = PetscFree(snes->data);CHKERRQ(ierr); 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 PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes) 100 { 101 SNES_KSPONLY *ksponly; 102 PetscErrorCode ierr; 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 = NULL; 109 snes->ops->view = NULL; 110 snes->ops->reset = NULL; 111 112 snes->usesksp = PETSC_TRUE; 113 snes->usesnpc = PETSC_FALSE; 114 115 snes->alwayscomputesfinalresidual = PETSC_FALSE; 116 117 ierr = PetscNewLog(snes,&ksponly);CHKERRQ(ierr); 118 snes->data = (void*)ksponly; 119 PetscFunctionReturn(0); 120 } 121 122 /*MC 123 SNESKSPTRANSPOSEONLY - Nonlinear solver that only performs one Newton step and does not compute any norms. 124 The main purpose of this solver is to solve transposed linear problems using the SNES interface, without 125 any additional overhead in the form of vector operations within adjoint solvers. 126 127 Level: beginner 128 129 .seealso: SNESCreate(), SNES, SNESSetType(), SNESKSPTRANSPOSEONLY, SNESNEWTONLS, SNESNEWTONTR 130 M*/ 131 PETSC_EXTERN PetscErrorCode SNESCreate_KSPTRANSPOSEONLY(SNES snes) 132 { 133 SNES_KSPONLY *kspo; 134 PetscErrorCode ierr; 135 136 PetscFunctionBegin; 137 ierr = SNESCreate_KSPONLY(snes);CHKERRQ(ierr); 138 ierr = PetscObjectChangeTypeName((PetscObject)snes,SNESKSPTRANSPOSEONLY);CHKERRQ(ierr); 139 kspo = (SNES_KSPONLY*)snes->data; 140 kspo->transpose_solve = PETSC_TRUE; 141 PetscFunctionReturn(0); 142 } 143