1: !
  2: !
  3: !  Description:  Illustrates the use of VecSetValues() to set
  4: !  multiple values at once; demonstrates VecGetArray().
  5: !
  6: ! -----------------------------------------------------------------------

  8: program main
  9: #include <petsc/finclude/petscvec.h>
 10:   use petscvec
 11:   implicit none

 13: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 14: !                 Beginning of program
 15: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 17:   PetscInt, parameter :: n = 6
 18:   PetscScalar, dimension(n) ::  xwork
 19:   PetscScalar, pointer, dimension(:) ::  xx_v, yy_v
 20:   PetscInt, dimension(n) :: loc
 21:   PetscInt i
 22:   PetscErrorCode ierr
 23:   Vec x, y

 25:   PetscCallA(PetscInitialize(ierr))

 27: !  Create initial vector and duplicate it

 29:   PetscCallA(VecCreateSeq(PETSC_COMM_SELF, n, x, ierr))
 30:   PetscCallA(VecDuplicate(x, y, ierr))

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

 36:   do 10 i = 1, n
 37:     loc(i) = i - 1
 38:     xwork(i) = 10.0*real(i)
 39: 10  continue

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

 46:     PetscCallA(VecSetValues(x, n, loc, xwork, INSERT_VALUES, ierr))

 48: !  Assemble vector

 50:     PetscCallA(VecAssemblyBegin(x, ierr))
 51:     PetscCallA(VecAssemblyEnd(x, ierr))

 53: !  View vector
 54:     PetscCallA(PetscObjectSetName(x, 'initial vector:', ierr))
 55:     PetscCallA(VecView(x, PETSC_VIEWER_STDOUT_SELF, ierr))
 56:     PetscCallA(VecCopy(x, y, ierr))

 58: !  Get a pointer to vector data.
 59: !    - For default PETSc vectors, VecGetArray() returns a pointer to
 60: !      the data array.  Otherwise, the routine is implementation dependent.
 61: !    - You MUST call VecRestoreArray() when you no longer need access to
 62: !      the array.

 64:     PetscCallA(VecGetArray(x, xx_v, ierr))
 65:     PetscCallA(VecGetArray(y, yy_v, ierr))

 67: !  Modify vector data

 69:     do 30 i = 1, n
 70:       xx_v(i) = 100.0*real(i)
 71:       yy_v(i) = 1000.0*real(i)
 72: 30    continue

 74: !  Restore vectors

 76:       PetscCallA(VecRestoreArray(x, xx_v, ierr))
 77:       PetscCallA(VecRestoreArray(y, yy_v, ierr))

 79: !  View vectors
 80:       PetscCallA(PetscObjectSetName(x, 'new vector 1:', ierr))
 81:       PetscCallA(VecView(x, PETSC_VIEWER_STDOUT_SELF, ierr))

 83:       PetscCallA(PetscObjectSetName(y, 'new vector 2:', ierr))
 84:       PetscCallA(VecView(y, PETSC_VIEWER_STDOUT_SELF, ierr))

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

 89:       PetscCallA(VecDestroy(x, ierr))
 90:       PetscCallA(VecDestroy(y, ierr))
 91:       PetscCallA(PetscFinalize(ierr))
 92:     end

 94: !
 95: !/*TEST
 96: !
 97: !     test:
 98: !
 99: !TEST*/