VectorNorm.F

c**********************************************************************
#include "author.inc"
c*    $Id: VectorNorm.F,v 1.8 1995/07/30 17:49:33 turner Exp $
c*
c*    Computes the norm of x.
c*
c*    WARNING: Note that JT_VectorNorm *must* be declared real in
c*             routines that use it.
c*
c*    <PARAMETER LIST>
c*     Input:
c*      norm - determines which vector norm to use
c*         0  ==>  infinity norm
c*         1  ==>  1-norm
c*         2  ==>  2-norm
c*      n - number of elements
c*      x - vector
c*
c*           infinity norm  ==>   max   | x |
c*                              1<=i<=n    i
c*
c*                        ---
c*                        \
c*           1-norm  ==>  /   | x |
c*                        ---    i
c*                         n
c*
c*                             ---
c*                             \    2          T
c*           2-norm  ==>  sqrt /   x  = sqrt (x x)
c*                             ---  i
c*                              n
c*
c*     Output:
c*      JT_VectorNorm - norm of x
c*      status - return status
c*
#include "copyright.inc"
c**********************************************************************
      real function JT_VectorNorm(norm, n, x, status)
      implicit none
c
c ... Input:
      integer norm, n
      real x(n)
c
c ... Output:
      integer status
c
c ... Local:
      integer i
      real zero
#ifdef use_blas
      integer index, ISAMAX
      real SASUM, SNRM2
#endif
      parameter (zero=0.0d0)
c
      JT_VectorNorm = zero
#ifdef use_blas
      if (norm .eq. 0) then
       index = ISAMAX(n, x, 1)
       JT_VectorNorm = ABS(x(index))
      elseif (norm .eq. 1) then
       JT_VectorNorm = SASUM(n, x, 1)
      else
       JT_VectorNorm = SNRM2(n, x, 1)
      endif
#else
      if (norm .eq. 0) then
       do i=1,n
        JT_VectorNorm = MAX( JT_VectorNorm , ABS(x(i)) )
       enddo
      elseif (norm .eq. 1) then
       do i=1,n
        JT_VectorNorm = JT_VectorNorm + ABS(x(i))
       enddo
      else
       do i=1,n
        JT_VectorNorm = JT_VectorNorm + x(i)*x(i)
       enddo
       JT_VectorNorm = SQRT(JT_VectorNorm)
      endif
#endif
c
      status = 0
      return
      end