#! /usr/bin/env python3
#
def fermat_boeck_test ( ):

#*****************************************************************************80
#
## fermat_boeck_test() tests fermat_boeck().
#
#  Licensing:
#
#    This code is distributed under the MIT license. 
#
#  Modified:
#
#    11 April 2025
#
#  Author:
#
#    John Burkardt
#
  import numpy as np
  import platform

  print ( '' )
  print ( 'fermat_boeck_test():' )
  print ( '  python version: ' + platform.python_version ( ) )
  print ( '  numpy version:  ' + np.version.version )
  print ( '  fermat_boeck() factors an integer using the Fermat method.' )

  for n in [ 39, 25951, 2025090001, 2027651281 ]:

    pq = fermat_boeck ( n )

    if ( not pq ):
      print ( '' )
      print ( '  Could not find factors in 100 tries.' )
    else:
      p = pq[0]
      q = pq[1]
      print ( '' )
      print ( '  Seeking factors of n = ', n )
      print ( '  ', p, ' * ', q, ' = ', p * q )
#
#  Terminate.
#
  print ( '' )
  print ( 'fermat_boeck_test():' )
  print ( '  Normal end of execution.' )

  return

def fermat_boeck ( n ):

#*****************************************************************************80
#
## fermat_boeck() uses Fermat factorization on an integer.
#
#  Discussion:
#
#    The method determines values A and B such that 
#      N = A^2 - B^2  = ( A + B ) * ( A - B ).
#
#    The method is most efficient when B is small, that is, the two
#    factors are close.  For encryption, the key N is often generated
#    by two primes that are roughly the same size, and this make this
#    method complete more quickly.
#
#    The gmpy2 library allows for multiple precision arithmetic.
#
#  Licensing:
#
#    This code is distributed under the MIT license. 
#
#  Modified:
#
#    11 April 2025
#
#  Author:
#
#    John Burkardt
#
#  Reference:
#
#    Manon Bischoff,
#    This more than 380-year-old trick can crack some modern encryption,
#    Spektrum der Wissenschaft,
#    09 April 2025.
#
#    Hanno Boeck,
#    Fermat factorization in the wild,
#    08 January 2023.
#
#  Input:
#
#    integer n: the number to be factored.
#
#  Output:
#
#    integer [ p, q ]: factors of n.
#
  import gmpy2
  from gmpy2 import is_square
  from gmpy2 import isqrt

  tries = 100

  a = gmpy2.isqrt ( n )

  c = 0
  while not gmpy2.is_square ( a**2 - n ):
    a = a + 1
    c = c + 1
    if ( tries < c ):
      return False
  bsq = a**2 - n
  b = gmpy2.isqrt ( bsq )
  p = a + b
  q = a - b

  return [ p, q ]

def timestamp ( ):

#*****************************************************************************80
#
## timestamp() prints the date as a timestamp.
#
#  Licensing:
#
#    This code is distributed under the MIT license. 
#
#  Modified:
#
#    21 August 2019
#
#  Author:
#
#    John Burkardt
#
  import time

  t = time.time ( )
  print ( time.ctime ( t ) )

  return

if ( __name__ == "__main__" ):
  timestamp ( )
  fermat_boeck_test ( )
  timestamp ( )