Line data    Source code

       1             : /*
       2             :    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
       3             :    SLEPc - Scalable Library for Eigenvalue Problem Computations
       4             :    Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
       5             : 
       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             : */
      13             : 
      14             : #include <slepc/private/mfnimpl.h>       /*I "slepcmfn.h" I*/
      15             : 
      16             : /*@
      17             :    MFNSetUp - Sets up all the internal data structures necessary for the
      18             :    execution of the matrix function solver.
      19             : 
      20             :    Collective
      21             : 
      22             :    Input Parameter:
      23             : .  mfn   - matrix function context
      24             : 
      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.
      29             : 
      30             :    Level: developer
      31             : 
      32             : .seealso: MFNCreate(), MFNSolve(), MFNDestroy()
      33             : @*/
      34         103 : PetscErrorCode MFNSetUp(MFN mfn)
      35             : {
      36         103 :   PetscInt       N;
      37             : 
      38         103 :   PetscFunctionBegin;
      39         103 :   PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
      40             : 
      41             :   /* reset the convergence flag from the previous solves */
      42         103 :   mfn->reason = MFN_CONVERGED_ITERATING;
      43             : 
      44         103 :   if (mfn->setupcalled) PetscFunctionReturn(PETSC_SUCCESS);
      45          17 :   PetscCall(PetscLogEventBegin(MFN_SetUp,mfn,0,0,0));
      46             : 
      47             :   /* Set default solver type (MFNSetFromOptions was not called) */
      48          17 :   if (!((PetscObject)mfn)->type_name) PetscCall(MFNSetType(mfn,MFNKRYLOV));
      49          17 :   if (!mfn->fn) PetscCall(MFNGetFN(mfn,&mfn->fn));
      50          17 :   if (!((PetscObject)mfn->fn)->type_name) PetscCall(FNSetFromOptions(mfn->fn));
      51             : 
      52             :   /* Check problem dimensions */
      53          17 :   PetscCheck(mfn->A,PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_WRONGSTATE,"MFNSetOperator must be called first");
      54          17 :   PetscCall(MatGetSize(mfn->A,&N,NULL));
      55          17 :   if (mfn->ncv > N) mfn->ncv = N;
      56             : 
      57             :   /* call specific solver setup */
      58          17 :   PetscUseTypeMethod(mfn,setup);
      59             : 
      60             :   /* set tolerance if not yet set */
      61          17 :   if (mfn->tol==(PetscReal)PETSC_DETERMINE) mfn->tol = SLEPC_DEFAULT_TOL;
      62             : 
      63          17 :   PetscCall(PetscLogEventEnd(MFN_SetUp,mfn,0,0,0));
      64          17 :   mfn->setupcalled = 1;
      65          17 :   PetscFunctionReturn(PETSC_SUCCESS);
      66             : }
      67             : 
      68             : /*@
      69             :    MFNSetOperator - Sets the matrix for which the matrix function is to be computed.
      70             : 
      71             :    Collective
      72             : 
      73             :    Input Parameters:
      74             : +  mfn - the matrix function context
      75             : -  A   - the problem matrix
      76             : 
      77             :    Notes:
      78             :    It must be called before MFNSetUp(). If it is called again after MFNSetUp() then
      79             :    the MFN object is reset.
      80             : 
      81             :    Level: beginner
      82             : 
      83             : .seealso: MFNSolve(), MFNSetUp(), MFNReset()
      84             : @*/
      85          17 : PetscErrorCode MFNSetOperator(MFN mfn,Mat A)
      86             : {
      87          17 :   PetscInt       m,n;
      88             : 
      89          17 :   PetscFunctionBegin;
      90          17 :   PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
      91          17 :   PetscValidHeaderSpecific(A,MAT_CLASSID,2);
      92          17 :   PetscCheckSameComm(mfn,1,A,2);
      93             : 
      94          17 :   PetscCall(MatGetSize(A,&m,&n));
      95          17 :   PetscCheck(m==n,PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_WRONG,"A is a non-square matrix");
      96          17 :   PetscCall(PetscObjectReference((PetscObject)A));
      97          17 :   if (mfn->setupcalled) PetscCall(MFNReset(mfn));
      98          17 :   else PetscCall(MatDestroy(&mfn->A));
      99          17 :   mfn->A = A;
     100          17 :   mfn->setupcalled = 0;
     101          17 :   PetscFunctionReturn(PETSC_SUCCESS);
     102             : }
     103             : 
     104             : /*@
     105             :    MFNGetOperator - Gets the matrix associated with the MFN object.
     106             : 
     107             :    Not Collective
     108             : 
     109             :    Input Parameter:
     110             : .  mfn - the MFN context
     111             : 
     112             :    Output Parameters:
     113             : .  A  - the matrix for which the matrix function is to be computed
     114             : 
     115             :    Level: intermediate
     116             : 
     117             : .seealso: MFNSolve(), MFNSetOperator()
     118             : @*/
     119           2 : PetscErrorCode MFNGetOperator(MFN mfn,Mat *A)
     120             : {
     121           2 :   PetscFunctionBegin;
     122           2 :   PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
     123           2 :   PetscAssertPointer(A,2);
     124           2 :   *A = mfn->A;
     125           2 :   PetscFunctionReturn(PETSC_SUCCESS);
     126             : }
     127             : 
     128             : /*@
     129             :    MFNAllocateSolution - Allocate memory storage for common variables such
     130             :    as the basis vectors.
     131             : 
     132             :    Collective
     133             : 
     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
     138             : 
     139             :    Developer Notes:
     140             :    This is SLEPC_EXTERN because it may be required by user plugin MFN
     141             :    implementations.
     142             : 
     143             :    Level: developer
     144             : 
     145             : .seealso: MFNSetUp()
     146             : @*/
     147          17 : PetscErrorCode MFNAllocateSolution(MFN mfn,PetscInt extra)
     148             : {
     149          17 :   PetscInt       oldsize,requested;
     150          17 :   Vec            t;
     151             : 
     152          17 :   PetscFunctionBegin;
     153          17 :   requested = mfn->ncv + extra;
     154             : 
     155             :   /* oldsize is zero if this is the first time setup is called */
     156          17 :   PetscCall(BVGetSizes(mfn->V,NULL,NULL,&oldsize));
     157             : 
     158             :   /* allocate basis vectors */
     159          17 :   if (!mfn->V) PetscCall(MFNGetBV(mfn,&mfn->V));
     160          17 :   if (!oldsize) {
     161          16 :     if (!((PetscObject)mfn->V)->type_name) PetscCall(BVSetType(mfn->V,BVMAT));
     162          16 :     PetscCall(MatCreateVecsEmpty(mfn->A,&t,NULL));
     163          16 :     PetscCall(BVSetSizesFromVec(mfn->V,t,requested));
     164          16 :     PetscCall(VecDestroy(&t));
     165           1 :   } else PetscCall(BVResize(mfn->V,requested,PETSC_FALSE));
     166          17 :   PetscFunctionReturn(PETSC_SUCCESS);
     167             : }