Actual source code: ex20f90.F90

      1: program main

      3: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      4: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      5: !
      6: !
      7: !     This examples uses Fortran 90 MODULES instead of include files
      8: !
      9: #include <petsc/finclude/petscvec.h>
     10:   use petscvec
     11:   implicit none

     13: !
     14: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     15: !                   Variable declarations
     16: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     17: !
     18: !  Variables:
     19: !     x, y, w - vectors
     20: !     z       - array of vectors
     21: !
     22:   type(tVec) x, y, w
     23:   type(tVec), pointer :: z(:)

     25:   PetscReal norm, v, v1, v2, tol
     26:   PetscInt n, ithree
     27:   PetscErrorCode ierr
     28:   PetscMPIInt rank
     29:   PetscBool flg
     30:   PetscScalar one, two, three
     31:   PetscScalar dots(3), dot
     32:   PetscReal nfloat

     34: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     35: !                 Beginning of program
     36: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

     38:   PetscCallA(PetscInitialize(ierr))
     39:   tol = 1.e-10_PETSC_REAL_KIND
     40:   one = 1.0
     41:   two = 2.0
     42:   three = 3.0
     43:   n = 20
     44:   ithree = 3

     46:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-n', n, flg, ierr))
     47:   nfloat = n
     48:   PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))

     50: !  Create a vector, specifying only its global dimension.
     51: !  When using VecCreate(), VecSetSizes() and VecSetFromOptions(),
     52: !  the vector format (currently parallel
     53: !  or sequential) is determined at runtime.  Also, the parallel
     54: !  partitioning of the vector is determined by PETSc at runtime.
     55: !
     56:   PetscCallA(VecCreate(PETSC_COMM_WORLD, x, ierr))
     57:   PetscCallA(VecSetSizes(x, PETSC_DECIDE, n, ierr))
     58:   PetscCallA(VecSetFromOptions(x, ierr))

     60: !  Duplicate some work vectors (of the same format and
     61: !  partitioning as the initial vector).

     63:   PetscCallA(VecDuplicate(x, y, ierr))
     64:   PetscCallA(VecDuplicate(x, w, ierr))

     66: !  Duplicate more work vectors (of the same format and
     67: !  partitioning as the initial vector).  Here we duplicate
     68: !  an array of vectors, which is often more convenient than
     69: !  duplicating individual ones.

     71:   PetscCallA(VecDuplicateVecs(x, ithree, z, ierr))

     73: !  Set the vectors to entries to a constant value.

     75:   PetscCallA(VecSet(x, one, ierr))
     76:   PetscCallA(VecSet(y, two, ierr))
     77:   PetscCallA(VecSet(z(1), one, ierr))
     78:   PetscCallA(VecSet(z(2), two, ierr))
     79:   PetscCallA(VecSet(z(3), three, ierr))

     81: !  Demonstrate various basic vector routines.

     83:   PetscCallA(VecDot(x, x, dot, ierr))
     84:   PetscCallA(VecMDot(x, ithree, z, dots, ierr))

     86: !  Note: If using a complex numbers version of PETSc, then
     87: !  PETSC_USE_COMPLEX is defined in the makefiles; otherwise,
     88: !  (when using real numbers) it is undefined.

     90:   if (rank == 0) then
     91: #if defined(PETSC_USE_COMPLEX)
     92:     write (6, 100) int(PetscRealPart(dot))
     93:     write (6, 110) int(PetscRealPart(dots(1))), int(PetscRealPart(dots(2))), int(PetscRealPart(dots(3)))
     94: #else
     95:     write (6, 100) int(dot)
     96:     write (6, 110) int(dots(1)), int(dots(2)), int(dots(3))
     97: #endif
     98:     write (6, 120)
     99:   end if
    100: 100 format('Vector length ', i6)
    101: 110 format('Vector length ', 3(i6))
    102: 120 format('All other values should be near zero')

    104:   PetscCallA(VecScale(x, two, ierr))
    105:   PetscCallA(VecNorm(x, NORM_2, norm, ierr))
    106:   v = abs(norm - 2.0*sqrt(nfloat))
    107:   if (v > -tol .and. v < tol) v = 0.0
    108:   if (rank == 0) write (6, 130) v
    109: 130 format('VecScale ', 1pe9.2)

    111:   PetscCallA(VecCopy(x, w, ierr))
    112:   PetscCallA(VecNorm(w, NORM_2, norm, ierr))
    113:   v = abs(norm - 2.0*sqrt(nfloat))
    114:   if (v > -tol .and. v < tol) v = 0.0
    115:   if (rank == 0) write (6, 140) v
    116: 140 format('VecCopy ', 1pe9.2)

    118:   PetscCallA(VecAXPY(y, three, x, ierr))
    119:   PetscCallA(VecNorm(y, NORM_2, norm, ierr))
    120:   v = abs(norm - 8.0*sqrt(nfloat))
    121:   if (v > -tol .and. v < tol) v = 0.0
    122:   if (rank == 0) write (6, 150) v
    123: 150 format('VecAXPY ', 1pe9.2)

    125:   PetscCallA(VecAYPX(y, two, x, ierr))
    126:   PetscCallA(VecNorm(y, NORM_2, norm, ierr))
    127:   v = abs(norm - 18.0*sqrt(nfloat))
    128:   if (v > -tol .and. v < tol) v = 0.0
    129:   if (rank == 0) write (6, 160) v
    130: 160 format('VecAYXP ', 1pe9.2)

    132:   PetscCallA(VecSwap(x, y, ierr))
    133:   PetscCallA(VecNorm(y, NORM_2, norm, ierr))
    134:   v = abs(norm - 2.0*sqrt(nfloat))
    135:   if (v > -tol .and. v < tol) v = 0.0
    136:   if (rank == 0) write (6, 170) v
    137: 170 format('VecSwap ', 1pe9.2)

    139:   PetscCallA(VecNorm(x, NORM_2, norm, ierr))
    140:   v = abs(norm - 18.0*sqrt(nfloat))
    141:   if (v > -tol .and. v < tol) v = 0.0
    142:   if (rank == 0) write (6, 180) v
    143: 180 format('VecSwap ', 1pe9.2)

    145:   PetscCallA(VecWAXPY(w, two, x, y, ierr))
    146:   PetscCallA(VecNorm(w, NORM_2, norm, ierr))
    147:   v = abs(norm - 38.0*sqrt(nfloat))
    148:   if (v > -tol .and. v < tol) v = 0.0
    149:   if (rank == 0) write (6, 190) v
    150: 190 format('VecWAXPY ', 1pe9.2)

    152:   PetscCallA(VecPointwiseMult(w, y, x, ierr))
    153:   PetscCallA(VecNorm(w, NORM_2, norm, ierr))
    154:   v = abs(norm - 36.0*sqrt(nfloat))
    155:   if (v > -tol .and. v < tol) v = 0.0
    156:   if (rank == 0) write (6, 200) v
    157: 200 format('VecPointwiseMult ', 1pe9.2)

    159:   PetscCallA(VecPointwiseDivide(w, x, y, ierr))
    160:   PetscCallA(VecNorm(w, NORM_2, norm, ierr))
    161:   v = abs(norm - 9.0*sqrt(nfloat))
    162:   if (v > -tol .and. v < tol) v = 0.0
    163:   if (rank == 0) write (6, 210) v
    164: 210 format('VecPointwiseDivide ', 1pe9.2)

    166:   dots(1) = one
    167:   dots(2) = three
    168:   dots(3) = two
    169:   PetscCallA(VecSet(x, one, ierr))
    170:   PetscCallA(VecMAXPY(x, ithree, dots, z, ierr))
    171:   PetscCallA(VecNorm(z(1), NORM_2, norm, ierr))
    172:   v = abs(norm - sqrt(nfloat))
    173:   if (v > -tol .and. v < tol) v = 0.0
    174:   PetscCallA(VecNorm(z(2), NORM_2, norm, ierr))
    175:   v1 = abs(norm - 2.0*sqrt(nfloat))
    176:   if (v1 > -tol .and. v1 < tol) v1 = 0.0
    177:   PetscCallA(VecNorm(z(3), NORM_2, norm, ierr))
    178:   v2 = abs(norm - 3.0*sqrt(nfloat))
    179:   if (v2 > -tol .and. v2 < tol) v2 = 0.0
    180:   if (rank == 0) write (6, 220) v, v1, v2
    181: 220 format('VecMAXPY ', 3(1pe9.2))

    183: !  Free work space.  All PETSc objects should be destroyed when they
    184: !  are no longer needed.

    186:   PetscCallA(VecDestroy(x, ierr))
    187:   PetscCallA(VecDestroy(y, ierr))
    188:   PetscCallA(VecDestroy(w, ierr))
    189:   PetscCallA(VecDestroyVecs(ithree, z, ierr))
    190:   PetscCallA(PetscFinalize(ierr))

    192: end

    194: !/*TEST
    195: !
    196: !     test:
    197: !
    198: !TEST*/