1: static char help[] = "Check TAOSNES convergence reasons on a quadratic problem.\n";

  3: #include <petsctao.h>

  5: static PetscErrorCode FormObjectiveGradient(Tao tao, Vec X, PetscReal *f, Vec G, PetscCtx ctx)
  6: {
  7:   PetscScalar dot;

  9:   PetscFunctionBeginUser;
 10:   PetscCall(VecDot(X, X, &dot));
 11:   *f = 0.5 * PetscRealPart(dot);
 12:   PetscCall(VecCopy(X, G));
 13:   PetscFunctionReturn(PETSC_SUCCESS);
 14: }

 16: static PetscErrorCode FormHessian(Tao tao, Vec X, Mat H, Mat P, PetscCtx ctx)
 17: {
 18:   PetscFunctionBeginUser;
 19:   PetscCall(MatZeroEntries(H));
 20:   PetscCall(MatShift(H, 1.0));
 21:   PetscFunctionReturn(PETSC_SUCCESS);
 22: }

 24: int main(int argc, char **argv)
 25: {
 26:   const SNESConvergedReason snes_reasons[] = {SNES_CONVERGED_FNORM_ABS, SNES_DIVERGED_MAX_IT, SNES_CONVERGED_FNORM_RELATIVE, SNES_DIVERGED_FUNCTION_COUNT, SNES_DIVERGED_LINEAR_SOLVE};
 27:   const TaoConvergedReason  tao_reasons[]  = {TAO_CONVERGED_GATOL, TAO_DIVERGED_MAXITS, TAO_CONVERGED_GTTOL, TAO_DIVERGED_MAXFCN, TAO_DIVERGED_USER};
 28:   Tao                       tao;
 29:   SNES                      snes;
 30:   KSP                       ksp;
 31:   PC                        pc;
 32:   Vec                       X;
 33:   Mat                       H;

 35:   PetscFunctionBeginUser;
 36:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 37:   PetscCall(VecCreateSeq(PETSC_COMM_SELF, 2, &X));
 38:   PetscCall(MatCreateSeqAIJ(PETSC_COMM_SELF, 2, 2, 1, NULL, &H));
 39:   for (PetscInt i = 0; i < 2; ++i) PetscCall(MatSetValue(H, i, i, 1.0, INSERT_VALUES));
 40:   PetscCall(MatAssemblyBegin(H, MAT_FINAL_ASSEMBLY));
 41:   PetscCall(MatAssemblyEnd(H, MAT_FINAL_ASSEMBLY));
 42:   PetscCall(TaoCreate(PETSC_COMM_SELF, &tao));
 43:   PetscCall(TaoSetType(tao, TAOSNES));
 44:   PetscCall(TaoSetSolution(tao, X));
 45:   PetscCall(TaoSetObjectiveAndGradient(tao, NULL, FormObjectiveGradient, NULL));
 46:   PetscCall(TaoSetHessian(tao, H, H, FormHessian, NULL));
 47:   PetscCall(TaoSNESGetSNES(tao, &snes));
 48:   PetscCall(TaoSetFromOptions(tao));
 49:   PetscCall(TaoSetUp(tao));
 50:   PetscCall(SNESGetKSP(snes, &ksp));
 51:   PetscCall(KSPSetType(ksp, KSPCG));
 52:   PetscCall(KSPGetPC(ksp, &pc));
 53:   PetscCall(PCSetType(pc, PCNONE));
 54:   for (size_t test = 0; test < PETSC_STATIC_ARRAY_LENGTH(snes_reasons); ++test) {
 55:     SNESConvergedReason snes_reason;
 56:     TaoConvergedReason  tao_reason;
 57:     PetscInt            snes_its, tao_its;

 59:     PetscCall(VecSet(X, 1.0));
 60:     PetscCall(SNESSetTolerances(snes, test == 2 ? 0.0 : 1.e-12, test == 2 ? 0.5 : 0.0, 0.0, test == 1 ? 0 : 10, test == 3 ? 1 : 100));
 61:     PetscCall(KSPSetTolerances(ksp, 1.e-12, PETSC_CURRENT, PETSC_CURRENT, test == 4 ? 0 : 10));
 62:     PetscCall(TaoSolve(tao));
 63:     PetscCall(SNESGetConvergedReason(snes, &snes_reason));
 64:     PetscCall(TaoGetConvergedReason(tao, &tao_reason));
 65:     PetscCall(SNESGetIterationNumber(snes, &snes_its));
 66:     PetscCall(TaoGetIterationNumber(tao, &tao_its));
 67:     PetscCheck(snes_reason == snes_reasons[test], PETSC_COMM_SELF, PETSC_ERR_PLIB, "Unexpected SNES reason %s", SNESConvergedReasons[snes_reason]);
 68:     PetscCheck(tao_reason == tao_reasons[test], PETSC_COMM_SELF, PETSC_ERR_PLIB, "SNES reason %s mapped to TAO reason %s", SNESConvergedReasons[snes_reason], TaoConvergedReasons[tao_reason]);
 69:     PetscCheck(tao_its == snes_its, PETSC_COMM_SELF, PETSC_ERR_PLIB, "TAO and SNES iteration counts differ");
 70:     PetscCall(PetscPrintf(PETSC_COMM_SELF, "SNES %s -> TAO %s\n", SNESConvergedReasons[snes_reason], TaoConvergedReasons[tao_reason]));
 71:   }
 72:   PetscCall(TaoDestroy(&tao));
 73:   PetscCall(MatDestroy(&H));
 74:   PetscCall(VecDestroy(&X));
 75:   PetscCall(PetscFinalize());
 76:   return 0;
 77: }

 79: /*TEST

 81:   test:
 82:     suffix: 1
 83:     requires: !single !complex
 84:     args: -snes_type {{newtonls newtontr}}

 86: TEST*/