//---------------------------------------------------------------------------
// ŠNikolai V. Shokhirev, 2004-2008  <nikolai@shokhirev.com> http://www.shokhirev.com/nikolai.html
// Variant 1 - "Classical"
//---------------------------------------------------------------------------
#include "complex.h"
#include "MathUtils.h"

//---------------------------------------------------------------------------

Complex::Complex(const Complex& c) {Re = c.Re; Im = c.Im; }

Complex::Complex( real r, real i) { Re = r; Im = i; }

//const Complex& Complex::operator=(const Complex& c)
Complex& Complex::operator=(const Complex& c)
{
  if(this != &c) {Re = c.Re; Im = c.Im;}
  return *this;
}

//const Complex& Complex::operator=(const real r)
Complex& Complex::operator=(const real r)
{
  Re = r; Im = 0.0;
  return *this;
}

Complex& Complex::operator+=(const Complex& c)
{
  Re += c.Re;  Im += c.Im;
  return *this;
}

Complex& Complex::operator-=(const Complex& c)
{
  Re -= c.Re;  Im -= c.Im;
  return *this;
}

Complex& Complex::operator*=(const Complex& c)
{
  Re *= c.Re;  Im *= c.Im;
  return *this;
}

Complex& Complex::operator/=(const Complex& c)
{
  Re /= c.Re;  Im /= c.Im;
  return *this;
}

/*  Can be used anywhere because of automatic conversion friend
  Complex operator +( const Complex &, const Complex & );

  Complex operator +(const Complex &c1, const Complex &c2)
  {
    real r = c1.Re + c2.Re, i = c1.Im + c2.Im;
    return Complex(r, i);
  } */

Complex Complex::operator+ (const Complex& c) { return Complex(Re + c.Re, Im + c.Im); }

Complex Complex::operator- (const Complex& c) { return Complex(Re - c.Re, Im - c.Im); }

Complex Complex::operator* (const Complex& c) { return Complex(Re*c.Re - Im*c.Im, Re*c.Im + Im*c.Re); }

Complex Complex::operator/ (const Complex& c)
{
  real cc = mod2(c);
  return Complex((Re*c.Re + Im*c.Im)/cc, (Im*c.Re - Re*c.Im)/cc);
}

//scalar math

Complex Complex::operator+ (real r) { return Complex(Re + r, Im); }

Complex Complex::operator- (real r) { return Complex(Re - r, Im); }

Complex Complex::operator* (real r) { return Complex(Re*r, Im*r); }

Complex Complex::operator/ (real r) { return Complex(Re/r , Im/r); }

//scalar math where scalars come first

Complex operator+ (real r, const Complex& c) { return Complex(r+c.Re, c.Im); }


Complex operator- (real r, const Complex& c) { return Complex(r-c.Re, -c.Im); }


Complex operator* (real r, const Complex& c) { return Complex(r*c.Re, r*c.Im); }


Complex operator/ (real r, const Complex& c)
{
  real cc = mod2(c);
  return Complex(r*c.Re/cc, -r*c.Im/cc);
}

// Complex library


Complex conjug(Complex c)
{
  return Complex(c.Re,-c.Im);
}

real mod2(const Complex& c) { return (c.Re*c.Re + c.Im*c.Im); }


real mod(const Complex& c) { return sqrt(mod2(c)); }


bool Equal(Complex c1, Complex c2, real eps)

{

  return (fabs(c1.Re-c2.Re) < eps && fabs(c1.Im-c2.Im) < eps);

}


bool Equal(Complex c1, real r, real i, real eps)

{

  return (fabs(c1.Re-r) < eps && fabs(c1.Im-i) < eps);

}


Complex exp(const Complex& c)

{
  return (exp(c.Re)*Complex(cos(c.Im),sin(c.Im)));
}

// result = sqrt(z), result.Re > 0
Complex Csqrt(Complex z)
{
  Complex result;
  result.Re = sqrt(abs(0.5*(mod(z) + z.Re)));
  result.Im = sqrt(abs(0.5*(mod(z) - z.Re)));
  return result;
};