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