Actual source code: mfnsolve.c
slepc-3.22.1 2024-10-28
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: }