Actual source code: test1f.F90

      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: #include <slepc/finclude/slepcbv.h>
     17: program test1f
     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
     32:   PetscScalar          :: z(KMAX), val
     33:   PetscScalar, pointer :: qq(:, :)
     34:   PetscReal            :: nrm
     35:   PetscBool            :: flg
     36:   PetscErrorCode       :: ierr
     37:   PetscScalar, parameter :: one = 1.0, mone = -1.0, two = 2.0, zero = 0.0

     39: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     40: ! Beginning of program
     41: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

     43:   n = 10
     44:   k = 5
     45:   l = 3
     46:   PetscCallA(SlepcInitialize(PETSC_NULL_CHARACTER, ierr))
     47:   PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))
     48:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-n', n, flg, ierr))
     49:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-k', k, flg, ierr))
     50:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-l', l, flg, ierr))
     51:   PetscCheckA(k <= KMAX, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, 'Program currently limited to k=35')
     52:   if (rank == 0) then
     53:     write (*, '(/a,i3,a,i3,a)') 'Test BV with', k, ' columns of length', n, ' (Fortran)'
     54:   end if

     56: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     57: ! Initialize data
     58: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

     60: ! ** Create template vector
     61:   PetscCallA(VecCreate(PETSC_COMM_WORLD, t, ierr))
     62:   PetscCallA(VecSetSizes(t, PETSC_DECIDE, n, ierr))
     63:   PetscCallA(VecSetFromOptions(t, ierr))

     65: ! ** Create BV object X
     66:   PetscCallA(BVCreate(PETSC_COMM_WORLD, X, ierr))
     67:   PetscCallA(BVSetSizesFromVec(X, t, k, ierr))
     68:   PetscCallA(BVSetFromOptions(X, ierr))

     70: ! ** Fill X entries
     71:   do j = 0, k - 1
     72:     PetscCallA(BVGetColumn(X, j, v, ierr))
     73:     PetscCallA(VecSet(v, zero, ierr))
     74:     do i = 0, 3
     75:       if (i + j < n) then
     76:         val = 3*i + j - 2
     77:         PetscCallA(VecSetValue(v, i + j, val, INSERT_VALUES, ierr))
     78:       end if
     79:     end do
     80:     PetscCallA(VecAssemblyBegin(v, ierr))
     81:     PetscCallA(VecAssemblyEnd(v, ierr))
     82:     PetscCallA(BVRestoreColumn(X, j, v, ierr))
     83:   end do

     85: ! ** Create BV object Y
     86:   PetscCallA(BVCreate(PETSC_COMM_WORLD, Y, ierr))
     87:   PetscCallA(BVSetSizesFromVec(Y, t, l, ierr))
     88:   PetscCallA(BVSetFromOptions(Y, ierr))

     90: ! ** Fill Y entries
     91:   do j = 0, l - 1
     92:     PetscCallA(BVGetColumn(Y, j, v, ierr))
     93:     val = real(j + 1, PETSC_REAL_KIND)/4.0
     94:     PetscCallA(VecSet(v, val, ierr))
     95:     PetscCallA(BVRestoreColumn(Y, j, v, ierr))
     96:   end do

     98: ! ** Create Mat
     99:   PetscCallA(MatCreateSeqDense(PETSC_COMM_SELF, k, l, PETSC_NULL_SCALAR_ARRAY, Q, ierr))
    100:   PetscCallA(MatDenseGetArray(Q, qq, ierr))
    101:   do i = 1, k
    102:     do j = 1, l
    103:       if (i < j) then
    104:         qq(i, j) = 2.0
    105:       else
    106:         qq(i, j) = -0.5
    107:       end if
    108:     end do
    109:   end do
    110:   PetscCallA(MatDenseRestoreArray(Q, qq, ierr))

    112: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    113: ! Test several operations
    114: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

    116: ! ** Test BVMult
    117:   PetscCallA(BVMult(Y, two, one, X, Q, ierr))

    119: ! ** Test BVMultVec
    120:   PetscCallA(BVGetColumn(Y, 0_PETSC_INT_KIND, v, ierr))
    121:   z(1) = 2.0
    122:   do i = 2, k
    123:     z(i) = -0.5*z(i - 1)
    124:   end do
    125:   PetscCallA(BVMultVec(X, mone, one, v, z, ierr))
    126:   PetscCallA(BVRestoreColumn(Y, 0_PETSC_INT_KIND, v, ierr))

    128: ! ** Test BVDot
    129:   PetscCallA(MatCreateSeqDense(PETSC_COMM_SELF, l, k, PETSC_NULL_SCALAR_ARRAY, M, ierr))
    130:   PetscCallA(BVDot(X, Y, M, ierr))

    132: ! ** Test BVDotVec
    133:   PetscCallA(BVGetColumn(Y, 0_PETSC_INT_KIND, v, ierr))
    134:   PetscCallA(BVDotVec(X, v, z, ierr))
    135:   PetscCallA(BVRestoreColumn(Y, 0_PETSC_INT_KIND, v, ierr))

    137: ! ** Test BVMultInPlace and BVScale
    138:   PetscCallA(BVMultInPlace(X, Q, 1_PETSC_INT_KIND, l, ierr))
    139:   PetscCallA(BVScale(X, two, ierr))

    141: ! ** Test BVNorm
    142:   PetscCallA(BVNormColumn(X, 0_PETSC_INT_KIND, NORM_2, nrm, ierr))
    143:   if (rank == 0) then
    144:     write (*, '(a,f8.4)') '2-Norm of X[0] = ', nrm
    145:   end if
    146:   PetscCallA(BVNorm(X, NORM_FROBENIUS, nrm, ierr))
    147:   if (rank == 0) then
    148:     write (*, '(a,f8.4)') 'Frobenius Norm of X = ', nrm
    149:   end if

    151: ! *** Clean up
    152:   PetscCallA(BVDestroy(X, ierr))
    153:   PetscCallA(BVDestroy(Y, ierr))
    154:   PetscCallA(VecDestroy(t, ierr))
    155:   PetscCallA(MatDestroy(Q, ierr))
    156:   PetscCallA(MatDestroy(M, ierr))
    157:   PetscCallA(SlepcFinalize(ierr))
    158: end program test1f

    160: !/*TEST
    161: !
    162: !   test:
    163: !      suffix: 1
    164: !      nsize: 1
    165: !      args: -bv_type {{vecs contiguous svec mat}separate output}
    166: !      output_file: output/test1f_1.out
    167: !
    168: !   test:
    169: !      suffix: 2
    170: !      nsize: 2
    171: !      args: -bv_type {{vecs contiguous svec mat}separate output}
    172: !      output_file: output/test1f_1.out
    173: !
    174: !TEST*/