1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    MFN routines related to the solution process
 12: */

 14: #include <slepc/private/mfnimpl.h>

 16: static PetscErrorCode MFNSolve_Private(MFN mfn,Vec b,Vec x)
 17: {
 18:   PetscFunctionBegin;
 19:   PetscCall(VecSetErrorIfLocked(x,3));

 21:   /* call setup */
 22:   PetscCall(MFNSetUp(mfn));
 23:   mfn->its = 0;

 25:   PetscCall(MFNViewFromOptions(mfn,NULL,"-mfn_view_pre"));

 27:   /* check nonzero right-hand side */
 28:   PetscCall(VecNorm(b,NORM_2,&mfn->bnorm));
 29:   PetscCheck(mfn->bnorm,PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_WRONG,"Cannot pass a zero b vector to MFNSolve()");

 31:   /* call solver */
 32:   PetscCall(PetscLogEventBegin(MFN_Solve,mfn,b,x,0));
 33:   if (b!=x) PetscCall(VecLockReadPush(b));
 34:   PetscUseTypeMethod(mfn,solve,b,x);
 35:   if (b!=x) PetscCall(VecLockReadPop(b));
 36:   PetscCall(PetscLogEventEnd(MFN_Solve,mfn,b,x,0));

 38:   PetscCheck(mfn->reason,PetscObjectComm((PetscObject)mfn),PETSC_ERR_PLIB,"Internal error, solver returned without setting converged reason");

 40:   PetscCheck(!mfn->errorifnotconverged || mfn->reason>=0,PetscObjectComm((PetscObject)mfn),PETSC_ERR_NOT_CONVERGED,"MFNSolve has not converged");

 42:   /* various viewers */
 43:   PetscCall(MFNViewFromOptions(mfn,NULL,"-mfn_view"));
 44:   PetscCall(MFNConvergedReasonViewFromOptions(mfn));
 45:   PetscCall(MatViewFromOptions(mfn->A,(PetscObject)mfn,"-mfn_view_mat"));
 46:   PetscCall(VecViewFromOptions(b,(PetscObject)mfn,"-mfn_view_rhs"));
 47:   PetscCall(VecViewFromOptions(x,(PetscObject)mfn,"-mfn_view_solution"));
 48:   PetscFunctionReturn(PETSC_SUCCESS);
 49: }

 51: /*@
 52:    MFNSolve - Solves the matrix function problem. Given a vector b, the
 53:    vector x = f(A)*b is returned.

 55:    Collective

 57:    Input Parameters:
 58: +  mfn - matrix function context obtained from MFNCreate()
 59: -  b   - the right hand side vector

 61:    Output Parameter:
 62: .  x   - the solution (this may be the same vector as b, then b will be
 63:          overwritten with the answer)

 65:    Options Database Keys:
 66: +  -mfn_view - print information about the solver used
 67: .  -mfn_view_mat binary - save the matrix to the default binary viewer
 68: .  -mfn_view_rhs binary - save right hand side vector to the default binary viewer
 69: .  -mfn_view_solution binary - save computed solution vector to the default binary viewer
 70: -  -mfn_converged_reason - print reason for convergence, and number of iterations

 72:    Notes:
 73:    The matrix A is specified with MFNSetOperator().
 74:    The function f is specified with MFNSetFN().

 76:    Level: beginner

 78: .seealso: MFNCreate(), MFNSetUp(), MFNDestroy(), MFNSetTolerances(),
 79:           MFNSetOperator(), MFNSetFN()
 80: @*/
 81: PetscErrorCode MFNSolve(MFN mfn,Vec b,Vec x)
 82: {
 83:   PetscFunctionBegin;
 86:   PetscCheckSameComm(mfn,1,b,2);
 88:   if (b!=x) PetscCheckSameComm(mfn,1,x,3);
 89:   mfn->transpose_solve = PETSC_FALSE;
 90:   PetscCall(MFNSolve_Private(mfn,b,x));
 91:   PetscFunctionReturn(PETSC_SUCCESS);
 92: }

 94: /*@
 95:    MFNSolveTranspose - Solves the transpose matrix function problem. Given a vector b,
 96:    the vector x = f(A^T)*b is returned.

 98:    Collective

100:    Input Parameters:
101: +  mfn - matrix function context obtained from MFNCreate()
102: -  b   - the right hand side vector

104:    Output Parameter:
105: .  x   - the solution (this may be the same vector as b, then b will be
106:          overwritten with the answer)

108:    Note:
109:    See available options at MFNSolve().

111:    Level: beginner

113: .seealso: MFNSolve()
114: @*/
115: PetscErrorCode MFNSolveTranspose(MFN mfn,Vec b,Vec x)
116: {
117:   PetscFunctionBegin;
120:   PetscCheckSameComm(mfn,1,b,2);
122:   if (b!=x) PetscCheckSameComm(mfn,1,x,3);
123:   mfn->transpose_solve = PETSC_TRUE;
124:   if (!mfn->AT) PetscCall(MatCreateTranspose(mfn->A,&mfn->AT));
125:   PetscCall(MFNSolve_Private(mfn,b,x));
126:   PetscFunctionReturn(PETSC_SUCCESS);
127: }

