function value = is_prime_fermat ( n, k )

%*****************************************************************************80
%
%% is_prime_fermat() applies Fermat's test for primality.
%
%  Discussion:
%
%    The test is only probabilistic.
%    If the input n is a prime, is_prime(n) will be true.
%    If n is not prime, then is_prime(n) is probably false.
%    Reliability can be increased by increasing the value of k, which
%    repeats the test with different bases.
%    However, there are some values of n which are not prime, but for
%    which is_prime(n) will always return true.
%
%  Licensing:
%
%    This code is distributed under the MIT license. 
%
%  Modified:
%
%    20 September 2024
%
%  Author:
%
%    Original C++ code on Geeks for Geeks.
%    This version by John Burkardt.
%
%  Reference:
%
%    Thomas Cormen, Charles Leiserson, Ronald Rivest, Clifford Stein,
%    Introduction to Algorithms,
%    Section 31.8: Primality testing,
%    MIT Press; McGraw-Hill, 2001.
%
%  Input:
%
%    integer n: the value to be tested.
%
%    integer k: the number of times the test is to be carried out.
%
%  Output:
%
%    boolean value: true if n passed the Fermat test.
%
  if ( n <= 1 || n == 4 )
    value = false;
    return
  end

  if ( n <= 3 )
    value = true;
    return
  end
%
%  Do the test k times.
%
  for i = 1 : k
%
%  Pick a random number in [2..n-2]        
%
    a = randi ( [ 2, n-2 ], 1, 1 ); 
% 
%  Checking that a and n are co-prime.
%
    if ( gcd ( n, a ) ~= 1 )
      value = false;
      return
    end
%
%  Check Fermat's little theorem.
%
    if ( prime_power ( a, n - 1, n ) ~= 1 )
      value = false;
      return
    end

  end
 
  value = true;

  return
end

