1 #include <petsc/private/taoimpl.h> /*I "petsctao.h" I*/ 2 3 typedef struct { 4 SNES snes; 5 } Tao_SNES; 6 7 static PetscErrorCode TaoSolve_SNES(Tao tao) 8 { 9 Tao_SNES *taosnes = (Tao_SNES *)tao->data; 10 PetscInt its; 11 12 PetscFunctionBegin; 13 /* TODO SNES fails if KSP reaches max_it, while TAO accepts whatever we got */ 14 PetscCall(SNESSolve(taosnes->snes, NULL, tao->solution)); 15 /* TODO REASONS */ 16 tao->reason = TAO_CONVERGED_USER; 17 PetscCall(SNESGetIterationNumber(taosnes->snes, &its)); 18 PetscCall(TaoSetIterationNumber(tao, its)); 19 PetscFunctionReturn(PETSC_SUCCESS); 20 } 21 22 static PetscErrorCode TaoDestroy_SNES(Tao tao) 23 { 24 Tao_SNES *taosnes = (Tao_SNES *)tao->data; 25 26 PetscFunctionBegin; 27 PetscCall(SNESDestroy(&taosnes->snes)); 28 PetscCall(PetscFree(tao->data)); 29 PetscFunctionReturn(PETSC_SUCCESS); 30 } 31 32 static PetscErrorCode TAOSNESObj(SNES snes, Vec X, PetscReal *f, void *ctx) 33 { 34 Tao tao = (Tao)ctx; 35 36 PetscFunctionBegin; 37 PetscCall(TaoComputeObjective(tao, X, f)); 38 PetscFunctionReturn(PETSC_SUCCESS); 39 } 40 41 static PetscErrorCode TAOSNESFunc(SNES snes, Vec X, Vec F, void *ctx) 42 { 43 Tao tao = (Tao)ctx; 44 45 PetscFunctionBegin; 46 PetscCall(TaoComputeGradient(tao, X, F)); 47 PetscFunctionReturn(PETSC_SUCCESS); 48 } 49 50 static PetscErrorCode TAOSNESJac(SNES snes, Vec X, Mat A, Mat P, void *ctx) 51 { 52 Tao tao = (Tao)ctx; 53 54 PetscFunctionBegin; 55 PetscCall(TaoComputeHessian(tao, X, A, P)); 56 PetscFunctionReturn(PETSC_SUCCESS); 57 } 58 59 static PetscErrorCode TAOSNESMonitor(SNES snes, PetscInt its, PetscReal fnorm, void *ctx) 60 { 61 Tao tao = (Tao)ctx; 62 PetscReal obj; 63 Vec X; 64 65 PetscFunctionBegin; 66 PetscCall(SNESGetSolution(snes, &X)); 67 PetscCall(TaoComputeObjective(tao, X, &obj)); 68 PetscCall(TaoSetIterationNumber(tao, its)); 69 PetscCall(TaoMonitor(tao, its, obj, fnorm, 0.0, 0.0)); 70 PetscFunctionReturn(PETSC_SUCCESS); 71 } 72 73 static PetscErrorCode TaoSetUp_SNES(Tao tao) 74 { 75 Tao_SNES *taosnes = (Tao_SNES *)tao->data; 76 Mat A, P; 77 78 PetscFunctionBegin; 79 PetscCall(SNESSetSolution(taosnes->snes, tao->solution)); 80 PetscCall(SNESSetObjective(taosnes->snes, TAOSNESObj, tao)); 81 PetscCall(SNESSetFunction(taosnes->snes, NULL, TAOSNESFunc, tao)); 82 PetscCall(SNESMonitorSet(taosnes->snes, TAOSNESMonitor, tao, NULL)); 83 PetscCall(TaoGetHessian(tao, &A, &P, NULL, NULL)); 84 if (A) PetscCall(SNESSetJacobian(taosnes->snes, A, P, TAOSNESJac, tao)); 85 /* TODO TYPES */ 86 PetscCall(SNESSetUp(taosnes->snes)); 87 PetscFunctionReturn(PETSC_SUCCESS); 88 } 89 90 static PetscErrorCode TaoSetFromOptions_SNES(Tao tao, PetscOptionItems *PetscOptionsObject) 91 { 92 Tao_SNES *taosnes = (Tao_SNES *)tao->data; 93 94 PetscFunctionBegin; 95 PetscCall(SNESSetFromOptions(taosnes->snes)); 96 PetscFunctionReturn(PETSC_SUCCESS); 97 } 98 99 static PetscErrorCode TaoView_SNES(Tao tao, PetscViewer viewer) 100 { 101 Tao_SNES *taosnes = (Tao_SNES *)tao->data; 102 103 PetscFunctionBegin; 104 PetscCall(SNESView(taosnes->snes, viewer)); 105 PetscFunctionReturn(PETSC_SUCCESS); 106 } 107 108 /*MC 109 TAOSNES - nonlinear solver using SNES 110 111 Level: advanced 112 113 .seealso: `TaoCreate()`, `Tao`, `TaoSetType()`, `TaoType` 114 M*/ 115 PETSC_EXTERN PetscErrorCode TaoCreate_SNES(Tao tao) 116 { 117 Tao_SNES *taosnes; 118 119 PetscFunctionBegin; 120 tao->ops->destroy = TaoDestroy_SNES; 121 tao->ops->setup = TaoSetUp_SNES; 122 tao->ops->setfromoptions = TaoSetFromOptions_SNES; 123 tao->ops->view = TaoView_SNES; 124 tao->ops->solve = TaoSolve_SNES; 125 126 PetscCall(PetscNew(&taosnes)); 127 tao->data = (void *)taosnes; 128 PetscCall(SNESCreate(PetscObjectComm((PetscObject)tao), &taosnes->snes)); 129 PetscCall(PetscObjectIncrementTabLevel((PetscObject)taosnes->snes, (PetscObject)tao, 1)); 130 PetscFunctionReturn(PETSC_SUCCESS); 131 } 132