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[] = "Illustrates passing a sparser matrix to build the preconditioner.\n\n"
      12             :   "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           2 : int main(int argc,char **argv)
      19             : {
      20           2 :   Mat            A,A0;            /* operator matrix */
      21           2 :   EPS            eps;             /* eigenproblem solver context */
      22           2 :   ST             st;
      23           2 :   PetscInt       N,n=24,m,Istart,Iend,II,i,j;
      24           2 :   PetscBool      flag,terse;
      25             : 
      26           2 :   PetscFunctionBeginUser;
      27           2 :   PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
      28             : 
      29           2 :   PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL));
      30           2 :   PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,&flag));
      31           2 :   if (!flag) m=n;
      32           2 :   N = n*m;
      33           2 :   PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nModified 2-D Laplacian Eigenproblem, N=%" PetscInt_FMT " (%" PetscInt_FMT "x%" PetscInt_FMT " grid)\n\n",N,n,m));
      34             : 
      35             :   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      36             :        Compute the operator matrix A and a sparser approximation A0
      37             :      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
      38             : 
      39           2 :   PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
      40           2 :   PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N));
      41           2 :   PetscCall(MatSetFromOptions(A));
      42           2 :   PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
      43        1154 :   for (II=Istart;II<Iend;II++) {
      44        1152 :     i = II/n; j = II-i*n;
      45        1152 :     if (i>0) PetscCall(MatSetValue(A,II,II-n,-0.2,INSERT_VALUES));
      46        1152 :     if (i<m-1) PetscCall(MatSetValue(A,II,II+n,-0.2,INSERT_VALUES));
      47        1152 :     if (j>0) PetscCall(MatSetValue(A,II,II-1,-3.0,INSERT_VALUES));
      48        1152 :     if (j<n-1) PetscCall(MatSetValue(A,II,II+1,-3.0,INSERT_VALUES));
      49        1152 :     PetscCall(MatSetValue(A,II,II,7.0,INSERT_VALUES));
      50             :   }
      51           2 :   PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
      52           2 :   PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
      53             : 
      54           2 :   PetscCall(MatCreate(PETSC_COMM_WORLD,&A0));
      55           2 :   PetscCall(MatSetSizes(A0,PETSC_DECIDE,PETSC_DECIDE,N,N));
      56           2 :   PetscCall(MatSetFromOptions(A0));
      57           2 :   PetscCall(MatGetOwnershipRange(A0,&Istart,&Iend));
      58        1154 :   for (II=Istart;II<Iend;II++) {
      59        1152 :     i = II/n; j = II-i*n;
      60        1152 :     if (j>0) PetscCall(MatSetValue(A0,II,II-1,-3.0,INSERT_VALUES));
      61        1152 :     if (j<n-1) PetscCall(MatSetValue(A0,II,II+1,-3.0,INSERT_VALUES));
      62        1152 :     PetscCall(MatSetValue(A0,II,II,7.0,INSERT_VALUES));
      63             :   }
      64           2 :   PetscCall(MatAssemblyBegin(A0,MAT_FINAL_ASSEMBLY));
      65           2 :   PetscCall(MatAssemblyEnd(A0,MAT_FINAL_ASSEMBLY));
      66             : 
      67             :   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      68             :                 Create the eigensolver and set various options
      69             :      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
      70             : 
      71           2 :   PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
      72           2 :   PetscCall(EPSSetOperators(eps,A,NULL));
      73           2 :   PetscCall(EPSSetProblemType(eps,EPS_HEP));
      74           2 :   PetscCall(EPSGetST(eps,&st));
      75           2 :   PetscCall(STSetType(st,STSINVERT));
      76           2 :   PetscCall(STSetPreconditionerMat(st,A0));
      77           2 :   PetscCall(EPSSetTarget(eps,0.0));
      78           2 :   PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE));
      79           2 :   PetscCall(EPSSetFromOptions(eps));
      80             : 
      81             :   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      82             :                  Solve the eigensystem and display solution
      83             :      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
      84             : 
      85           2 :   PetscCall(EPSSolve(eps));
      86             : 
      87             :   /* show detailed info unless -terse option is given by user */
      88           2 :   PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse));
      89           2 :   if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL));
      90             :   else {
      91           0 :     PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL));
      92           0 :     PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD));
      93           0 :     PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD));
      94           0 :     PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
      95             :   }
      96           2 :   PetscCall(EPSDestroy(&eps));
      97           2 :   PetscCall(MatDestroy(&A));
      98           2 :   PetscCall(MatDestroy(&A0));
      99           2 :   PetscCall(SlepcFinalize());
     100             :   return 0;
     101             : }
     102             : 
     103             : /*TEST
     104             : 
     105             :    testset:
     106             :       args: -eps_nev 4 -terse
     107             :       output_file: output/ex46_1.out
     108             :       requires: !single
     109             :       test:
     110             :          suffix: 1
     111             :       test:
     112             :          suffix: 2
     113             :          args: -st_ksp_type bcgs -st_pc_type bjacobi
     114             : 
     115             : TEST*/