package scientific.C4LIB

import std.math.*
// import std.util.*

@C
foreign func atan2(y: Float64, x: Float64): Float64

class FloatComplex {
  var real: Float64 = 0.0
  var imag: Float64 = 0.0

  init (a: Float64, b: Float64) {
    real = a
    imag = b
  }
  operator func -(): FloatComplex {
    FloatComplex(-real, -imag)
  }
  operator func +(right: FloatComplex): FloatComplex {
    FloatComplex(this.real + right.real, this.imag + right.imag)
  }
  operator func -(right: FloatComplex): FloatComplex {
    FloatComplex(this.real - right.real, this.imag - right.imag)
  }
  operator func *(right: FloatComplex): FloatComplex {
    FloatComplex(this.real * right.real - this.imag * right.imag, this.real * right.imag + this.imag * right.real)
  }
  operator func *(right: Float64): FloatComplex {
    FloatComplex(this.real * right, this.imag * right)
  }
  operator func /(right: FloatComplex): FloatComplex {
    var d_real = (this.real * right.real + this.imag * right.imag) / (right.real * right.real + right.imag * right.imag)
    var d_imag = (this.imag * right.real - this.real * right.imag) / (right.real * right.real + right.imag * right.imag)
    FloatComplex(d_real, d_imag)
  }
  operator func /(right: Float64): FloatComplex {
    FloatComplex(this.real / right, this.imag / right)
  }

}

func creal(c: FloatComplex): Float64 {
  return c.real
}

func cimag(c: FloatComplex): Float64 {
  return c.imag
}

func c4_abs(c: FloatComplex): Float64 {

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

    C4_ABS returns the absolute value of a C4.

  Parameters:

    Input, float complex C, the argument.

    Output, float C4_ABS, the function value.
*/

  var value: Float64
  value = sqrt ( pow ( creal ( c ), Float64(2) ) + pow ( cimag ( c ), Float64(2) ) )
  
  return value
}

func c4_acos(c1: FloatComplex): FloatComplex {


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

    C4_ACOS evaluates the inverse cosine of a C4.

  Discussion:

    Here we use the relationship:

      C4_ACOS ( Z ) = pi/2 - C4_ASIN ( Z ).

  Parameters:

    Input, float complex C1, the argument.

    Output, float complex C4_ACOS, the function value.
*/

  var c2: FloatComplex
  var c3_imag: Float64
  var c3_real: Float64
  var pi_half: Float64 = 1.5707963

  c2 = c4_asin ( c1 )

  c3_real = pi_half - creal ( c2 )
  c3_imag = - cimag ( c2 )

  var c3 = FloatComplex(c3_real, c3_imag)
  return c3
}

func c4_acosh(c1: FloatComplex): FloatComplex {


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

    C4_ACOSH evaluates the inverse hyperbolic cosine of a C4.

  Discussion:

    Here we use the relationship:

      C4_ACOSH ( Z ) = i * C4_ACOS ( Z ).

  Parameters:

    Input, float complex C1, the argument.

    Output, float complex C4_ACOSH, the function value.
*/
  var c2: FloatComplex
  c2 = c4_I ( ) * c4_acos ( c1 )

  return c2
}

func c4_asin(c1: FloatComplex): FloatComplex {


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

    C4_ASIN evaluates the inverse sine of a C4.

  Discussion:

    Here we use the relationship:

      C4_ASIN ( Z ) = - i * log ( i * z + sqrt ( 1 - z * z ) )


  Parameters:

    Input, float complex C1, the argument.

    Output, float complex C4_ASIN, the function value.
*/

  var c2: FloatComplex
  var c3: FloatComplex
  var c4: FloatComplex
  var ce: FloatComplex

  c2 = c4_I()
  c3 = c4_sqrt ( c4_one() - c1 * c1 )
  c4 = c4_log ( c3 + c2 * c1 )
  ce = - c2 * c4

  return ce
}

func c4_asinh(c1: FloatComplex): FloatComplex {

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

    C4_ASINH evaluates the inverse hyperbolic sine of a C4.

  Discussion:

    Here we use the relationship:

      C4_ASINH ( Z ) = - i * C4_ASIN ( i * Z ).

  Parameters:

    Input, float complex C1, the argument.

    Output, float complex C4_ASINH, the function value.
*/

  var c2: FloatComplex
  var c3: FloatComplex
  var c4: FloatComplex
  var c5: FloatComplex
  var c6: FloatComplex

  c2 = c4_I ( )
  c3 = c2 * c1
  c4 = c4_asin ( c3 )
  c5 = c2 * c4
  c6 = - c5

  return c6
}

func c4_atan(c1: FloatComplex): FloatComplex {


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

    C4_ATAN evaluates the inverse tangent of a C4.

  Discussion:

    Here we use the relationship:

      C4_ATAN ( Z ) = ( i / 2 ) * log ( ( 1 - i * z ) / ( 1 + i * z ) )

  Parameters:

    Input, float complex C1, the argument.

    Output, float complex C4_ATAN, the function value.
*/

  var c2: FloatComplex
  var c3: FloatComplex
  var c4: FloatComplex
  var c5: FloatComplex
  var c6: FloatComplex
  var c7: FloatComplex
  var c8: FloatComplex
  var c9: FloatComplex
  var cx: FloatComplex
  c2 = c4_I ( )
  c3 = c4_one ( )
  c4 = c2 * c1
  c5 = c3 - c4 
  c6 = c3 + c4 
  c7 = c5 / c6 

  c8 = c4_log ( c7 )
  c9 =  c2 * c8
  cx = c9 / 2.0

  return cx
}

