Actual source code: ex11f90.F90

  1: program main
  2: #include <petsc/finclude/petscvec.h>
  3:   use petscvec
  4:   implicit none

  6:   Vec        ::   x
  7:   PetscReal  :: norm
  8:   PetscMPIInt :: rank
  9:   PetscInt, parameter :: n = 20
 10:   PetscErrorCode :: ierr
 11:   PetscScalar, parameter :: sone = 1.0
 12:   PetscBool :: flg
 13:   character(len=PETSC_MAX_PATH_LEN) :: outputString
 14:   PetscInt, parameter :: zero = 0, one = 1, two = 2

 16:   PetscCallA(PetscInitialize(ierr))

 18:   PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))

 20:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-n', n, flg, ierr))

 22:   !Create a vector, specifying only its global dimension.
 23:   !When using VecCreate(), VecSetSizes() and VecSetFromOptions(),
 24:   !the vector format (currently parallel,
 25:   !shared, or sequential) is determined at runtime.  Also, the parallel
 26:   !partitioning of the vector is determined by PETSc at runtime.

 28:   !Routines for creating particular vector types directly are:
 29:   !VecCreateSeq() - uniprocessor vector
 30:   !VecCreateMPI() - distributed vector, where the user can
 31:   !determine the parallel partitioning
 32:   !VecCreateShared() - parallel vector that uses shared memory
 33:   !(available only on the SGI) otherwise,
 34:   !is the same as VecCreateMPI()

 36:   !With VecCreate(), VecSetSizes() and VecSetFromOptions() the option
 37:   !-vec_type mpi or -vec_type shared causes the
 38:   !particular type of vector to be formed.

 40:   PetscCallA(VecCreate(PETSC_COMM_WORLD, x, ierr))

 42:   PetscCallA(VecSetSizes(x, PETSC_DECIDE, n, ierr))
 43:   !
 44:   PetscCallA(VecSetBlockSize(x, two, ierr))
 45:   PetscCallA(VecSetFromOptions(x, ierr))

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

 49:   PetscCallA(VecSet(x, sone, ierr))

 51:   PetscCallA(VecNorm(x, NORM_2, norm, ierr))
 52:   write (outputString, *) norm
 53:   PetscCallA(PetscPrintf(PETSC_COMM_WORLD, 'L_2 Norm of entire vector: '//trim(outputString)//'\n', ierr))

 55:   PetscCallA(VecNorm(x, NORM_1, norm, ierr))
 56:   write (outputString, *) norm
 57:   PetscCallA(PetscPrintf(PETSC_COMM_WORLD, 'L_1 Norm of entire vector: '//trim(outputString)//'\n', ierr))

 59:   PetscCallA(VecNorm(x, NORM_INFINITY, norm, ierr))
 60:   write (outputString, *) norm
 61:   PetscCallA(PetscPrintf(PETSC_COMM_WORLD, 'L_inf Norm of entire vector: '//trim(outputString)//'\n', ierr))

 63:   PetscCallA(VecStrideNorm(x, zero, NORM_2, norm, ierr))
 64:   write (outputString, *) norm
 65:   PetscCallA(PetscPrintf(PETSC_COMM_WORLD, 'L_2 Norm of sub-vector 0: '//trim(outputString)//'\n', ierr))

 67:   PetscCallA(VecStrideNorm(x, zero, NORM_1, norm, ierr))
 68:   write (outputString, *) norm
 69:   PetscCallA(PetscPrintf(PETSC_COMM_WORLD, 'L_1 Norm of sub-vector 0: '//trim(outputString)//'\n', ierr))

 71:   PetscCallA(VecStrideNorm(x, zero, NORM_INFINITY, norm, ierr))
 72:   write (outputString, *) norm
 73:   PetscCallA(PetscPrintf(PETSC_COMM_WORLD, 'L_inf Norm of sub-vector 0: '//trim(outputString)//'\n', ierr))

 75:   PetscCallA(VecStrideNorm(x, one, NORM_2, norm, ierr))
 76:   write (outputString, *) norm
 77:   PetscCallA(PetscPrintf(PETSC_COMM_WORLD, 'L_2 Norm of sub-vector 1: '//trim(outputString)//'\n', ierr))

 79:   PetscCallA(VecStrideNorm(x, one, NORM_1, norm, ierr))
 80:   write (outputString, *) norm
 81:   PetscCallA(PetscPrintf(PETSC_COMM_WORLD, 'L_1 Norm of sub-vector 1: '//trim(outputString)//'\n', ierr))

 83:   PetscCallA(VecStrideNorm(x, one, NORM_INFINITY, norm, ierr))
 84:   write (outputString, *) norm
 85:   PetscCallA(PetscPrintf(PETSC_COMM_WORLD, 'L_inf Norm of sub-vector 1: '//trim(outputString)//'\n', ierr))

 87:   !Free work space.  All PETSc objects should be destroyed when they
 88:   !are no longer needed.
 89:   PetscCallA(VecDestroy(x, ierr))
 90:   PetscCallA(PetscFinalize(ierr))

 92: end program

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