1: !
  2: !
  3: !  Description:  Illustrates the use of VecSetValues() to set
  4: !  multiple values at once; demonstrates VecGetArray().
  5: !
  6: !/*T
  7: !   Concepts: vectors^assembling;
  8: !   Concepts: vectors^arrays of vectors;
  9: !   Processors: 1
 10: !T*/
 11: ! -----------------------------------------------------------------------

 13:       program main
 14:       implicit none

 16: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 17: !                    Include files
 18: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 19: !
 20: !  The following include statements are required for Fortran programs
 21: !  that use PETSc vectors:
 22: !     petscsys.h       - base PETSc routines
 23: !     petscvec.h    - vectors

 25: #include <finclude/petscsys.h>
 26: #include <finclude/petscvec.h>

 28: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 29: !                   Macro definitions
 30: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 31: !
 32: !  Macros to make clearer the process of setting values in vectors and
 33: !  getting values from vectors.
 34: !
 35: !   - The element xx_a(ib) is element ib+1 in the vector x
 36: !   - Here we add 1 to the base array index to facilitate the use of
 37: !     conventional Fortran 1-based array indexing.
 38: !
 39: #define xx_a(ib)  xx_v(xx_i + (ib))
 40: #define yy_a(ib)  yy_v(yy_i + (ib))

 42: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 43: !                 Beginning of program
 44: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 46:        PetscScalar xwork(6)
 47:        PetscScalar xx_v(1),yy_v(1)
 48:        PetscInt     i,n,loc(6),isix
 49:        PetscErrorCode ierr
 50:        PetscOffset xx_i,yy_i
 51:        Vec         x,y

 53:        call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 54:        n = 6
 55:        isix = 6

 57: !  Create initial vector and duplicate it

 59:        call VecCreateSeq(PETSC_COMM_SELF,n,x,ierr)
 60:        call VecDuplicate(x,y,ierr)

 62: !  Fill work arrays with vector entries and locations.  Note that
 63: !  the vector indices are 0-based in PETSc (for both Fortran and
 64: !  C vectors)

 66:        do 10 i=1,n
 67:           loc(i) = i-1
 68:           xwork(i) = 10.0*i
 69:   10   continue

 71: !  Set vector values.  Note that we set multiple entries at once.
 72: !  Of course, usually one would create a work array that is the
 73: !  natural size for a particular problem (not one that is as long
 74: !  as the full vector).

 76:        call VecSetValues(x,isix,loc,xwork,INSERT_VALUES,ierr)

 78: !  Assemble vector

 80:        call VecAssemblyBegin(x,ierr)
 81:        call VecAssemblyEnd(x,ierr)

 83: !  View vector
 84:        call PetscObjectSetName(x, "initial vector:",ierr)
 85:        call VecView(x,PETSC_VIEWER_STDOUT_SELF,ierr)
 86:        call VecCopy(x,y,ierr)

 88: !  Get a pointer to vector data.
 89: !    - For default PETSc vectors, VecGetArray() returns a pointer to
 90: !      the data array.  Otherwise, the routine is implementation dependent.
 91: !    - You MUST call VecRestoreArray() when you no longer need access to
 92: !      the array.
 93: !    - Note that the Fortran interface to VecGetArray() differs from the
 94: !      C version.  See the users manual for details.

 96:        call VecGetArray(x,xx_v,xx_i,ierr)
 97:        call VecGetArray(y,yy_v,yy_i,ierr)

 99: !  Modify vector data

101:        do 30 i=1,n
102:           xx_a(i) = 100.0*i
103:           yy_a(i) = 1000.0*i
104:   30   continue

106: !  Restore vectors

108:        call VecRestoreArray(x,xx_v,xx_i,ierr)
109:        call VecRestoreArray(y,yy_v,yy_i,ierr)

111: !  View vectors
112:        call PetscObjectSetName(x, "new vector 1:",ierr)
113:        call VecView(x,PETSC_VIEWER_STDOUT_SELF,ierr)

115:        call PetscObjectSetName(y, "new vector 2:",ierr)
116:        call VecView(y,PETSC_VIEWER_STDOUT_SELF,ierr)

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

121:        call VecDestroy(x,ierr)
122:        call VecDestroy(y,ierr)
123:        call PetscFinalize(ierr)
124:        end