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 problem setup
 12: */

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

 16: /*@
 17:    MFNSetUp - Sets up all the internal data structures necessary for the
 18:    execution of the matrix function solver.

 20:    Collective

 22:    Input Parameter:
 23: .  mfn   - matrix function context

 25:    Notes:
 26:    This function need not be called explicitly in most cases, since MFNSolve()
 27:    calls it. It can be useful when one wants to measure the set-up time
 28:    separately from the solve time.

 30:    Level: developer

 32: .seealso: MFNCreate(), MFNSolve(), MFNDestroy()
 33: @*/
 34: PetscErrorCode MFNSetUp(MFN mfn)
 35: {
 36:   PetscInt       N;

 38:   PetscFunctionBegin;

 41:   /* reset the convergence flag from the previous solves */
 42:   mfn->reason = MFN_CONVERGED_ITERATING;

 44:   if (mfn->setupcalled) PetscFunctionReturn(PETSC_SUCCESS);
 45:   PetscCall(PetscLogEventBegin(MFN_SetUp,mfn,0,0,0));

 47:   /* Set default solver type (MFNSetFromOptions was not called) */
 48:   if (!((PetscObject)mfn)->type_name) PetscCall(MFNSetType(mfn,MFNKRYLOV));
 49:   if (!mfn->fn) PetscCall(MFNGetFN(mfn,&mfn->fn));
 50:   if (!((PetscObject)mfn->fn)->type_name) PetscCall(FNSetFromOptions(mfn->fn));

 52:   /* Check problem dimensions */
 53:   PetscCheck(mfn->A,PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_WRONGSTATE,"MFNSetOperator must be called first");
 54:   PetscCall(MatGetSize(mfn->A,&N,NULL));
 55:   if (mfn->ncv > N) mfn->ncv = N;

 57:   /* call specific solver setup */
 58:   PetscUseTypeMethod(mfn,setup);

 60:   /* set tolerance if not yet set */
 61:   if (mfn->tol==(PetscReal)PETSC_DETERMINE) mfn->tol = SLEPC_DEFAULT_TOL;

 63:   PetscCall(PetscLogEventEnd(MFN_SetUp,mfn,0,0,0));
 64:   mfn->setupcalled = 1;
 65:   PetscFunctionReturn(PETSC_SUCCESS);
 66: }

 68: /*@
 69:    MFNSetOperator - Sets the matrix for which the matrix function is to be computed.

 71:    Collective

 73:    Input Parameters:
 74: +  mfn - the matrix function context
 75: -  A   - the problem matrix

 77:    Notes:
 78:    It must be called before MFNSetUp(). If it is called again after MFNSetUp() then
 79:    the MFN object is reset.

 81:    Level: beginner

 83: .seealso: MFNSolve(), MFNSetUp(), MFNReset()
 84: @*/
 85: PetscErrorCode MFNSetOperator(MFN mfn,Mat A)
 86: {
 87:   PetscInt       m,n;

 89:   PetscFunctionBegin;
 92:   PetscCheckSameComm(mfn,1,A,2);

 94:   PetscCall(MatGetSize(A,&m,&n));
 95:   PetscCheck(m==n,PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_WRONG,"A is a non-square matrix");
 96:   PetscCall(PetscObjectReference((PetscObject)A));
 97:   if (mfn->setupcalled) PetscCall(MFNReset(mfn));
 98:   else PetscCall(MatDestroy(&mfn->A));
 99:   mfn->A = A;
100:   mfn->setupcalled = 0;
101:   PetscFunctionReturn(PETSC_SUCCESS);
102: }

104: /*@
105:    MFNGetOperator - Gets the matrix associated with the MFN object.

107:    Not Collective

109:    Input Parameter:
110: .  mfn - the MFN context

112:    Output Parameters:
113: .  A  - the matrix for which the matrix function is to be computed

115:    Level: intermediate

117: .seealso: MFNSolve(), MFNSetOperator()
118: @*/
119: PetscErrorCode MFNGetOperator(MFN mfn,Mat *A)
120: {
121:   PetscFunctionBegin;
123:   PetscAssertPointer(A,2);
124:   *A = mfn->A;
125:   PetscFunctionReturn(PETSC_SUCCESS);
126: }

128: /*@
129:    MFNAllocateSolution - Allocate memory storage for common variables such
130:    as the basis vectors.

132:    Collective

134:    Input Parameters:
135: +  mfn   - matrix function context
136: -  extra - number of additional positions, used for methods that require a
137:            working basis slightly larger than ncv

139:    Developer Notes:
140:    This is SLEPC_EXTERN because it may be required by user plugin MFN
141:    implementations.

143:    Level: developer

145: .seealso: MFNSetUp()
146: @*/
147: PetscErrorCode MFNAllocateSolution(MFN mfn,PetscInt extra)
148: {
149:   PetscInt       oldsize,requested;
150:   Vec            t;

152:   PetscFunctionBegin;
153:   requested = mfn->ncv + extra;

155:   /* oldsize is zero if this is the first time setup is called */
156:   PetscCall(BVGetSizes(mfn->V,NULL,NULL,&oldsize));

158:   /* allocate basis vectors */
159:   if (!mfn->V) PetscCall(MFNGetBV(mfn,&mfn->V));
160:   if (!oldsize) {
161:     if (!((PetscObject)mfn->V)->type_name) PetscCall(BVSetType(mfn->V,BVMAT));
162:     PetscCall(MatCreateVecsEmpty(mfn->A,&t,NULL));
163:     PetscCall(BVSetSizesFromVec(mfn->V,t,requested));
164:     PetscCall(VecDestroy(&t));
165:   } else PetscCall(BVResize(mfn->V,requested,PETSC_FALSE));
166:   PetscFunctionReturn(PETSC_SUCCESS);
167: }