#! /usr/bin/env python
#
def drot ( n, x, incx, y, incy, c, s ):

#*****************************************************************************80
#
## DROT applies a plane rotation.
#
#  Discussion:
#
#    Note that the FORTRAN version of this function allowed users to pass in
#    X and Y data that was noncontiguous, (such as rows of a FORTRAN matrix).
#    The rotated data overwrote the input data, and so it might therefore
#    also be noncontiguous.
#
#    This function does not assume that the output overwrites the input,
#    and treats the output vectors as new items of length exactly N.
#
#    Note, moreover, that Python does NOT allow a matrix to be treated as a 
#    vector in quite the simple way that FORTRAN does.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    06 September 2016
#
#  Author:
#
#    Python version by John Burkardt.
#
#  Reference:
#
#    Jack Dongarra, Cleve Moler, Jim Bunch and Pete Stewart,
#    LINPACK User's Guide,
#    SIAM, 1979.
#
#    Charles Lawson, Richard Hanson, David Kincaid, Fred Krogh,
#    Basic Linear Algebra Subprograms for Fortran Usage,
#    Algorithm 539,
#    ACM Transactions on Mathematical Software,
#    Volume 5, Number 3, September 1979, pages 308-323.
#
#  Parameters:
#
#    Input, integer N, the number of entries in the vectors.
#
#    Input, real X(INCX*N), one of the vectors to be rotated.
#
#    Input, integer INCX, the increment between successive entries of X.
#
#    Input, real Y(INCX*N), one of the vectors to be rotated.
#
#    Input, integer INCY, the increment between successive elements of Y.
#
#    Input, real C, S, parameters (presumably the cosine and
#    sine of some angle) that define a plane rotation.
#
#    Output, real XR(N), the rotated vector.
#
#    Output, real YR(N), the rotated vector.
#
  import numpy as np

  xr = np.zeros ( n )
  yr = np.zeros ( n )

  if ( n <= 0 ):

    pass

  elif ( incx == 1 and incy == 1 ):

    xr[0:n] = c * x[0:n] + s * y[0:n]
    yr[0:n] = c * y[0:n] - s * x[0:n]

  else:

    if ( 0 <= incx ):
      ix = 0
    else:
      ix = ( - n + 1 ) * incx

    if ( 0 <= incy ):
      iy = 0
    else:
      iy = ( - n + 1 ) * incy

    for i in range ( 0, n ):
      xr[ix] = c * x[ix] + s * y[iy]
      yr[iy] = c * y[iy] - s * x[ix]
      ix = ix + incx
      iy = iy + incy

  return xr, yr

def drot_test ( ):

#*****************************************************************************80
#
## DROT_TEST tests DROT.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    06 September 2016
#
#  Author:
#
#    John Burkardt
#
  import numpy as np
  import platform

  n = 6
  x = np.zeros ( n )
  y = np.zeros ( n )

  for i in range ( 0, n ):
    x[i] = float ( i + 1 )
    y[i] = x[i] * x[i] - 12.0
 
  print ( '' )
  print ( 'DROT_TEST' )
  print ( '  Python version: %s' % ( platform.python_version ( ) ) )
  print ( '  DROT carries out a Givens rotation.' )
  print ( '' )
  print ( '  Vectors X and Y' )
  print ( '' )
  for i in range ( 0, n ):
    print ( '  %6d  %14.6g  %14.6g' % ( i, x[i], y[i] ) )

  c = 0.5
  s = np.sqrt ( 1.0 - c * c )
  xr, yr = drot ( n, x, 1, y, 1, c, s )
  print ( '' )
  print ( '  xr, yr = drot ( n, x, 1, y, 1, %g, %g )' % ( c, s ) )
  print ( '' )
  print ( '  Rotated vectors XR and YR' )
  print ( '' )
  for i in range ( 0, n ):
    print ( '  %6d  %14.6g  %14.6g' % ( i, xr[i], yr[i] ) )

  t = np.arctan2 ( y[0], x[0] )
  c = np.cos ( t )
  s = np.sin ( t )
  xr, yr = drot ( n, x, 1, y, 1, c, s )
  print ( '' )
  print ( '  xr, yr = drot ( n, x, 1, y, 1, %g, %g )' % ( c, s ) )
  print ( '' )
  print ( '  Rotated vectors XR and YR' )
  print ( '' )
  for i in range ( 0, n ):
    print ( '  %6d  %14.6g  %14.6g' % ( i, xr[i], yr[i] ) )
#
#  Terminate.
#
  print ( '' )
  print ( 'DROT_TEST' )
  print ( '  Normal end of execution.' )
  return

if ( __name__ == '__main__' ):
  from timestamp import timestamp
  timestamp ( )
  drot_test ( )
  timestamp ( )