#include "library.h"


struct complex
{ double real, imag; };

complex make_complex(double r, double i)
{ complex c;
  c.real = r;
  c.imag = i;
  return c; }

void print(complex c)
{ cout <<  c.real << " + " << c.imag << "j\n"; }

complex add(complex a, complex b)
{ complex c;
  c.real = a.real + b.real;
  c.imag = a.imag + b.imag;
  return c; }

complex sub(complex a, complex b)
{ complex c;
  c.real = a.real - b.real;
  c.imag = a.imag - b.imag;
  return c; }

complex mul(complex a, complex b)
{ complex c;
  c.real = a.real * b.real - a.imag * b.imag;
  c.imag = a.real * b.imag + a.imag * b.real;
  return c; }

complex div(complex a, complex b)
{ complex c;
  double bottom = b.real * b.real + b.imag * b.imag;
  if (bottom==0.0)
    cout << "Error, division by zero!\n";
  else
  { c.real = (a.real * b.real + a.imag * b.imag) / bottom;
    c.imag = (b.real * a.imag - a.real * b.imag) / bottom; }
  return c; }

double magnitude(complex a)
{ return sqrt(a.real*a.real + a.imag*a.imag); }

const int winw = 1000, winh = 800;
const double realmax = 2.0, imagmax = 1.5;

double xtoreal(double x)
{ return x/winw*(realmax-realmin)-realmax; }

double ytoimag(double y)
{ return y/winh*(2*imagmax)-imagmax; }

double realtox(double r)
{ return (r+realmax)/(2*realmax)*winw; }

double imagtoy(double i)
{ return (i+imagmax)/(2*imagmax)*winh; }

void move_to(complex c)
{ move_to(realtox(c.real), imagtoy(c.imag)); }

void draw_to(complex c)
{ draw_to(realtox(c.real), imagtoy(c.imag)); }

void draw_an_x()
{ move_relative(-10, -10);
  draw_relative(20, 20);
  move_relative(-20, 0);
  draw_relative(20, -20); }

void draw_axes()
{ double real0 = realtox(0);
  move_to(real0, 0.0);
  draw_to(real0, (double)winh);
  double imag0 = imagtoy(0);
  move_to(0.0, imag0);
  draw_to((double)winw, imag0); }

complex click_to_complex()
{ wait_for_mouse_click();
  double x = get_click_x();
  double y = get_click_y();
  double r = xtoreal(x);
  double i = ytoimag(y);
  complex c = make_complex(r, i);
  return c; }

void plot(complex c)
{ complex z = make_complex(0, 0);
  for (int i = 0; i < 100; i += 1)
  { z = add(mul(z, z), c);
    if (i == 0)
      move_to(z);
    else
      draw_to(z); } }

bool is_stable(complex c)
{ complex z = make_complex(0, 0);
  for (int i = 0; i < 100; i += 1)
  { z = add(mul(z, z), c);
    if (magnitude(z) >= 10)
      return false; }
  return true; }

void main() // last version: plots the Mandelbrot set
{ make_window(winw, winh);
  draw_axes();
  set_pen_width(3);
  complex prev = make_complex(0, 0);
  for (int x = 0; x < winw; x +=1)
    for (int y = 0; y < winh; y += 1)
    { complex c = make_complex(xtoreal(x), ytoimag(y));
      if (is_stable(c))
        set_pen_color(color::red);
      else
        set_pen_color(color::blue);
      move_to(c);
      draw_point(); } }

/*
void main() // third version: indicate whether we are in a stable area or not
{ make_window(winw, winh);
  draw_axes();
  set_pen_width(3);
  complex prev = make_complex(0, 0);
  while (true)
  { complex c = click_to_complex();
    if (is_stable(c))
      set_pen_color(color::red);
    else
      set_pen_color(color::blue);
    move_to(c);
    draw_point(); } }
*/

/*
void main() // second version: plot successive values of z = z*z+c
{ make_window(winw, winh);
  draw_axes();
  complex prev = make_complex(0, 0);
  while (true)
  { complex c = click_to_complex();
    set_pen_color(color::grey);
    plot(prev);
    set_pen_color(color::black);
    plot(c);
    prev = c; } }
*/



/*

void main() // first version: just a calculator
{ make_window(winw, winh);
  draw_axes();
  while (true)
  { complex c = click_to_complex();
    move_to(c); draw_an_x();
    print(c);
    complex d = click_to_complex();
    move_to(d); draw_an_x();
    print(d);
    move_to(add(c, d)); draw_an_x(); write_string("+");
    move_to(sub(c, d)); draw_an_x(); write_string("-");
    move_to(mul(c, d)); draw_an_x(); write_string("*");
    move_to(div(c, d)); draw_an_x(); write_string("/");
 } }
*/