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 SNES_KSPONLY *ksponly = (SNES_KSPONLY *)snes->data; 9 PetscInt lits; 10 Vec Y, X, F; 11 12 PetscFunctionBegin; 13 PetscCheck(!snes->xl && !snes->xu && !snes->ops->computevariablebounds, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name); 14 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 if (!snes->vec_func_init_set) { 26 PetscCall(SNESComputeFunction(snes, X, F)); 27 } else snes->vec_func_init_set = PETSC_FALSE; 28 29 if (snes->numbermonitors) { 30 PetscReal fnorm; 31 PetscCall(VecNorm(F, NORM_2, &fnorm)); 32 PetscCall(SNESMonitor(snes, 0, fnorm)); 33 } 34 35 /* Call general purpose update function */ 36 PetscTryTypeMethod(snes, update, 0); 37 38 /* Solve J Y = F, where J is Jacobian matrix */ 39 PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre)); 40 41 SNESCheckJacobianDomainerror(snes); 42 43 PetscCall(KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre)); 44 if (ksponly->transpose_solve) { 45 PetscCall(KSPSolveTranspose(snes->ksp, F, Y)); 46 } else { 47 PetscCall(KSPSolve(snes->ksp, F, Y)); 48 } 49 snes->reason = SNES_CONVERGED_ITS; 50 SNESCheckKSPSolve(snes); 51 52 PetscCall(KSPGetIterationNumber(snes->ksp, &lits)); 53 PetscCall(PetscInfo(snes, "iter=%" PetscInt_FMT ", linear solve iterations=%" PetscInt_FMT "\n", snes->iter, lits)); 54 snes->iter++; 55 56 /* Take the computed step. */ 57 PetscCall(VecAXPY(X, -1.0, Y)); 58 if (snes->numbermonitors) { 59 PetscReal fnorm; 60 PetscCall(SNESComputeFunction(snes, X, F)); 61 PetscCall(VecNorm(F, NORM_2, &fnorm)); 62 PetscCall(SNESMonitor(snes, 1, fnorm)); 63 } 64 PetscFunctionReturn(0); 65 } 66 67 static PetscErrorCode SNESSetUp_KSPONLY(SNES snes) { 68 PetscFunctionBegin; 69 PetscCall(SNESSetUpMatrices(snes)); 70 PetscFunctionReturn(0); 71 } 72 73 static PetscErrorCode SNESDestroy_KSPONLY(SNES snes) { 74 PetscFunctionBegin; 75 PetscCall(PetscFree(snes->data)); 76 PetscFunctionReturn(0); 77 } 78 79 /* -------------------------------------------------------------------------- */ 80 /*MC 81 SNESKSPONLY - Nonlinear solver that only performs one Newton step and does not compute any norms. 82 The main purpose of this solver is to solve linear problems using the SNES interface, without 83 any additional overhead in the form of vector operations. 84 85 Level: beginner 86 87 .seealso: `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESNEWTONTR` 88 M*/ 89 PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes) { 90 SNES_KSPONLY *ksponly; 91 92 PetscFunctionBegin; 93 snes->ops->setup = SNESSetUp_KSPONLY; 94 snes->ops->solve = SNESSolve_KSPONLY; 95 snes->ops->destroy = SNESDestroy_KSPONLY; 96 snes->ops->setfromoptions = NULL; 97 snes->ops->view = NULL; 98 snes->ops->reset = NULL; 99 100 snes->usesksp = PETSC_TRUE; 101 snes->usesnpc = PETSC_FALSE; 102 103 snes->alwayscomputesfinalresidual = PETSC_FALSE; 104 105 PetscCall(PetscNewLog(snes, &ksponly)); 106 snes->data = (void *)ksponly; 107 PetscFunctionReturn(0); 108 } 109 110 /*MC 111 SNESKSPTRANSPOSEONLY - Nonlinear solver that only performs one Newton step and does not compute any norms. 112 The main purpose of this solver is to solve transposed linear problems using the SNES interface, without 113 any additional overhead in the form of vector operations within adjoint solvers. 114 115 Level: beginner 116 117 .seealso: `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESKSPTRANSPOSEONLY`, `SNESNEWTONLS`, `SNESNEWTONTR` 118 M*/ 119 PETSC_EXTERN PetscErrorCode SNESCreate_KSPTRANSPOSEONLY(SNES snes) { 120 SNES_KSPONLY *kspo; 121 122 PetscFunctionBegin; 123 PetscCall(SNESCreate_KSPONLY(snes)); 124 PetscCall(PetscObjectChangeTypeName((PetscObject)snes, SNESKSPTRANSPOSEONLY)); 125 kspo = (SNES_KSPONLY *)snes->data; 126 kspo->transpose_solve = PETSC_TRUE; 127 PetscFunctionReturn(0); 128 } 129