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             : static char help[] = "Test interface to external package PRIMME.\n\n"
      12             :   "This is based on ex3.c. The command line options are:\n"
      13             :   "  -n <n>, where <n> = number of grid subdivisions in x dimension.\n"
      14             :   "  -m <m>, where <m> = number of grid subdivisions in y dimension.\n\n";
      15             : 
      16             : #include <slepceps.h>
      17             : 
      18           4 : int main(int argc,char **argv)
      19             : {
      20           4 :   Mat             A;           /* matrix */
      21           4 :   EPS             eps;         /* eigenproblem solver context */
      22           4 :   ST              st;          /* spectral transformation context */
      23           4 :   KSP             ksp;
      24           4 :   PC              pc;
      25           4 :   PetscInt        N,n=35,m,Istart,Iend,II,i,j,bs;
      26           4 :   PetscBool       flag;
      27           4 :   EPSPRIMMEMethod meth;
      28             : 
      29           4 :   PetscFunctionBeginUser;
      30           4 :   PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
      31             : 
      32           4 :   PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL));
      33           4 :   PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,&flag));
      34           4 :   if (!flag) m=n;
      35           4 :   N = n*m;
      36           4 :   PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nStandard eigenproblem with PRIMME, N=%" PetscInt_FMT " (%" PetscInt_FMT "x%" PetscInt_FMT " grid)\n\n",N,n,m));
      37             : 
      38             :   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      39             :          Compute the matrices that define the eigensystem, Ax=kBx
      40             :      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
      41             : 
      42           4 :   PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
      43           4 :   PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N));
      44           4 :   PetscCall(MatSetFromOptions(A));
      45             : 
      46           4 :   PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
      47        3679 :   for (II=Istart;II<Iend;II++) {
      48        3675 :     i = II/n; j = II-i*n;
      49        3675 :     if (i>0) PetscCall(MatSetValue(A,II,II-n,-1.0,INSERT_VALUES));
      50        3675 :     if (i<m-1) PetscCall(MatSetValue(A,II,II+n,-1.0,INSERT_VALUES));
      51        3675 :     if (j>0) PetscCall(MatSetValue(A,II,II-1,-1.0,INSERT_VALUES));
      52        3675 :     if (j<n-1) PetscCall(MatSetValue(A,II,II+1,-1.0,INSERT_VALUES));
      53        3675 :     PetscCall(MatSetValue(A,II,II,4.0,INSERT_VALUES));
      54             :   }
      55             : 
      56           4 :   PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
      57           4 :   PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
      58             : 
      59             :   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      60             :                 Create the eigensolver and set various options
      61             :      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
      62             : 
      63           4 :   PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
      64           4 :   PetscCall(EPSSetOperators(eps,A,NULL));
      65           4 :   PetscCall(EPSSetProblemType(eps,EPS_HEP));
      66           4 :   PetscCall(EPSSetType(eps,EPSPRIMME));
      67             : 
      68             :   /*
      69             :      Set several options
      70             :   */
      71           4 :   PetscCall(EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL));
      72           4 :   PetscCall(EPSGetST(eps,&st));
      73           4 :   PetscCall(STSetType(st,STPRECOND));
      74           4 :   PetscCall(STGetKSP(st,&ksp));
      75           4 :   PetscCall(KSPGetPC(ksp,&pc));
      76           4 :   PetscCall(KSPSetType(ksp,KSPPREONLY));
      77           4 :   PetscCall(PCSetType(pc,PCICC));
      78             : 
      79           4 :   PetscCall(EPSPRIMMESetBlockSize(eps,4));
      80           4 :   PetscCall(EPSPRIMMESetMethod(eps,EPS_PRIMME_GD_OLSEN_PLUSK));
      81           4 :   PetscCall(EPSSetFromOptions(eps));
      82             : 
      83             :   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      84             :                  Compute eigenvalues and display info
      85             :      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
      86             : 
      87           4 :   PetscCall(EPSSolve(eps));
      88           4 :   PetscCall(EPSPRIMMEGetBlockSize(eps,&bs));
      89           4 :   PetscCall(EPSPRIMMEGetMethod(eps,&meth));
      90           4 :   PetscCall(PetscPrintf(PETSC_COMM_WORLD," PRIMME: using block size %" PetscInt_FMT ", method %s\n",bs,EPSPRIMMEMethods[meth]));
      91             : 
      92           4 :   PetscCall(EPSErrorView(eps,EPS_ERROR_ABSOLUTE,NULL));
      93             : 
      94           4 :   PetscCall(EPSDestroy(&eps));
      95           4 :   PetscCall(MatDestroy(&A));
      96           4 :   PetscCall(SlepcFinalize());
      97             :   return 0;
      98             : }
      99             : 
     100             : /*TEST
     101             : 
     102             :    build:
     103             :       requires: primme
     104             : 
     105             :    testset:
     106             :       args: -eps_nev 4
     107             :       requires: primme
     108             :       output_file: output/test34_1.out
     109             :       test:
     110             :          suffix: 1
     111             :       test:
     112             :          suffix: 2
     113             :          args: -st_pc_type bjacobi -eps_target 0.01 -eps_target_real -eps_refined
     114             :          nsize: 2
     115             :       test:
     116             :          suffix: 3
     117             :          args: -eps_smallest_magnitude -eps_harmonic
     118             : 
     119             : TEST*/