Actual source code: ex1f90.F90
slepc-3.6.3 2016-03-29
1: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2: ! SLEPc - Scalable Library for Eigenvalue Problem Computations
3: ! Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain
4: !
5: ! This file is part of SLEPc.
6: !
7: ! SLEPc is free software: you can redistribute it and/or modify it under the
8: ! terms of version 3 of the GNU Lesser General Public License as published by
9: ! the Free Software Foundation.
10: !
11: ! SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
12: ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
13: ! FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
14: ! more details.
15: !
16: ! You should have received a copy of the GNU Lesser General Public License
17: ! along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
18: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
19: !
20: ! Program usage: mpirun -np n ex1f90 [-help] [-n <n>] [all SLEPc options]
21: !
22: ! Description: Simple example that solves an eigensystem with the EPS object.
23: ! The standard symmetric eigenvalue problem to be solved corresponds to the
24: ! Laplacian operator in 1 dimension.
25: !
26: ! The command line options are:
27: ! -n <n>, where <n> = number of grid points = matrix size
28: !
29: ! ----------------------------------------------------------------------
30: !
31: program main
33: #include <slepc/finclude/slepcepsdef.h>
34: use slepceps
36: implicit none
38: ! For usage without modules, uncomment the following lines and remove
39: ! the previous lines between 'program main' and 'implicit none'
40: !
41: !#include <petsc/finclude/petsc.h>
42: !#include <slepc/finclude/slepc.h>
44: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
45: ! Declarations
46: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
47: !
48: ! Variables:
49: ! A operator matrix
50: ! eps eigenproblem solver context
52: #if defined(PETSC_USE_FORTRAN_DATATYPES)
53: type(Mat) A
54: type(EPS) eps
55: #else
56: Mat A
57: EPS eps
58: #endif
59: EPSType tname
60: PetscInt n, i, Istart, Iend, one, two, three
61: PetscInt nev
62: PetscInt row(1), col(3)
63: PetscMPIInt rank
64: PetscErrorCode ierr
65: PetscBool flg, terse
66: PetscScalar value(3)
68: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
69: ! Beginning of program
70: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
72: one = 1
73: two = 2
74: three = 3
75: call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
76: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
77: n = 30
78: call PetscOptionsGetInt(PETSC_NULL_CHARACTER,'-n',n,flg,ierr)
80: if (rank .eq. 0) then
81: write(*,100) n
82: endif
83: 100 format (/'1-D Laplacian Eigenproblem, n =',I4,' (Fortran)')
85: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
86: ! Compute the operator matrix that defines the eigensystem, Ax=kx
87: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
89: call MatCreate(PETSC_COMM_WORLD,A,ierr)
90: call MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)
91: call MatSetFromOptions(A,ierr)
92: call MatSetUp(A,ierr)
94: call MatGetOwnershipRange(A,Istart,Iend,ierr)
95: if (Istart .eq. 0) then
96: row(1) = 0
97: col(1) = 0
98: col(2) = 1
99: value(1) = 2.0
100: value(2) = -1.0
101: call MatSetValues(A,one,row,two,col,value,INSERT_VALUES,ierr)
102: Istart = Istart+1
103: endif
104: if (Iend .eq. n) then
105: row(1) = n-1
106: col(1) = n-2
107: col(2) = n-1
108: value(1) = -1.0
109: value(2) = 2.0
110: call MatSetValues(A,one,row,two,col,value,INSERT_VALUES,ierr)
111: Iend = Iend-1
112: endif
113: value(1) = -1.0
114: value(2) = 2.0
115: value(3) = -1.0
116: do i=Istart,Iend-1
117: row(1) = i
118: col(1) = i-1
119: col(2) = i
120: col(3) = i+1
121: call MatSetValues(A,one,row,three,col,value,INSERT_VALUES,ierr)
122: enddo
124: call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
125: call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
127: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
128: ! Create the eigensolver and display info
129: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
131: ! ** Create eigensolver context
132: call EPSCreate(PETSC_COMM_WORLD,eps,ierr)
134: ! ** Set operators. In this case, it is a standard eigenvalue problem
135: call EPSSetOperators(eps,A,PETSC_NULL_OBJECT,ierr)
136: call EPSSetProblemType(eps,EPS_HEP,ierr)
138: ! ** Set solver parameters at runtime
139: call EPSSetFromOptions(eps,ierr)
141: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
142: ! Solve the eigensystem
143: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
145: call EPSSolve(eps,ierr)
147: ! ** Optional: Get some information from the solver and display it
148: call EPSGetType(eps,tname,ierr)
149: if (rank .eq. 0) then
150: write(*,120) tname
151: endif
152: 120 format (' Solution method: ',A)
153: call EPSGetDimensions(eps,nev,PETSC_NULL_INTEGER, &
154: & PETSC_NULL_INTEGER,ierr)
155: if (rank .eq. 0) then
156: write(*,130) nev
157: endif
158: 130 format (' Number of requested eigenvalues:',I4)
160: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
161: ! Display solution and clean up
162: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
164: ! ** show detailed info unless -terse option is given by user
165: call PetscOptionsHasName(PETSC_NULL_CHARACTER,'-terse',terse,ierr)
166: if (terse) then
167: call EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_NULL_OBJECT,ierr)
168: else
169: call PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, &
170: & PETSC_VIEWER_ASCII_INFO_DETAIL,ierr)
171: call EPSReasonView(eps,PETSC_VIEWER_STDOUT_WORLD,ierr)
172: call EPSErrorView(eps,EPS_ERROR_RELATIVE, &
173: & PETSC_VIEWER_STDOUT_WORLD,ierr)
174: call PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD,ierr)
175: endif
176: call EPSDestroy(eps,ierr)
177: call MatDestroy(A,ierr)
179: call SlepcFinalize(ierr)
180: end