Actual source code: test1f.F90

slepc-3.11.2 2019-07-30
Report Typos and Errors
  1: !
  2: !  - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3: !  SLEPc - Scalable Library for Eigenvalue Problem Computations
  4: !  Copyright (c) 2002-2019, 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 PetscObjectSetName(X,'X',ierr);CHKERRA(ierr)
 79:       call BVSetSizesFromVec(X,t,k,ierr);CHKERRA(ierr)
 80:       call BVSetFromOptions(X,ierr);CHKERRA(ierr)

 82: !     ** Fill X entries
 83:       do j=0,k-1
 84:         call BVGetColumn(X,j,v,ierr);CHKERRA(ierr)
 85:         call VecSet(v,zero,ierr);CHKERRA(ierr)
 86:         do i=0,3
 87:           if (i+j<n) then
 88:             val = 3*i+j-2
 89:             call VecSetValue(v,i+j,val,INSERT_VALUES,ierr);CHKERRA(ierr)
 90:           end if
 91:         end do
 92:         call VecAssemblyBegin(v,ierr);CHKERRA(ierr)
 93:         call VecAssemblyEnd(v,ierr);CHKERRA(ierr)
 94:         call BVRestoreColumn(X,j,v,ierr);CHKERRA(ierr)
 95:       end do

 97: !     ** Create BV object Y
 98:       call BVCreate(PETSC_COMM_WORLD,Y,ierr);CHKERRA(ierr)
 99:       call PetscObjectSetName(Y,'Y',ierr);CHKERRA(ierr)
100:       call BVSetSizesFromVec(Y,t,l,ierr);CHKERRA(ierr)
101:       call BVSetFromOptions(Y,ierr);CHKERRA(ierr)

103: !     ** Fill Y entries
104:       do j=0,l-1
105:         call BVGetColumn(Y,j,v,ierr);CHKERRA(ierr)
106:         val = real(j+1)/4.0
107:         call VecSet(v,val,ierr);CHKERRA(ierr)
108:         call BVRestoreColumn(Y,j,v,ierr);CHKERRA(ierr)
109:       end do

111: !     ** Create Mat
112:       call MatCreateSeqDense(PETSC_COMM_SELF,k,l,PETSC_NULL_SCALAR,Q,ierr);CHKERRA(ierr)
113:       call PetscObjectSetName(Q,'Q',ierr);CHKERRA(ierr)
114:       call MatDenseGetArray(Q,qq,iq,ierr);CHKERRA(ierr)
115:       do i=0,k-1
116:         do j=0,l-1
117:           if (i<j) then
118:             qq(iq+1+i+j*k) = 2.0
119:           else
120:             qq(iq+1+i+j*k) = -0.5
121:           end if
122:         end do
123:       end do
124:       call MatDenseRestoreArray(Q,qq,iq,ierr);CHKERRA(ierr)

126: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
127: !     Test several operations
128: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

130: !     ** Test BVMult
131:       call BVMult(Y,two,one,X,Q,ierr);CHKERRA(ierr)

133: !     ** Test BVMultVec
134:       call BVGetColumn(Y,izero,v,ierr);CHKERRA(ierr)
135:       z(1) = 2.0
136:       do i=2,k
137:         z(i) = -0.5*z(i-1)
138:       end do
139:       call BVMultVec(X,mone,one,v,z,ierr);CHKERRA(ierr)
140:       call BVRestoreColumn(Y,izero,v,ierr);CHKERRA(ierr)

142: !     ** Test BVDot
143:       call MatCreateSeqDense(PETSC_COMM_SELF,l,k,PETSC_NULL_SCALAR,M,ierr);CHKERRA(ierr)
144:       call PetscObjectSetName(M,'M',ierr);CHKERRA(ierr)
145:       call BVDot(X,Y,M,ierr);CHKERRA(ierr)

147: !     ** Test BVDotVec
148:       call BVGetColumn(Y,izero,v,ierr);CHKERRA(ierr)
149:       call BVDotVec(X,v,z,ierr);CHKERRA(ierr)
150:       call BVRestoreColumn(Y,izero,v,ierr);CHKERRA(ierr)

152: !     ** Test BVMultInPlace and BVScale
153:       call BVMultInPlace(X,Q,ione,l,ierr);CHKERRA(ierr)
154:       call BVScale(X,two,ierr);CHKERRA(ierr)

156: !     ** Test BVNorm
157:       call BVNormColumn(X,izero,NORM_2,nrm,ierr);CHKERRA(ierr)
158:       if (rank .eq. 0) then
159:         write(*,120) nrm
160:       endif
161:  120  format ('2-Norm of X[0] = ',f8.4)
162:       call BVNorm(X,NORM_FROBENIUS,nrm,ierr);CHKERRA(ierr)
163:       if (rank .eq. 0) then
164:         write(*,130) nrm
165:       endif
166:  130  format ('Frobenius Norm of X = ',f8.4)

168: !     *** Clean up
169:       call BVDestroy(X,ierr);CHKERRA(ierr)
170:       call BVDestroy(Y,ierr);CHKERRA(ierr)
171:       call VecDestroy(t,ierr);CHKERRA(ierr)
172:       call MatDestroy(Q,ierr);CHKERRA(ierr)
173:       call MatDestroy(M,ierr);CHKERRA(ierr)
174:       call SlepcFinalize(ierr)
175:       end

177: !/*TEST
178: !
179: !   test:
180: !      suffix: 1
181: !      nsize: 1
182: !      args: -bv_type {{vecs contiguous svec mat}separate output}
183: !      output_file: output/test1f_1.out
184: !
185: !   test:
186: !      suffix: 2
187: !      nsize: 2
188: !      args: -bv_type {{vecs contiguous svec mat}separate output}
189: !      output_file: output/test1f_1.out
190: !
191: !TEST*/