# include <complex.h>
# include <math.h>
# include <stdbool.h>
# include <stdio.h>
# include <stdlib.h>
# include <string.h>

# include "c8poly.h"

/******************************************************************************/

void c8poly_print ( int n, double complex c[], char *title )

/******************************************************************************/
/*
  Purpose:

    c8poly_print() prints out a polynomial.

  Licensing:

    This code is distributed under the MIT license.

  Modified:

    07 December 2023

  Author:

    John Burkardt

  Input:

    int N, the dimension of C.

    double complex C[N+1], the polynomial coefficients.
    C(0) is the constant term and
    C(N) is the coefficient of X^N.

    char *TITLE, a title.
*/
{
  int i;

  printf ( "\n" );

  if ( n < 0 )
  {
    printf ( "  p(z) = 0\n" );
    return;
  }

  if ( 2 <= n )
  {
    printf ( "  p(z) = (%g,%g) * z^%d\n", creal ( c[n] ), cimag ( c[n] ), n );
  }
  else if ( n == 1 )
  {
    printf ( "  p(z) = (%g,%g) * z\n", creal ( c[n] ), cimag ( c[n] ) );
  }
  else if ( n == 0 )
  {
    printf ( "  p(z) = (%g,%g)\n", creal ( c[n] ), cimag ( c[n] ) );
  }

  for ( i = n - 1; 0 <= i; i-- )
  {
    if ( 2 <= i )
    {
      printf ( "         (%g,%g) * z^%d\n", creal ( c[i] ), cimag ( c[i] ), i );
    }
    else if ( i == 1 )
    {
      printf ( "         (%g,%g) * z\n", creal ( c[i] ), cimag ( c[i] ) );
    }
    else if ( i == 0 )
    {
      printf ( "         (%g,%g)\n", creal ( c[i] ), cimag ( c[i] ) );
    }
  }
  return;
}
/******************************************************************************/

double complex *roots_to_c8poly ( int n, double complex r[] )

/******************************************************************************/
/*
  Purpose:

    roots_to_c8poly() converts polynomial roots to polynomial coefficients.

  Licensing:

    This code is distributed under the MIT license.

  Modified:

    07 December 2023

  Author:

    John Burkardt

  Input:

    int N, the number of roots.

    double complex R[N], the roots.

  Output:

    double complex ROOTS_TO_C8POLY[N+1], the coefficients of the polynomial.
*/
{
  double complex *c;
  int i;
  int j;

  c = ( double complex * ) malloc ( ( n + 1 ) * sizeof ( double complex ) );
/*
  Initialize C to (0, 0, ..., 0, 1).
  Essentially, we are setting up a divided difference table.
*/
  for ( i = 0; i < n; i++ )
  {
    c[i] = 0.0;
  }
  c[n] = 1.0;
/*
  Convert to standard polynomial form by shifting the abscissas
  of the divided difference table to 0.
*/
  for ( j = 1; j <= n; j++ )
  {
    for ( i = 1; i <= n + 1 - j; i++ )
    {
      c[n-i] = c[n-i] - r[n+1-i-j] * c[n-i+1];
    }
  }
  return c;
}