11ef27442SStefano Zampini #include <petsc/private/snesimpl.h> /*I "petscsnes.h" I*/ 2b79b07cfSJed Brown 31ef27442SStefano Zampini typedef struct { 41ef27442SStefano Zampini PetscBool transpose_solve; 51ef27442SStefano Zampini } SNES_KSPONLY; 6b79b07cfSJed Brown 7d71ae5a4SJacob Faibussowitsch static PetscErrorCode SNESSolve_KSPONLY(SNES snes) 8d71ae5a4SJacob Faibussowitsch { 91ef27442SStefano Zampini SNES_KSPONLY *ksponly = (SNES_KSPONLY *)snes->data; 10b79b07cfSJed Brown PetscInt lits; 11b79b07cfSJed Brown Vec Y, X, F; 12b79b07cfSJed Brown 13b79b07cfSJed Brown PetscFunctionBegin; 140b121fc5SBarry Smith 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); 15c579b300SPatrick Farrell 16b79b07cfSJed Brown snes->numFailures = 0; 17b79b07cfSJed Brown snes->numLinearSolveFailures = 0; 18b79b07cfSJed Brown snes->reason = SNES_CONVERGED_ITERATING; 19b79b07cfSJed Brown snes->iter = 0; 20b79b07cfSJed Brown snes->norm = 0.0; 21b79b07cfSJed Brown 22b79b07cfSJed Brown X = snes->vec_sol; 23b79b07cfSJed Brown F = snes->vec_func; 24b79b07cfSJed Brown Y = snes->vec_sol_update; 25b79b07cfSJed Brown 26451891c6SStefano Zampini if (!snes->vec_func_init_set) { 279566063dSJacob Faibussowitsch PetscCall(SNESComputeFunction(snes, X, F)); 28451891c6SStefano Zampini } else snes->vec_func_init_set = PETSC_FALSE; 29451891c6SStefano Zampini 30a5eaddd9SBarry Smith if (snes->numbermonitors) { 31a5eaddd9SBarry Smith PetscReal fnorm; 329566063dSJacob Faibussowitsch PetscCall(VecNorm(F, NORM_2, &fnorm)); 335f3c5e7aSBarry Smith SNESCheckFunctionNorm(snes, fnorm); 349566063dSJacob Faibussowitsch PetscCall(SNESMonitor(snes, 0, fnorm)); 35a5eaddd9SBarry Smith } 36b79b07cfSJed Brown 375341784dSBarry Smith /* Call general purpose update function */ 38dbbe0bcdSBarry Smith PetscTryTypeMethod(snes, update, 0); 395341784dSBarry Smith 40b79b07cfSJed Brown /* Solve J Y = F, where J is Jacobian matrix */ 419566063dSJacob Faibussowitsch PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre)); 4207b62357SFande Kong 4307b62357SFande Kong SNESCheckJacobianDomainerror(snes); 4407b62357SFande Kong 459566063dSJacob Faibussowitsch PetscCall(KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre)); 461ef27442SStefano Zampini if (ksponly->transpose_solve) { 479566063dSJacob Faibussowitsch PetscCall(KSPSolveTranspose(snes->ksp, F, Y)); 481ef27442SStefano Zampini } else { 499566063dSJacob Faibussowitsch PetscCall(KSPSolve(snes->ksp, F, Y)); 501ef27442SStefano Zampini } 51422a814eSBarry Smith snes->reason = SNES_CONVERGED_ITS; 52422a814eSBarry Smith SNESCheckKSPSolve(snes); 531aa26658SKarl Rupp 549566063dSJacob Faibussowitsch PetscCall(KSPGetIterationNumber(snes->ksp, &lits)); 5563a3b9bcSJacob Faibussowitsch PetscCall(PetscInfo(snes, "iter=%" PetscInt_FMT ", linear solve iterations=%" PetscInt_FMT "\n", snes->iter, lits)); 56b79b07cfSJed Brown snes->iter++; 57b79b07cfSJed Brown 58b79b07cfSJed Brown /* Take the computed step. */ 599566063dSJacob Faibussowitsch PetscCall(VecAXPY(X, -1.0, Y)); 60a5eaddd9SBarry Smith if (snes->numbermonitors) { 61a5eaddd9SBarry Smith PetscReal fnorm; 629566063dSJacob Faibussowitsch PetscCall(SNESComputeFunction(snes, X, F)); 639566063dSJacob Faibussowitsch PetscCall(VecNorm(F, NORM_2, &fnorm)); 645f3c5e7aSBarry Smith SNESCheckFunctionNorm(snes, fnorm); 659566063dSJacob Faibussowitsch PetscCall(SNESMonitor(snes, 1, fnorm)); 66a5eaddd9SBarry Smith } 673ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 68b79b07cfSJed Brown } 69b79b07cfSJed Brown 70d71ae5a4SJacob Faibussowitsch static PetscErrorCode SNESSetUp_KSPONLY(SNES snes) 71d71ae5a4SJacob Faibussowitsch { 72b79b07cfSJed Brown PetscFunctionBegin; 739566063dSJacob Faibussowitsch PetscCall(SNESSetUpMatrices(snes)); 743ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 75b79b07cfSJed Brown } 76b79b07cfSJed Brown 77d71ae5a4SJacob Faibussowitsch static PetscErrorCode SNESDestroy_KSPONLY(SNES snes) 78d71ae5a4SJacob Faibussowitsch { 79b79b07cfSJed Brown PetscFunctionBegin; 809566063dSJacob Faibussowitsch PetscCall(PetscFree(snes->data)); 813ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 82b79b07cfSJed Brown } 83b79b07cfSJed Brown 84b79b07cfSJed Brown /*MC 85420bcc1bSBarry Smith SNESKSPONLY - Nonlinear solver that performs one Newton step with `KSPSolve()` and does not compute any norms. 86b79b07cfSJed Brown 87b79b07cfSJed Brown Level: beginner 88b79b07cfSJed Brown 89420bcc1bSBarry Smith Note: 90420bcc1bSBarry Smith The main purpose of this solver is to solve linear problems using the `SNES` interface, without 91420bcc1bSBarry Smith any additional overhead in the form of vector norm operations. 92420bcc1bSBarry Smith 93420bcc1bSBarry Smith .seealso: [](ch_snes), `SNES`, `SNESType`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESNEWTONTR`, `SNESKSPTRANSPOSEONLY` 94b79b07cfSJed Brown M*/ 95d71ae5a4SJacob Faibussowitsch PETSC_EXTERN PetscErrorCode SNESCreate_KSPONLY(SNES snes) 96d71ae5a4SJacob Faibussowitsch { 971ef27442SStefano Zampini SNES_KSPONLY *ksponly; 98b79b07cfSJed Brown 99b79b07cfSJed Brown PetscFunctionBegin; 100b79b07cfSJed Brown snes->ops->setup = SNESSetUp_KSPONLY; 101b79b07cfSJed Brown snes->ops->solve = SNESSolve_KSPONLY; 102b79b07cfSJed Brown snes->ops->destroy = SNESDestroy_KSPONLY; 1039e5d0892SLisandro Dalcin snes->ops->setfromoptions = NULL; 1049e5d0892SLisandro Dalcin snes->ops->view = NULL; 1059e5d0892SLisandro Dalcin snes->ops->reset = NULL; 106b79b07cfSJed Brown 10742f4f86dSBarry Smith snes->usesksp = PETSC_TRUE; 108efd4aadfSBarry Smith snes->usesnpc = PETSC_FALSE; 10942f4f86dSBarry Smith 1104fc747eaSLawrence Mitchell snes->alwayscomputesfinalresidual = PETSC_FALSE; 1114fc747eaSLawrence Mitchell 112*77e5a1f9SBarry Smith PetscCall(SNESParametersInitialize(snes)); 113*77e5a1f9SBarry Smith 1144dfa11a4SJacob Faibussowitsch PetscCall(PetscNew(&ksponly)); 1151ef27442SStefano Zampini snes->data = (void *)ksponly; 1163ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 1171ef27442SStefano Zampini } 1181ef27442SStefano Zampini 1191ef27442SStefano Zampini /*MC 120420bcc1bSBarry Smith SNESKSPTRANSPOSEONLY - Nonlinear solver that performs one Newton step with `KSPSolveTranspose()` and does not compute any norms. 1211ef27442SStefano Zampini 1221ef27442SStefano Zampini Level: beginner 1231ef27442SStefano Zampini 124420bcc1bSBarry Smith Note: 125420bcc1bSBarry Smith The main purpose of this solver is to solve transposed linear problems using the `SNES` interface, without 126420bcc1bSBarry Smith any additional overhead in the form of vector operations within adjoint solvers. 127420bcc1bSBarry Smith 128420bcc1bSBarry Smith .seealso: [](ch_snes), `SNES`, `SNESType`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESKS`, `SNESNEWTONLS`, `SNESNEWTONTR` 1291ef27442SStefano Zampini M*/ 130d71ae5a4SJacob Faibussowitsch PETSC_EXTERN PetscErrorCode SNESCreate_KSPTRANSPOSEONLY(SNES snes) 131d71ae5a4SJacob Faibussowitsch { 1321ef27442SStefano Zampini SNES_KSPONLY *kspo; 1331ef27442SStefano Zampini 1341ef27442SStefano Zampini PetscFunctionBegin; 1359566063dSJacob Faibussowitsch PetscCall(SNESCreate_KSPONLY(snes)); 1369566063dSJacob Faibussowitsch PetscCall(PetscObjectChangeTypeName((PetscObject)snes, SNESKSPTRANSPOSEONLY)); 1371ef27442SStefano Zampini kspo = (SNES_KSPONLY *)snes->data; 1381ef27442SStefano Zampini kspo->transpose_solve = PETSC_TRUE; 1393ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 140b79b07cfSJed Brown } 141