function p2 = collatz_polynomial_next ( p1 )

%*****************************************************************************80
%
%% collatz_polynomial_next() computes the next Collatz polynomial.
%
%  Discussion:
%
%    The Collatz polynomial sequence is defined for polynomials with
%    coeficients mod 2.
%
%    If P1(x) is divisible by x, then
%      P2(x) = P1(x) / x
%    Else
%      P2(x) = P1(x) * (x+1) + 1
%
%    A polynomial of degree N is defined by a vector of length N+1,
%    whose first entry is the constant term, and whose last entry
%    is the coefficient of X^N.
%
%    For this function, it is assumed that all coefficients are 
%    either 0 or 1.
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    15 September 2022
%
%  Author:
%
%    John Burkardt
%
%  Reference:
%
%    John D Cook,
%    Polynomial analog of the Collatz conjecture,
%    27 July 2021,
%    https://www.johndcook.com/blog/2021/07/27/polynomial-collatz/
%    
%  Input:
%
%    integer P1(), the current polynomial.
%
%  Output:
%
%    integer P2(), the next polynomial.
%
  if ( ~ i4vec_is_binary ( p1 ) )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'collatz_polyomial_next(): Fatal error!\n' );
    fprintf ( 1, '  Input polynomial has a least one coefficient that is not 0 or 1.\n' );
    error ( 'collatz_polyomial_next(): Fatal error!\n' );
  end

  n = polynomial_degree ( p1 );

  if ( n == 0 )
    p2 = p1;
  elseif ( p1(1) == 0 )
    p2 = zeros ( 1, n );
    p2(1:n) = p1(2:n+1);
  else
    p2 = zeros ( 1, n + 2 );
    p2(1:n+1) = p1(1:n+1);
    p2(2:n+2) = p2(2:n+2) + p1(1:n+1);
    p2(1) = p2(1) + 1;
    p2 = mod ( p2, 2 );
  end

  return
end
