//  ============================================================
//
//  inverseBinary.h
//
//  Copyright 2002, Dennis Meilicke and Rene Peralta
//
//  ============================================================
//
//  Description:
//
//  NTTL Implementation of:
//      Inverse
//
//  While the NTTL function ExtendedGCD can also be used to 
//  calculate the modular inverse, this function can be
//  optimized to do fewer operations than the full ExtendedGCD.
//
//  This algorithm is based on the binary GCD algorithm.  I
//  remember exactly where this algorithm came from.  Try and
//  find the source...
//
//  Note that this algorithm works fastest when the two
//  arguments are about the same size.
//
//  ============================================================
//
//  To Do:
//
//  There has to be a better way to handle negative arguments.
//
//  As currently coded, this breaks the Return Value
//  Optimization.  It should be fixed so that it always returns
//  the same variable.
//
//  ============================================================
#ifndef __nttlHeader_inverseBinary__
#define __nttlHeader_inverseBinary__


#include <nttl/traits.h>


template< typename T >
T
InverseBinary( T y, T x )
//
//  Compute i = y^-1 mod x
//
//  Assumes x (modulus) is odd.
//  Works fastest when |y| ~= |x|
//
{
  nttlTraits<T> tr;

  long flag = 0;
  if( tr.IsNegative( y ) ) 
  {
    flag = 1;
    y = x-y;
  } 
  
  while( tr.IsEven( x ) && tr.IsEven( y ) )
  {
    x >>= 1;
    y >>= 1;
  }

  T u = x;
  T v = y;
  T B = 0;
  T D = 1;

  while( u != 0 )
  {
    while( tr.IsEven( u ) )
    {
      u >>= 1;
      if( tr.IsEven( B ) )
      {
        B >>= 1;
      }
      else
      {
        B -= x;  B >>= 1;
      }
    }
    while( tr.IsEven( v ) )
    {
      v >>= 1;
      if( tr.IsEven( D ) )
      {
        D >>= 1;
      }
      else
      {
        D -= x;  D >>= 1;
      }
    }

    if( u >= v )
    {
      u -=v;
      B -= D;
    }
    else
    {
      v -= u;
      D -= B;
    }
  }

  //
  //  GCD = v * 2^g, where g is the number of common
  //  factors of 2 removed in the first while loop
  //  of the algorithm.
  //
 
  if( tr.IsNegative( D ) )
    D += x;

  //  cout << "exiting  bInv " << endl;
  if (flag) return (x-D);
  return D;
}


//  ============================================================


#ifndef __nttlAlgorithm_Inverse__
#define __nttlAlgorithm_Inverse__


template< typename T >
inline
T
Inverse( const T& x, const T& modulus )
{
  return InverseBinary( x, modulus );
}


#endif  //  __nttlAlgorithm_Inverse__


//  ============================================================


#endif  //  __nttlHeader_inverseBinary__