Actual source code: test14f.F90

slepc-3.8.3 2018-04-03
Report Typos and Errors
  1: !
  2: !  - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3: !  SLEPc - Scalable Library for Eigenvalue Problem Computations
  4: !  Copyright (c) 2002-2017, 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> ./test14f [-help] [-n <n>] [all SLEPc options]
 11: !
 12: !  Description: Simple example that tests solving a DSNHEP problem.
 13: !
 14: !  The command line options are:
 15: !    -n <n>, where <n> = matrix size
 16: !
 17: ! ----------------------------------------------------------------------
 18: !
 19:       program main
 20: #include <slepc/finclude/slepcds.h>
 21:       use slepcds
 22:       implicit none

 24: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 25: !     Declarations
 26: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 27: !
 28: !  Variables:
 29: !     A     problem matrix
 30: !     ds    dense solver context

 32:       Mat            A
 33:       DS             ds
 34:       PetscInt       n, i, ld, zero
 35:       PetscMPIInt    rank
 36:       PetscErrorCode ierr
 37:       PetscBool      flg
 38:       PetscScalar    aa(1), wr(100), wi(100)
 39:       PetscReal      re, im
 40:       PetscOffset    ia

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

 46:       zero = 0
 47:       call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
 48:       if (ierr .ne. 0) then
 49:         print*,'SlepcInitialize failed'
 50:         stop
 51:       endif
 52:       call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr);CHKERRA(ierr)
 53:       n = 10
 54:       call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-n',n,flg,ierr);CHKERRA(ierr)
 55:       if (n .gt. 100) then
 56:         if (rank .eq. 0) then
 57:           write(*,100) n
 58:         endif
 59:         SETERRA(PETSC_COMM_SELF,1,' ')
 60:       endif
 61:  100  format (/'Program currently limited to n=100, you set n=',I3)

 63:       if (rank .eq. 0) then
 64:         write(*,110) n
 65:       endif
 66:  110  format (/'Solve a Dense System of type NHEP, n =',I3,' (Fortran)')

 68: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 69: !     Create DS object
 70: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 72:       call DSCreate(PETSC_COMM_WORLD,ds,ierr);CHKERRA(ierr)
 73:       call DSSetType(ds,DSNHEP,ierr);CHKERRA(ierr)
 74:       call DSSetFromOptions(ds,ierr);CHKERRA(ierr)
 75:       ld = n
 76:       call DSAllocate(ds,ld,ierr);CHKERRA(ierr)
 77:       call DSSetDimensions(ds,n,zero,zero,zero,ierr);CHKERRA(ierr)

 79: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 80: !     Fill with Grcar matrix
 81: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 83:       call DSGetMat(ds,DS_MAT_A,A,ierr);CHKERRA(ierr)
 84:       call MatDenseGetArray(A,aa,ia,ierr);CHKERRA(ierr)
 85:       call FillUpMatrix(n,aa(ia+1))
 86:       call MatDenseRestoreArray(A,aa,ia,ierr);CHKERRA(ierr)
 87:       call DSRestoreMat(ds,DS_MAT_A,A,ierr);CHKERRA(ierr)
 88:       call DSSetState(ds,DS_STATE_INTERMEDIATE,ierr);CHKERRA(ierr)

 90: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 91: !     Solve the problem and show eigenvalues
 92: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 94:       call DSSolve(ds,wr,wi,ierr);CHKERRA(ierr)
 95: !     call DSSort(ds,wr,wi,PETSC_NULL_SCALAR,PETSC_NULL_SCALAR,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)

 97:       if (rank .eq. 0) then
 98:         write(*,*) 'Computed eigenvalues ='
 99:         do i=1,n
100: #if defined(PETSC_USE_COMPLEX)
101:           re = PetscRealPart(wr(i))
102:           im = PetscImaginaryPart(wr(i))
103: #else
104:           re = wr(i)
105:           im = wi(i)
106: #endif
107:           if (abs(im).lt.1.d-10) then
108:             write(*,120) re
109:           else
110:             write(*,130) re, im
111:           endif
112:         end do
113:       endif
114:  120  format ('  ',F8.5)
115:  130  format ('  ',F8.5,SP,F8.5,'i')

117: !     *** Clean up
118:       call DSDestroy(ds,ierr);CHKERRA(ierr)
119:       call SlepcFinalize(ierr)
120:       end

122: ! -----------------------------------------------------------------

124:       subroutine FillUpMatrix(n,X)
125:       PetscInt    n,i,j
126:       PetscScalar X(n,n)

128:       do i=2,n
129:         X(i,i-1) = -1.d0
130:       end do
131:       do j=0,3
132:         do i=1,n-j
133:           X(i,i+j) = 1.d0
134:         end do
135:       end do
136:       return
137:       end