Actual source code: test1f.F90
slepc-3.8.3 2018-04-03
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> ./test1f [-help]
11: !
12: ! Description: Simple example that tests BV interface functions.
13: !
14: ! ----------------------------------------------------------------------
15: !
16: program main
17: #include <slepc/finclude/slepcbv.h>
18: use slepcbv
19: implicit none
21: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: ! Declarations
23: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
25: #define KMAX 35
27: Vec t,v
28: Mat Q,M
29: BV X,Y;
30: PetscMPIInt rank
31: PetscInt i,j,n,k,l,izero,ione
32: PetscScalar qq(1),z(KMAX),val
33: PetscScalar one,mone,two,zero
34: PetscOffset iq
35: PetscReal nrm
36: PetscBool flg
37: PetscErrorCode ierr
39: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
40: ! Beginning of program
41: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
43: n = 10
44: k = 5
45: l = 3
46: one = 1.0
47: mone = -1.0
48: two = 2.0
49: zero = 0.0
50: izero = 0
51: ione = 1
52: call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
53: if (ierr .ne. 0) then
54: print*,'SlepcInitialize failed'
55: stop
56: endif
57: call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr);CHKERRA(ierr)
58: call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-n',n,flg,ierr);CHKERRA(ierr)
59: call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-k',k,flg,ierr);CHKERRA(ierr)
60: call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-l',l,flg,ierr);CHKERRA(ierr)
61: if (k .gt. KMAX) then
62: if (rank .eq. 0) then
63: write(*,100) KMAX,k
64: endif
65: SETERRA(PETSC_COMM_SELF,1,' ')
66: endif
67: 100 format (/'Program currently limited to k=',I2,', you set k=',I3)
68: if (rank .eq. 0) then
69: write(*,110) k,n
70: endif
71: 110 format (/'Test BV with',I3,' columns of length',I3,' (Fortran)')
73: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
74: ! Initialize data
75: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
77: ! ** Create template vector
78: call VecCreate(PETSC_COMM_WORLD,t,ierr);CHKERRA(ierr)
79: call VecSetSizes(t,PETSC_DECIDE,n,ierr);CHKERRA(ierr)
80: call VecSetFromOptions(t,ierr);CHKERRA(ierr)
82: ! ** Create BV object X
83: call BVCreate(PETSC_COMM_WORLD,X,ierr);CHKERRA(ierr)
84: call PetscObjectSetName(X,'X',ierr);CHKERRA(ierr)
85: call BVSetSizesFromVec(X,t,k,ierr);CHKERRA(ierr)
86: call BVSetFromOptions(X,ierr);CHKERRA(ierr)
88: ! ** Fill X entries
89: do j=0,k-1
90: call BVGetColumn(X,j,v,ierr);CHKERRA(ierr)
91: call VecSet(v,zero,ierr);CHKERRA(ierr)
92: do i=0,3
93: if (i+j<n) then
94: val = 3*i+j-2
95: call VecSetValue(v,i+j,val,INSERT_VALUES,ierr);CHKERRA(ierr)
96: end if
97: end do
98: call VecAssemblyBegin(v,ierr);CHKERRA(ierr)
99: call VecAssemblyEnd(v,ierr);CHKERRA(ierr)
100: call BVRestoreColumn(X,j,v,ierr);CHKERRA(ierr)
101: end do
103: ! ** Create BV object Y
104: call BVCreate(PETSC_COMM_WORLD,Y,ierr);CHKERRA(ierr)
105: call PetscObjectSetName(Y,'Y',ierr);CHKERRA(ierr)
106: call BVSetSizesFromVec(Y,t,l,ierr);CHKERRA(ierr)
107: call BVSetFromOptions(Y,ierr);CHKERRA(ierr)
109: ! ** Fill Y entries
110: do j=0,l-1
111: call BVGetColumn(Y,j,v,ierr);CHKERRA(ierr)
112: val = real(j+1)/4.0
113: call VecSet(v,val,ierr);CHKERRA(ierr)
114: call BVRestoreColumn(Y,j,v,ierr);CHKERRA(ierr)
115: end do
117: ! ** Create Mat
118: call MatCreateSeqDense(PETSC_COMM_SELF,k,l,PETSC_NULL_SCALAR,Q,ierr);CHKERRA(ierr)
119: call PetscObjectSetName(Q,'Q',ierr);CHKERRA(ierr)
120: call MatDenseGetArray(Q,qq,iq,ierr);CHKERRA(ierr)
121: do i=0,k-1
122: do j=0,l-1
123: if (i<j) then
124: qq(iq+1+i+j*k) = 2.0
125: else
126: qq(iq+1+i+j*k) = -0.5
127: end if
128: end do
129: end do
130: call MatDenseRestoreArray(Q,qq,iq,ierr);CHKERRA(ierr)
132: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
133: ! Test several operations
134: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
136: ! ** Test BVMult
137: call BVMult(Y,two,one,X,Q,ierr);CHKERRA(ierr)
139: ! ** Test BVMultVec
140: call BVGetColumn(Y,izero,v,ierr);CHKERRA(ierr)
141: z(1) = 2.0
142: do i=2,k
143: z(i) = -0.5*z(i-1)
144: end do
145: call BVMultVec(X,mone,one,v,z,ierr);CHKERRA(ierr)
146: call BVRestoreColumn(Y,izero,v,ierr);CHKERRA(ierr)
148: ! ** Test BVDot
149: call MatCreateSeqDense(PETSC_COMM_SELF,l,k,PETSC_NULL_SCALAR,M,ierr);CHKERRA(ierr)
150: call PetscObjectSetName(M,'M',ierr);CHKERRA(ierr)
151: call BVDot(X,Y,M,ierr);CHKERRA(ierr)
153: ! ** Test BVDotVec
154: call BVGetColumn(Y,izero,v,ierr);CHKERRA(ierr)
155: call BVDotVec(X,v,z,ierr);CHKERRA(ierr)
156: call BVRestoreColumn(Y,izero,v,ierr);CHKERRA(ierr)
158: ! ** Test BVMultInPlace and BVScale
159: call BVMultInPlace(X,Q,ione,l,ierr);CHKERRA(ierr)
160: call BVScale(X,two,ierr);CHKERRA(ierr)
162: ! ** Test BVNorm
163: call BVNormColumn(X,izero,NORM_2,nrm,ierr);CHKERRA(ierr)
164: if (rank .eq. 0) then
165: write(*,120) nrm
166: endif
167: 120 format ('2-Norm of X[0] = ',f8.4)
168: call BVNorm(X,NORM_FROBENIUS,nrm,ierr);CHKERRA(ierr)
169: if (rank .eq. 0) then
170: write(*,130) nrm
171: endif
172: 130 format ('Frobenius Norm of X = ',f8.4)
174: ! *** Clean up
175: call BVDestroy(X,ierr);CHKERRA(ierr)
176: call BVDestroy(Y,ierr);CHKERRA(ierr)
177: call VecDestroy(t,ierr);CHKERRA(ierr)
178: call MatDestroy(Q,ierr);CHKERRA(ierr)
179: call MatDestroy(M,ierr);CHKERRA(ierr)
180: call SlepcFinalize(ierr)
181: end