1: !
  2: !      "$Id: ex4f.F,v 1.28 2001/08/07 03:02:34 balay Exp $";
  3: !
  4: !  Description:  Illustrates the use of VecSetValues() to set
  5: !  multiple values at once; demonstrates VecGetArray().
  6: !
  7: !/*T
  8: !   Concepts: vectors^assembling;
  9: !   Concepts: vectors^arrays of vectors;
 10: !   Processors: 1
 11: !T*/
 12: ! -----------------------------------------------------------------------

 14:       program main
 15:       implicit none

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

 26:  #include include/finclude/petsc.h
 27:  #include include/finclude/petscvec.h

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

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

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

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

 56: !  Create initial vector and duplicate it

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

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

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

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

 75:        call VecSetValues(x,6,loc,xwork,INSERT_VALUES,ierr)

 77: !  Assemble vector

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

 82: !  View vector

 84:        write(6,20)
 85:   20   format('initial vector:')
 86:        call VecView(x,PETSC_VIEWER_STDOUT_SELF,ierr)
 87:        call VecCopy(x,y,ierr)

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

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

100: !  Modify vector data

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

107: !  Restore vectors

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

112: !  View vectors

114:        write(6,40)
115:   40   format('new vector 1:')
116:        call VecView(x,PETSC_VIEWER_STDOUT_SELF,ierr)

118:        write(6,50)
119:   50   format('new vector 2:')
120:        call VecView(y,PETSC_VIEWER_STDOUT_SELF,ierr)

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

125:        call VecDestroy(x,ierr)
126:        call VecDestroy(y,ierr)
127:        call PetscFinalize(ierr)
128:        end
129: