Actual source code: test4.c

      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: */

     11: static char help[] = "Test the solution of a HEP without calling EPSSetFromOptions (based on ex1.c).\n\n"
     12:   "The command line options are:\n"
     13:   "  -n <n>, where <n> = number of grid subdivisions = matrix dimension.\n"
     14:   "  -type <eps_type> = eps type to test.\n\n";

     16: #include <slepceps.h>

     18: int main(int argc,char **argv)
     19: {
     20:   Mat            A;           /* problem matrix */
     21:   EPS            eps;         /* eigenproblem solver context */
     22:   PetscReal      tol=1000*PETSC_MACHINE_EPSILON;
     23:   PetscInt       n=30,i,Istart,Iend,nev;
     24:   PetscBool      flg,gd2;
     25:   char           epstype[30] = "krylovschur";

     27:   PetscFunctionBeginUser;
     28:   PetscCall(SlepcInitialize(&argc,&argv,NULL,help));

     30:   PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL));
     31:   PetscCall(PetscOptionsGetString(NULL,NULL,"-type",epstype,sizeof(epstype),NULL));
     32:   PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\n1-D Laplacian Eigenproblem, n=%" PetscInt_FMT "\n\n",n));

     34:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     35:      Compute the operator matrix that defines the eigensystem, Ax=kx
     36:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

     38:   PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
     39:   PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n));
     40:   PetscCall(MatSetFromOptions(A));

     42:   PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
     43:   for (i=Istart;i<Iend;i++) {
     44:     if (i>0) PetscCall(MatSetValue(A,i,i-1,-1.0,INSERT_VALUES));
     45:     if (i<n-1) PetscCall(MatSetValue(A,i,i+1,-1.0,INSERT_VALUES));
     46:     PetscCall(MatSetValue(A,i,i,2.0,INSERT_VALUES));
     47:   }
     48:   PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
     49:   PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));

     51:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     52:                 Create the eigensolver and set various options
     53:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
     54:   /*
     55:      Create eigensolver context
     56:   */
     57:   PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));

     59:   /*
     60:      Set operators. In this case, it is a standard eigenvalue problem
     61:   */
     62:   PetscCall(EPSSetOperators(eps,A,NULL));
     63:   PetscCall(EPSSetProblemType(eps,EPS_HEP));
     64:   PetscCall(EPSSetDimensions(eps,4,PETSC_DETERMINE,PETSC_DETERMINE));
     65:   PetscCall(EPSSetTolerances(eps,tol,PETSC_CURRENT));

     67:   /*
     68:      Set solver parameters at runtime
     69:   */
     70:   PetscCall(PetscStrcmp(epstype,"gd2",&flg));
     71:   if (flg) {
     72:     PetscCall(EPSSetType(eps,EPSGD));
     73:     PetscCall(EPSGDSetDoubleExpansion(eps,PETSC_TRUE));
     74:     PetscCall(EPSGDGetDoubleExpansion(eps,&gd2));  /* not used */
     75:   } else PetscCall(EPSSetType(eps,epstype));
     76:   PetscCall(PetscStrcmp(epstype,"jd",&flg));
     77:   if (flg) {
     78:     PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE));
     79:     PetscCall(EPSSetTarget(eps,4.0));
     80:   }
     81:   PetscCall(PetscStrcmp(epstype,"lanczos",&flg));
     82:   if (flg) PetscCall(EPSLanczosSetReorthog(eps,EPS_LANCZOS_REORTHOG_LOCAL));
     83:   PetscCall(PetscObjectTypeCompareAny((PetscObject)eps,&flg,EPSRQCG,EPSLOBPCG,""));
     84:   if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL));

     86:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     87:                       Solve the eigensystem
     88:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

     90:   PetscCall(EPSSolve(eps));
     91:   PetscCall(EPSGetDimensions(eps,&nev,NULL,NULL));
     92:   PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of requested eigenvalues: %" PetscInt_FMT "\n",nev));

     94:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     95:                     Display solution and clean up
     96:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

     98:   PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL));
     99:   PetscCall(EPSDestroy(&eps));
    100:   PetscCall(MatDestroy(&A));
    101:   PetscCall(SlepcFinalize());
    102:   return 0;
    103: }

    105: /*TEST

    107:    testset:
    108:       filter: sed -e "s/3.95905/3.95906/"
    109:       output_file: output/test4_1.out
    110:       test:
    111:          suffix: 1
    112:          args: -type {{krylovschur subspace arnoldi lanczos gd jd gd2 lapack}}
    113:       test:
    114:          suffix: 1_arpack
    115:          args: -type arpack
    116:          requires: arpack
    117:       test:
    118:          suffix: 1_primme
    119:          args: -type primme
    120:          requires: primme

    122:    test:
    123:       suffix: 2
    124:       args: -type {{rqcg lobpcg}}

    126: TEST*/