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> ./ex1f [-help] [-n <n>] [all SLEPc options]
 11: !
 12: !  Description: Simple example that solves an eigensystem with the EPS object.
 13: !  The standard symmetric eigenvalue problem to be solved corresponds to the
 14: !  Laplacian operator in 1 dimension.
 15: !
 16: !  The command line options are:
 17: !    -n <n>, where <n> = number of grid points = matrix size
 18: !
 19: ! ----------------------------------------------------------------------
 20: !
 21: #include <slepc/finclude/slepceps.h>
 22: program ex1f
 23:   use slepceps
 24:   implicit none

 26: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 27: !     Declarations
 28: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 29: !
 30:   Mat            :: A    ! operator matrix
 31:   EPS            :: eps  ! eigenproblem solver context
 32:   EPSType        :: tname
 33:   PetscInt       :: n, i, Istart, Iend
 34:   PetscInt       :: nev
 35:   PetscInt       :: row(1), col(3)
 36:   PetscMPIInt    :: rank
 37:   PetscErrorCode :: ierr
 38:   PetscBool      :: flg, terse
 39:   PetscScalar    :: val(3)

 41: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 42: ! Beginning of program
 43: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 45:   PetscCallA(SlepcInitialize(PETSC_NULL_CHARACTER, "ex1f test"//c_new_line, ierr))
 46:   PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))
 47:   n = 30
 48:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-n', n, flg, ierr))

 50:   if (rank == 0) then
 51:     write (*, '(/a,i4,a)') '1-D Laplacian Eigenproblem, n =', n, ' (Fortran)'
 52:   end if

 54: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 55: ! Compute the operator matrix that defines the eigensystem, Ax=kx
 56: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 58:   PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr))
 59:   PetscCallA(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, n, n, ierr))
 60:   PetscCallA(MatSetFromOptions(A, ierr))

 62:   PetscCallA(MatGetOwnershipRange(A, Istart, Iend, ierr))
 63:   if (Istart == 0) then
 64:     row(1) = 0
 65:     col(1) = 0
 66:     col(2) = 1
 67:     val(1) = 2.0
 68:     val(2) = -1.0
 69:     PetscCallA(MatSetValues(A, 1_PETSC_INT_KIND, row, 2_PETSC_INT_KIND, col, val, INSERT_VALUES, ierr))
 70:     Istart = Istart + 1
 71:   end if
 72:   if (Iend == n) then
 73:     row(1) = n - 1
 74:     col(1) = n - 2
 75:     col(2) = n - 1
 76:     val(1) = -1.0
 77:     val(2) = 2.0
 78:     PetscCallA(MatSetValues(A, 1_PETSC_INT_KIND, row, 2_PETSC_INT_KIND, col, val, INSERT_VALUES, ierr))
 79:     Iend = Iend - 1
 80:   end if
 81:   val(1) = -1.0
 82:   val(2) = 2.0
 83:   val(3) = -1.0
 84:   do i = Istart, Iend - 1
 85:     row(1) = i
 86:     col(1) = i - 1
 87:     col(2) = i
 88:     col(3) = i + 1
 89:     PetscCallA(MatSetValues(A, 1_PETSC_INT_KIND, row, 3_PETSC_INT_KIND, col, val, INSERT_VALUES, ierr))
 90:   end do

 92:   PetscCallA(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY, ierr))
 93:   PetscCallA(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY, ierr))

 95: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 96: ! Create the eigensolver and display info
 97: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 99: ! ** Create eigensolver context
100:   PetscCallA(EPSCreate(PETSC_COMM_WORLD, eps, ierr))

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

106: ! ** Set solver parameters at runtime
107:   PetscCallA(EPSSetFromOptions(eps, ierr))

109: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
110: ! Solve the eigensystem
111: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

113:   PetscCallA(EPSSolve(eps, ierr))

115: ! ** Optional: Get some information from the solver and display it
116:   PetscCallA(EPSGetType(eps, tname, ierr))
117:   if (rank == 0) then
118:     write (*, '(a,a)') ' Solution method: ', tname
119:   end if
120:   PetscCallA(EPSGetDimensions(eps, nev, PETSC_NULL_INTEGER, PETSC_NULL_INTEGER, ierr))
121:   if (rank == 0) then
122:     write (*, '(a,i4)') ' Number of requested eigenvalues:', nev
123:   end if

125: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
126: ! Display solution and clean up
127: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

129: ! ** show detailed info unless -terse option is given by user
130:   PetscCallA(PetscOptionsHasName(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-terse', terse, ierr))
131:   if (terse) then
132:     PetscCallA(EPSErrorView(eps, EPS_ERROR_RELATIVE, PETSC_NULL_VIEWER, ierr))
133:   else
134:     PetscCallA(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_INFO_DETAIL, ierr))
135:     PetscCallA(EPSConvergedReasonView(eps, PETSC_VIEWER_STDOUT_WORLD, ierr))
136:     PetscCallA(EPSErrorView(eps, EPS_ERROR_RELATIVE, PETSC_VIEWER_STDOUT_WORLD, ierr))
137:     PetscCallA(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD, ierr))
138:   end if
139:   PetscCallA(EPSDestroy(eps, ierr))
140:   PetscCallA(MatDestroy(A, ierr))

142:   PetscCallA(SlepcFinalize(ierr))
143: end program ex1f

145: !/*TEST
146: !
147: !   test:
148: !      args: -eps_nev 4 -terse
149: !      filter: sed -e "s/3.83791/3.83792/"
150: !
151: !TEST*/