func c4_atanh(c1: FloatComplex): FloatComplex {


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

    C4_ATANH evaluates the inverse hyperbolic tangent of a C4.

  Discussion:

    Here we use the relationship:

      C4_ATANH ( Z ) = - i * C4_ATAN ( i * Z ).

  Parameters:

    Input, float complex C1, the argument.

    Output, float complex C4_ATANH, the function value.
*/

  var c2: FloatComplex
  var c3: FloatComplex
  var c4: FloatComplex
  var c5: FloatComplex
  var c6: FloatComplex

  c2 = c4_I ( )

  c3 = c2 * c1
  c4 = c4_atan ( c3 )
  c5 = c2 * c4
  c6 = - c5

  return c6;
}
func c4_I(): FloatComplex {

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

    C4_I returns the value of I as a C4

  Parameters:

    Output, float complex C4_I, the value of I.
*/

  var c_I = FloatComplex(0.0, 1.0)
  return c_I
}

func c4_one(): FloatComplex {

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

    C4_ONE returns the value of 1 as a C4

  Parameters:

    Output, float complex C4_ONE, the value of 1.
*/

  var c_one = FloatComplex(1.0, 0.0)
  return c_one
}

func polar_to_c4(r: Float64, theta:Float64): FloatComplex {

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

    POLAR_TO_C4 converts a polar form to a C4.

  Parameters:

    Input, float R, THETA, the polar form.

    Output, float complex POLAR_TO_C4, the complex number.
*/

  var c = FloatComplex(r * cos ( theta ), r * sin ( theta ))
  return c
}

func c4_log(c1: FloatComplex): FloatComplex {

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

    C4_LOG evaluates the logarithm of a C4.

  Discussion:

    Here we use the relationship:

      C4_LOG ( Z ) = LOG ( MAG ( Z ) ) + i * ARG ( Z )

  Parameters:

    Input, float complex C1, the argument.

    Output, float complex C4_ATAN, the function value.
*/

  var arg: Float64
  var mag: Float64

  arg = c4_arg ( c1 )
  mag = c4_mag ( c1 )

  var c2 = FloatComplex(log ( mag ),  arg)

  return c2
}

func c4_mag(c: FloatComplex): Float64 {

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

    C4_MAG returns the magnitude of a C4.

  Parameters:

    Input, float complex C, the argument.

    Output, float C4_MAG, the function value.
*/

	var value: Float64

	value = sqrt ( pow ( creal ( c ), 2.0 ) + pow ( cimag ( c ), 2.0 ) )

	return value
}

func c4_sqrt(c: FloatComplex): FloatComplex {

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

    C4_SQRT computes a square root of a C4.

  Parameters:

    Input, float complex C1, the argument.

    Output, float complex C4_SQRT, the function value.
*/

  var c2: FloatComplex
  var r: Float64 = c4_abs( c )
  var t: Float64 = c4_arg( c )

  r = sqrt ( r )
  t = t / 2.0
  
  c2 = polar_to_c4 ( r, t )

  return c2
}

func c4_arg(c: FloatComplex): Float64 {


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

    C4_ARG returns the argument of a C4.

  Parameters:

    Input, float complex C, the complex number.

    Output, float C4_ARG, the function value.
*/

  var value: Float64

  if ( cimag ( c ) == 0.0 && creal ( c ) == 0.0 ) {
    value = 0.0
  } else {
    value = unsafe{ atan2 ( cimag ( c ), creal ( c ) ) }
  }
  return value
}


public func testC4lib() {
  let c1 = FloatComplex(1.0, 2.0)
  let c2 = FloatComplex(2.0, -1.0)

  var c: FloatComplex
  print("c1 = creal: ${ creal(c1) } cimag: ${ cimag(c1) }\n")
  print("c2 = creal: ${ creal(c2) } cimag: ${ cimag(c2) }\n")
  //c1 + c2
  c = c1 + c2
  print("c1 + c2 = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")
  
  //c1 - c2
  c = c1 - c2
  print("c1 - c2 = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")

  //c1 * c2
  c = c1 * c2
  print("c1 * c2 = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")

  //c1 / c2
  c = c1 / c2
  print("c1 / c2 = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")

  //c4_abs function
  print("abs(c1) = ${ c4_abs(c1) }\n")

  //c4_arg function
  print("arg(c1) = ${ c4_arg(c1) }\n")

  //c4_sqrt function
  c = c4_sqrt(c1)
  print("sqrt(c1) = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")

  //c4_log function
  c = c4_log(c1)
  print("log(c1) = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")

  //c4_acos function
  c = c4_acos(c1)
  print("acos(c1) = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")

  //c4_acosh function
  c = c4_acosh(c1)
  print("acosh(c1) = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")

  //c4_asin function
  c = c4_asin(c1)
  print("asin(c1) = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")

  //c4_asinh function
  c = c4_asinh(c1)
  print("asinh(c1) = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")

  //c4_atan function
  c = c4_atan(c1)
  print("atan(c1) = creal: ${unsafe { creal(c) } }  ciamg: ${unsafe { cimag(c) }}\n")

  //c4_atanh function
  c = c4_atanh(c1)
  print("atanh(c1) = creal: ${ creal(c) }  ciamg: ${ cimag(c) }\n")
}