129: /*@
130:    MFNGetIterationNumber - Gets the current iteration number. If the
131:    call to MFNSolve() is complete, then it returns the number of iterations
132:    carried out by the solution method.

134:    Not Collective

136:    Input Parameter:
137: .  mfn - the matrix function context

139:    Output Parameter:
140: .  its - number of iterations

142:    Note:
143:    During the i-th iteration this call returns i-1. If MFNSolve() is
144:    complete, then parameter "its" contains either the iteration number at
145:    which convergence was successfully reached, or failure was detected.
146:    Call MFNGetConvergedReason() to determine if the solver converged or
147:    failed and why.

149:    Level: intermediate

151: .seealso: MFNGetConvergedReason(), MFNSetTolerances()
152: @*/
153: PetscErrorCode MFNGetIterationNumber(MFN mfn,PetscInt *its)
154: {
155:   PetscFunctionBegin;
157:   PetscAssertPointer(its,2);
158:   *its = mfn->its;
159:   PetscFunctionReturn(PETSC_SUCCESS);
160: }

162: /*@
163:    MFNGetConvergedReason - Gets the reason why the MFNSolve() iteration was
164:    stopped.

166:    Not Collective

168:    Input Parameter:
169: .  mfn - the matrix function context

171:    Output Parameter:
172: .  reason - negative value indicates diverged, positive value converged

174:    Notes:

176:    Possible values for reason are
177: +  MFN_CONVERGED_TOL - converged up to tolerance
178: .  MFN_CONVERGED_ITS - solver completed the requested number of steps
179: .  MFN_DIVERGED_ITS - required more than max_it iterations to reach convergence
180: -  MFN_DIVERGED_BREAKDOWN - generic breakdown in method

182:    Can only be called after the call to MFNSolve() is complete.

184:    Basic solvers (e.g. unrestarted Krylov iterations) cannot determine if the
185:    computation is accurate up to the requested tolerance. In that case, the
186:    converged reason is set to MFN_CONVERGED_ITS if the requested number of steps
187:    (for instance, the ncv value in unrestarted Krylov methods) have been
188:    completed successfully.

190:    Level: intermediate

192: .seealso: MFNSetTolerances(), MFNSolve(), MFNConvergedReason, MFNSetErrorIfNotConverged()
193: @*/
194: PetscErrorCode MFNGetConvergedReason(MFN mfn,MFNConvergedReason *reason)
195: {
196:   PetscFunctionBegin;
198:   PetscAssertPointer(reason,2);
199:   *reason = mfn->reason;
200:   PetscFunctionReturn(PETSC_SUCCESS);
201: }