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: !  Program usage: mpiexec -n <np> ./test7f [-help] [-n <n>] [all SLEPc options]
 11: !
 12: !  Description: Simple example that solves an eigensystem with the EPS object.
 13: !  Same problem as ex1f but with simplified output.
 14: !
 15: !  The command line options are:
 16: !    -n <n>, where <n> = number of grid points = matrix size
 17: !
 18: ! ----------------------------------------------------------------------
 19: !
 20: #include <slepc/finclude/slepceps.h>
 21: program test7f
 22:   use slepceps
 23:   implicit none

 25: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 26: ! Declarations
 27: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 29:   Mat            :: A     ! operator matrix
 30:   EPS            :: eps   ! eigenproblem solver context
 31:   EPSType        :: tname
 32:   PetscInt       :: n, i, Istart, Iend
 33:   PetscInt       :: nev, nini
 34:   PetscInt       :: col(3)
 35:   PetscInt       :: i0, i1, i2, i3
 36:   PetscMPIInt    :: rank
 37:   PetscErrorCode :: ierr
 38:   PetscBool      :: flg
 39:   PetscScalar    :: val(3)
 40:   Vec            :: v
 41:   PetscScalar, parameter :: one = 1.0

 43: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 44: ! Beginning of program
 45: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 47:   PetscCallA(SlepcInitialize(PETSC_NULL_CHARACTER, ierr))
 48:   PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))
 49:   n = 30
 50:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-n', n, flg, ierr))

 52:   if (rank == 0) then
 53:     write (*, '(/a,i3,a)') '1-D Laplacian Eigenproblem, n =', n, ' (Fortran)'
 54:   end if

 56: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 57: ! Compute the operator matrix that defines the eigensystem, Ax=kx
 58: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 60:   PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr))
 61:   PetscCallA(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, n, n, ierr))
 62:   PetscCallA(MatSetFromOptions(A, ierr))

 64:   i0 = 0
 65:   i1 = 1
 66:   i2 = 2
 67:   i3 = 3
 68:   PetscCallA(MatGetOwnershipRange(A, Istart, Iend, ierr))
 69:   if (Istart == 0) then
 70:     i = 0
 71:     col(1) = 0
 72:     col(2) = 1
 73:     val(1) = 2.0
 74:     val(2) = -1.0
 75:     PetscCallA(MatSetValues(A, i1, [i], i2, col, val, INSERT_VALUES, ierr))
 76:     Istart = Istart + 1
 77:   end if
 78:   if (Iend == n) then
 79:     i = n - 1
 80:     col(1) = n - 2
 81:     col(2) = n - 1
 82:     val(1) = -1.0
 83:     val(2) = 2.0
 84:     PetscCallA(MatSetValues(A, i1, [i], i2, col, val, INSERT_VALUES, ierr))
 85:     Iend = Iend - 1
 86:   end if
 87:   val(1) = -1.0
 88:   val(2) = 2.0
 89:   val(3) = -1.0
 90:   do i = Istart, Iend - 1
 91:     col(1) = i - 1
 92:     col(2) = i
 93:     col(3) = i + 1
 94:     PetscCallA(MatSetValues(A, i1, [i], i3, col, val, INSERT_VALUES, ierr))
 95:   end do

 97:   PetscCallA(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY, ierr))
 98:   PetscCallA(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY, ierr))

100:   PetscCallA(MatCreateVecs(A, v, PETSC_NULL_VEC, ierr))
101:   if (Istart == 0) then
102:     PetscCallA(VecSetValue(v, i0, one, INSERT_VALUES, ierr))
103:   end if
104:   PetscCallA(VecAssemblyBegin(v, ierr))
105:   PetscCallA(VecAssemblyEnd(v, ierr))

107: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
108: ! Create the eigensolver and display info
109: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

111: ! ** Create eigensolver context
112:   PetscCallA(EPSCreate(PETSC_COMM_WORLD, eps, ierr))

114: ! ** Set operators. In this case, it is a standard eigenvalue problem
115:   PetscCallA(EPSSetOperators(eps, A, PETSC_NULL_MAT, ierr))
116:   PetscCallA(EPSSetProblemType(eps, EPS_HEP, ierr))

118: ! ** Set solver parameters at runtime
119:   PetscCallA(EPSSetFromOptions(eps, ierr))

121: ! ** Set initial vectors
122:   nini = 1
123:   PetscCallA(EPSSetInitialSpace(eps, nini, [v], ierr))

125: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
126: !     Solve the eigensystem
127: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

129:   PetscCallA(EPSSolve(eps, ierr))

131: ! ** Optional: Get some information from the solver and display it
132:   PetscCallA(EPSGetType(eps, tname, ierr))
133:   if (rank == 0) then
134:     write (*, '(a,a)') ' Solution method: ', tname
135:   end if
136:   PetscCallA(EPSGetDimensions(eps, nev, PETSC_NULL_INTEGER, PETSC_NULL_INTEGER, ierr))
137:   if (rank == 0) then
138:     write (*, '(a,i2)') ' Number of requested eigenvalues:', nev
139:   end if

141: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
142: ! Display solution and clean up
143: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

145:   PetscCallA(EPSErrorView(eps, EPS_ERROR_RELATIVE, PETSC_NULL_VIEWER, ierr))
146:   PetscCallA(EPSDestroy(eps, ierr))
147:   PetscCallA(MatDestroy(A, ierr))
148:   PetscCallA(VecDestroy(v, ierr))

150:   PetscCallA(SlepcFinalize(ierr))
151: end program test7f

153: !/*TEST
154: !
155: !   test:
156: !      suffix: 1
157: !      args: -eps_nev 4 -eps_ncv 19
158: !      filter: sed -e "s/83791/83792/"
159: !
160: !TEST*/