Actual source code: test1f.F90
slepc-3.17.0 2022-03-31
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> ./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; SETERRA(PETSC_COMM_SELF,1,'Program currently limited to k=35'); endif
62: if (rank .eq. 0) then
63: write(*,110) k,n
64: endif
65: 110 format (/'Test BV with',I3,' columns of length',I3,' (Fortran)')
67: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
68: ! Initialize data
69: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
71: ! ** Create template vector
72: call VecCreate(PETSC_COMM_WORLD,t,ierr);CHKERRA(ierr)
73: call VecSetSizes(t,PETSC_DECIDE,n,ierr);CHKERRA(ierr)
74: call VecSetFromOptions(t,ierr);CHKERRA(ierr)
76: ! ** Create BV object X
77: call BVCreate(PETSC_COMM_WORLD,X,ierr);CHKERRA(ierr)
78: call BVSetSizesFromVec(X,t,k,ierr);CHKERRA(ierr)
79: call BVSetFromOptions(X,ierr);CHKERRA(ierr)
81: ! ** Fill X entries
82: do j=0,k-1
83: call BVGetColumn(X,j,v,ierr);CHKERRA(ierr)
84: call VecSet(v,zero,ierr);CHKERRA(ierr)
85: do i=0,3
86: if (i+j<n) then
87: val = 3*i+j-2
88: call VecSetValue(v,i+j,val,INSERT_VALUES,ierr);CHKERRA(ierr)
89: end if
90: end do
91: call VecAssemblyBegin(v,ierr);CHKERRA(ierr)
92: call VecAssemblyEnd(v,ierr);CHKERRA(ierr)
93: call BVRestoreColumn(X,j,v,ierr);CHKERRA(ierr)
94: end do
96: ! ** Create BV object Y
97: call BVCreate(PETSC_COMM_WORLD,Y,ierr);CHKERRA(ierr)
98: call BVSetSizesFromVec(Y,t,l,ierr);CHKERRA(ierr)
99: call BVSetFromOptions(Y,ierr);CHKERRA(ierr)
101: ! ** Fill Y entries
102: do j=0,l-1
103: call BVGetColumn(Y,j,v,ierr);CHKERRA(ierr)
104: val = real(j+1)/4.0
105: call VecSet(v,val,ierr);CHKERRA(ierr)
106: call BVRestoreColumn(Y,j,v,ierr);CHKERRA(ierr)
107: end do
109: ! ** Create Mat
110: call MatCreateSeqDense(PETSC_COMM_SELF,k,l,PETSC_NULL_SCALAR,Q,ierr);CHKERRA(ierr)
111: call MatDenseGetArray(Q,qq,iq,ierr);CHKERRA(ierr)
112: do i=0,k-1
113: do j=0,l-1
114: if (i<j) then
115: qq(iq+1+i+j*k) = 2.0
116: else
117: qq(iq+1+i+j*k) = -0.5
118: end if
119: end do
120: end do
121: call MatDenseRestoreArray(Q,qq,iq,ierr);CHKERRA(ierr)
123: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
124: ! Test several operations
125: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
127: ! ** Test BVMult
128: call BVMult(Y,two,one,X,Q,ierr);CHKERRA(ierr)
130: ! ** Test BVMultVec
131: call BVGetColumn(Y,izero,v,ierr);CHKERRA(ierr)
132: z(1) = 2.0
133: do i=2,k
134: z(i) = -0.5*z(i-1)
135: end do
136: call BVMultVec(X,mone,one,v,z,ierr);CHKERRA(ierr)
137: call BVRestoreColumn(Y,izero,v,ierr);CHKERRA(ierr)
139: ! ** Test BVDot
140: call MatCreateSeqDense(PETSC_COMM_SELF,l,k,PETSC_NULL_SCALAR,M,ierr);CHKERRA(ierr)
141: call BVDot(X,Y,M,ierr);CHKERRA(ierr)
143: ! ** Test BVDotVec
144: call BVGetColumn(Y,izero,v,ierr);CHKERRA(ierr)
145: call BVDotVec(X,v,z,ierr);CHKERRA(ierr)
146: call BVRestoreColumn(Y,izero,v,ierr);CHKERRA(ierr)
148: ! ** Test BVMultInPlace and BVScale
149: call BVMultInPlace(X,Q,ione,l,ierr);CHKERRA(ierr)
150: call BVScale(X,two,ierr);CHKERRA(ierr)
152: ! ** Test BVNorm
153: call BVNormColumn(X,izero,NORM_2,nrm,ierr);CHKERRA(ierr)
154: if (rank .eq. 0) then
155: write(*,120) nrm
156: endif
157: 120 format ('2-Norm of X[0] = ',f8.4)
158: call BVNorm(X,NORM_FROBENIUS,nrm,ierr);CHKERRA(ierr)
159: if (rank .eq. 0) then
160: write(*,130) nrm
161: endif
162: 130 format ('Frobenius Norm of X = ',f8.4)
164: ! *** Clean up
165: call BVDestroy(X,ierr);CHKERRA(ierr)
166: call BVDestroy(Y,ierr);CHKERRA(ierr)
167: call VecDestroy(t,ierr);CHKERRA(ierr)
168: call MatDestroy(Q,ierr);CHKERRA(ierr)
169: call MatDestroy(M,ierr);CHKERRA(ierr)
170: call SlepcFinalize(ierr)
171: end
173: !/*TEST
174: !
175: ! test:
176: ! suffix: 1
177: ! nsize: 1
178: ! args: -bv_type {{vecs contiguous svec mat}separate output}
179: ! output_file: output/test1f_1.out
180: !
181: ! test:
182: ! suffix: 2
183: ! nsize: 2
184: ! args: -bv_type {{vecs contiguous svec mat}separate output}
185: ! output_file: output/test1f_1.out
186: !
187: !TEST*/