#ifndef __COMPLEX_H_INCLUDE
#define __COMPLEX_H_INCLUDE

#include <math.h>

//----  複素数クラスの宣言

struct Complex
{
  double real, imag;
  inline Complex& operator+=( const Complex& );
  inline Complex& operator-=( const Complex& );
  inline Complex& operator*=( const double& );
  inline Complex& operator*=( const Complex& );

  inline Complex( void ){}
  inline Complex( const double& _real, const double& _imag=0.0 ) : real(_real), imag(_imag) {}
};


//----  複素数クラスのメンバ関数
//----  複素数 += 複素数の定義
inline Complex& Complex::operator+=( const Complex& c )
{
  real += c.real;    imag += c.imag;
  return *this;
}
//----  複素数 -= 複素数の定義
inline Complex& Complex::operator-=( const Complex& c )
{
  real -= c.real;    imag -= c.imag;
  return *this;
}
inline Complex& Complex::operator*=( const Complex& c )
{
  const double _real = real*c.real - imag*c.imag;
  imag = real*c.imag + imag*c.real;
  real = _real;
  return *this;  
}
inline Complex& Complex::operator*=( const double& d )
{
  real *= d;
  imag *= d;
  return *this;  
}

//----  複素数クラスの関連関数

//----  単項演算子の定義
inline Complex  operator - ( const Complex& c )
{
  Complex t;
  t.real = -c.real;  t.imag = -c.imag;
  return t;
}


//----  複素数 + 実数の定義
inline Complex  operator+( const Complex& c, const double& real )
{
  Complex t;
  t.real = c.real+real;  t.imag = c.imag;
  return t;
}

//----  実数 + 複素数の定義
inline Complex  operator+( const double& real, const Complex& c )
{
  Complex t;
  t.real = c.real+real;  t.imag = c.imag;
  return t;
}

//----  複素数 + 複素数の定義
inline Complex  operator+( const Complex& c1, const Complex& c2 )
{
  Complex t;
  t.real = c1.real+c2.real;  t.imag = c1.imag+c2.imag;
  return t;
}

//----  複素数 - 実数の定義
inline Complex  operator-( const Complex& c, const double& real )
{
  Complex t;
  t.real = c.real-real;  t.imag = c.imag;
  return t;
}

//----  実数 - 複素数の定義
inline Complex  operator-( const double& real, const Complex& c )
{
  Complex t;
  t.real = real-c.real;  t.imag = -c.imag;
  return t;
}

//----  複素数 - 複素数の定義
inline Complex  operator-( const Complex& c1, const Complex& c2 )
{
  Complex t;
  t.real = c1.real-c2.real;  t.imag = c1.imag-c2.imag;
  return t;
}

//----  実数 * 複素数の定義
inline Complex  operator*( double k, const Complex& c )
{
  Complex t;
  t.real = k*c.real;  t.imag = k*c.imag;
  return t;
}

//----  複素数 * 実数の定義
inline Complex  operator*( const Complex& c, double k )
{
  Complex t;
  t.real = k*c.real;  t.imag = k*c.imag;
  return t;
}

//----  複素数 / 実数の定義
inline Complex  operator/( const Complex& c, double k )
{
  const double inv_k = 1.0/k;
  Complex t;
  t.real = inv_k*c.real;  t.imag = inv_k*c.imag;
  return t;
}

//----  実数 / 複素数の定義
inline Complex  operator/( double k, const Complex& c )
{
  const double fact   = k/(c.real*c.real + c.imag*c.imag);
  Complex t;
  t.real = +fact*c.real;
  t.imag = -fact*c.imag;
  return t;
}

//----  複素数 * 複素数の定義
inline Complex  operator*( const Complex& c1, const Complex& c2 )
{
  Complex t;
  t.real = c1.real*c2.real - c1.imag*c2.imag;
  t.imag = c1.real*c2.imag + c1.imag*c2.real;
  return t;
}

//----  複素数 / 複素数の定義
inline Complex  operator/( const Complex& c1, const Complex& c2 )
{
  const double fact   = 1.0/(c2.real*c2.real + c2.imag*c2.imag);
  Complex t;
  t.real = (c1.real*c2.real + c1.imag*c2.imag) * fact;
  t.imag = (c1.imag*c2.real - c1.real*c2.imag) * fact;
  return t;
}

inline double real( const Complex& c )
{
  return c.real;
}
inline double imag( const Complex& c )
{
  return c.imag;
}

//----  複素共役の定義
inline Complex conj( const Complex& c )
{
  Complex t;
  t.real = +c.real;
  t.imag = -c.imag;
  return t;
}

//----  絶対値の定義
inline double abs( const Complex& c )
{
  return hypot( c.real, c.imag );
}

//----  絶対値の2乗の定義
inline double norm( const Complex& c )
{
  return c.real*c.real+c.imag*c.imag;
}

//----  偏角の定義
inline double arg( const Complex& c )
{
  return atan2( c.imag, c.real );
}

//----  指数関数の定義
inline Complex  exp( const Complex& c )
{
  const double R = exp(c.real);
  Complex t;
  t.real = R*cos(c.imag);
  t.imag = R*sin(c.imag);
  return t;
}

//----  極座標形式で指定
inline Complex  polar( const double& r, const double& d )
{
  Complex t;
  t.real = r*cos(d);
  t.imag = r*sin(d);
  return t;
}

//---- 実数の単位
const Complex C1(1.0, 0.0);

//---- 虚数の単位
const Complex CI(0.0, 1.0);


//---- square function
#ifndef __TEMPLATE_SQR
#define __TEMPLATE_SQR
template <class T>
inline T sqr( const T& c )
{
  return c*c;
}
#endif // __TEMPLATE_SQR

#endif  // __COMPLEX_H_INCLUDE