/*---------------------------------------------------------------------------
  Copyright 2017-2020 Microsoft Research, Daan Leijen.

  This is free software; you can redistribute it and/or modify it under the
  terms of the Apache License, Version 2.0. A copy of the License can be
  found in the LICENSE file at the root of this distribution.
---------------------------------------------------------------------------*/

/* 64-bit IEEE double floating point numbers.

*/
module double

import parse
import int32

extern import
  c file "double-inline.h"
  js file "double-inline.js"


//-----------------------------------------
// Constants
//-----------------------------------------

// π
pub val pistd/num/double/pi: double          = 0x1.9_21FB_5444_2D18p+1 //3.141592653589793238462643383279502884

// π
pub val dbl-pistd/num/double/dbl-pi: double      = pistd/num/double/pi: double

// 2π
pub val dbl-twopistd/num/double/dbl-twopi: double   = 0x1.9_21FB_5444_2D18p+2

// π/2
pub val dbl-pi2std/num/double/dbl-pi2: double     = 0x1.9_21FB_5444_2D18p+0

// π/4
pub val dbl-pi4std/num/double/dbl-pi4: double     = 0x1.9_21FB_5444_2D18p-1

// 3π/4
pub val dbl-pi34std/num/double/dbl-pi34: double    =0x1.2_D97C_7F33_21D2p+1

// The [_e_](https://en.wikipedia.org/wiki/E_(mathematical_constant)) constant.
pub val dbl-estd/num/double/dbl-e: double       = 0x1.5_BF0A_8B14_5769p+1 // 2.718281828459045235360287471352662498

// The natural logarithm of 2
pub val dbl-log2std/num/double/dbl-log2: double    = 0x1.6_2E42_FEFA_39EFp-1 // 0.693147180559945309417232121458176568

// The natural logarithm of 10
pub val dbl-log10std/num/double/dbl-log10: double   = 0x1.2_6BB1_BBB5_5516p+1 // 2.302585092994045684017991454684364208

// The base-2 logarithm of _e_.
pub val dbl-log2estd/num/double/dbl-log2e: double   =  0x1.7_1547_652B_82FEp+0 // 1.442695040888963407359924681001892137

// The base-10 logarithm of _e_.
pub val dbl-log10estd/num/double/dbl-log10e: double  = 0x1.B_CB7B_1526_E50Ep-2 // 0.434294481903251827651128918916605082

// The square-root of 2
pub val dbl-sqrt2std/num/double/dbl-sqrt2: double   = 0x1.6_A09E_667F_3BCDp+0 // 1.414213562373095048801688724209698079

// `1.0 / sqrt(2.0)`  (= `sqrt(0.5)`)
pub val dbl-sqrt12std/num/double/dbl-sqrt12: double = 0x1.6_A09E_667F_3BCDp-1 // 0.707106781186547524400844362104849039

// [Euler's constant](https://en.wikipedia.org/wiki/Euler%E2%80%93Mascheroni_constant)
pub val dbl-eulerstd/num/double/dbl-euler: double   = 0x1.2_788C_FC6F_B619p-1 // 0.577215664901532860606512090082402431

// Maximum double value
pub val dbl-maxstd/num/double/dbl-max: double     = 0x1.F_FFFF_FFFF_FFFFp+1023 // 1.79769313486231570815e+308

// Smallest positive normalized double value
pub val dbl-minstd/num/double/dbl-min: double     = 0x1.0p-1022 // 2.22507385850720138309e-308

// Smallest positive subnormal value (i.e. [``DBL_TRUE_MIN``](http://en.cppreference.com/w/cpp/types/climits))
pub val dbl-true-minstd/num/double/dbl-true-min: double= 0x1.0p-1074 // 4.9406564584124654418e-324

// Machine epsilon: the difference between 1.0 and the next representable `:double` value.
pub val dbl-epsilonstd/num/double/dbl-epsilon: double = 0x1.0p-52 // 2.2204460492503130808e-16

// maximal precision in decimal digits of a `:double`.
pub val dbl-max-precstd/num/double/dbl-max-prec: int= 15

// Represents a value that is _not a number_ (NaN)
pub val nanstd/num/double/nan: double : doublestd/core/types/double: V    = make-nanstd/num/double/make-nan: () -> double()

// Positive infinity
pub val posinfstd/num/double/posinf: double : doublestd/core/types/double: V = make-posinfstd/num/double/make-posinf: () -> double()

// Negative infinity
pub val neginfstd/num/double/neginf: double : doublestd/core/types/double: V = make-neginfstd/num/double/make-neginf: () -> double()

//-----------------------------------------
// NaN, Infinity
//-----------------------------------------

extern make-nanstd/num/double/make-nan: () -> double() : doublestd/core/types/double: V
  c  inline "(double)NAN"
  cs inline "double.NaN"
  js inline "NaN"

extern make-posinfstd/num/double/make-posinf: () -> double() : doublestd/core/types/double: V
  c  inline "HUGE_VAL"
  cs inline "double.PositiveInfinity"
  js inline "Infinity"

extern make-neginfstd/num/double/make-neginf: () -> double() : doublestd/core/types/double: V
  c  inline "-HUGE_VAL"
  cs inline "double.NegativeInfinity"
  js inline "-Infinity"


// Is this value equal to NaN ?
pub inline extern is-nanstd/num/double/is-nan: (d : double) -> bool(d:doublestd/core/types/double: V) : boolstd/core/types/bool: V
  c  inline "isnan(#1)"
  cs "double.IsNaN"
  js "isNaN"

// Is this value equal to negative or positive infinity ?
pub extern is-infstd/num/double/is-inf: (d : double) -> bool(dd: double:doublestd/core/types/double: V) : boolstd/core/types/bool: V
  c  inline "isinf(#1)"
  cs "double.IsInfinity"
  js inline "((#1) === Infinity || (#1) === -Infinity)"

// Is this value equal to positive infinity ?
pub inline extern is-posinfstd/num/double/is-posinf: (d : double) -> bool(d:doublestd/core/types/double: V) : boolstd/core/types/bool: V
  c  inline "(isinf(#1) && !signbit(#1))"
  cs "double.IsPositiveInfinity"
  js inline "((#1) === Infinity)"

// Is this value equal to negative infinity ?
pub inline extern is-neginfstd/num/double/is-neginf: (d : double) -> bool(d:doublestd/core/types/double: V) : boolstd/core/types/bool: V
  c inline "(isinf(#1) && signbit(#1))"
  cs "double.IsNegativeInfinity"
  js inline "((#1) === -Infinity)"

// Is this a finite number (i.e. not `nan` or infinity)
pub inline extern is-finitestd/num/double/is-finite: (d : double) -> bool( d : doublestd/core/types/double: V ) : boolstd/core/types/bool: V
  c  inline "isfinite(#1)"
  cs inline "(!double.IsNaN(#1) && !double.IsInfinity(#1))"
  js "isFinite"

// Is this a negative zero value?
pub fun is-negzerostd/num/double/is-negzero: (d : double) -> bool( dd: double : doublestd/core/types/double: V ) : boolstd/core/types/bool: V
  (dd: double==std/core/(==).2: (double, double) -> bool0.0 &&std/core/types/(&&): (x : bool, y : bool) -> bool is-neginfstd/num/double/is-neginf: (d : double) -> bool(1.0 /std/core/(/).1: (double, double) -> double dd: double))

// Is this a [subnormal](https://en.wikipedia.org/wiki/Denormal_number) value?
// (i.e. `0 < d.abs < dbl-min`)
pub fun is-subnormalstd/num/double/is-subnormal: (d : double) -> bool( dd: double :doublestd/core/types/double: V ) : boolstd/core/types/bool: V
  (dd: double !=std/core/(!=).2: (double, double) -> bool 0.0 &&std/core/types/(&&): (x : bool, y : bool) -> bool dd: double.absstd/core/abs.1: (d : double) -> double <std/core/(<).4: (double, double) -> bool dbl-minstd/num/double/dbl-min: double)


//-----------------------------------------
// Rounding
//-----------------------------------------

// Round a double to its nearest integral value.
// If the value is halfway between two integers, the value is rounded to the even one.
pub inline extern roundstd/num/double/round: (d : double) -> double(d:doublestd/core/types/double: V) : doublestd/core/types/double: V
  c  "round"              // assume the rounding mode is set correctly by kklib
  cs "Math.Round"
  js "$std_core._double_round"

// Return the largest integer equal or less than `d`
pub inline extern floorstd/num/double/floor: (d : double) -> double(d:doublestd/core/types/double: V) : doublestd/core/types/double: V
  c  "floor"
  cs "Math.Floor"
  js "Math.floor"

// Return the smallest integer equal or larger than `d`
pub inline extern ceilingstd/num/double/ceiling: (d : double) -> double : (d:doublestd/core/types/double: V) -> doublestd/core/types/double: V
  c  "ceil"
  cs "Math.Ceiling"
  js "Math.ceil"

// Return the integral part of a `:double` `d` .
// If `d >= 0.0` , return the largest integer equal or less to `d` ,
// If `d < 0.0` , return the smallest integer equal or larger to `d` .
pub fun truncatestd/num/double/truncate: (d : double) -> double(dd: double : doublestd/core/types/double: V) : doublestd/core/types/double: V
  if dd: double >=std/core/(>=).3: (double, double) -> bool 0.0 then floorstd/num/double/floor: (d : double) -> double(dd: double) else ceilingstd/num/double/ceiling: (d : double) -> double(dd: double)

// Return the fractional part of a `:double` `d`.\
// `d.truncate + d.fraction === d`\
// `(-2.4).fraction === -0.4`
pub fun fractionstd/num/double/fraction: (d : double) -> double(dd: double : doublestd/core/types/double: V) : doublestd/core/types/double: V
  dd: double -std/core/(-).2: (double, double) -> double dd: double.truncatestd/num/double/truncate: (d : double) -> double

// Return the 'floored' fraction of `d`, always greater or equal to zero.\
// `d.floor + d.ffraction === d`\
// `(-2.4).ffraction === 0.6`
pub fun ffractionstd/num/double/ffraction: (d : double) -> double(dd: double : doublestd/core/types/double: V) : doublestd/core/types/double: V
  dd: double -std/core/(-).2: (double, double) -> double dd: double.floorstd/num/double/floor: (d : double) -> double


// Round a double to a specified precision. Rounds to the  even number in case of a tie.\
// `123.456.round-to-prec(2) == 123.46`\
// `123.456.round-to-prec(-1) == 120.0`\
pub fun round-to-precstd/num/double/round-to-prec: (d : double, prec : int) -> double( dd: double : doublestd/core/types/double: V, precprec: int : intstd/core/types/int: V  ) : doublestd/core/types/double: V
  if precprec: int <=std/core/(<=).1: (x : int, y : int) -> bool 0 then dd: double.roundstd/num/double/round: (d : double) -> double else
    val pp: double  = exp10std/num/double/exp10: (p : double) -> double(precprec: int.doublestd/core/double: (i : int) -> double)
    (dd: double *std/core/(*).1: (double, double) -> double pp: double).roundstd/num/double/round: (d : double) -> double /std/core/(/).1: (double, double) -> double pp: double


// fused multiply-add. Computes `(x*y)+z` as if to infinite precision
// with only the final result rounded back to a `:double`.
pub extern fmaddstd/num/double/fmadd: (x : double, y : double, z : double) -> double( xx: double : doublestd/core/types/double: V, yy: double : doublestd/core/types/double: V, zz: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "fma"
  js "_fmadd"
  cs "Math.FusedMultiplyAdd"

//-----------------------------------------
// Powers
//-----------------------------------------

// Return the square root of a value `d`
// Returns `nan` if `d == nan`  or if `d`  is negative.
// Returns `inf` if `d == inf` .
pub inline extern sqrtstd/num/double/sqrt: (d : double) -> double(d:doublestd/core/types/double: V) : doublestd/core/types/double: V
  c  "sqrt"
  cs "Math.Sqrt"
  js "Math.sqrt"


// Return the `d` raised to the power of `p`.
pub inline extern powstd/num/double/pow: (d : double, p : double) -> double(d:doublestd/core/types/double: V, p:doublestd/core/types/double: V) : doublestd/core/types/double: V
  c  "pow"
  cs "Math.Pow"
  js "Math.pow"

// Return the natural logarithm (in base _e_)  of a `:double` `d`
pub inline extern logstd/num/double/log: (d : double) -> double(d:doublestd/core/types/double: V) : doublestd/core/types/double: V
  c  "log"
  cs "Math.Log"
  js "Math.log"

// Return the logarithm in base 10 of a `:double` `d`.
pub fun log10std/num/double/log10: (d : double) -> double( dd: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  logstd/num/double/log: (d : double) -> double(dd: double) /std/core/(/).1: (double, double) -> double dbl-log10std/num/double/dbl-log10: double

// Return the logarithm in base 2 of a `:double` `d`.
pub fun log2std/num/double/log2: (d : double) -> double( dd: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  logstd/num/double/log: (d : double) -> double(dd: double) /std/core/(/).1: (double, double) -> double dbl-log2std/num/double/dbl-log2: double

// Return _e_ to the power of `p`.
pub inline extern expstd/num/double/exp: (p : double) -> double( p : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "exp"
  cs "Math.Exp"
  js "Math.exp"

// Return the 10 to the power of `p`.
pub fun exp10std/num/double/exp10: (p : double) -> double( pp: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  powstd/num/double/pow: (d : double, p : double) -> double(10.0,pp: double)

// Return the 2 to the power of `p`.
pub fun exp2std/num/double/exp2: (p : double) -> double( pp: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  powstd/num/double/pow: (d : double, p : double) -> double(2.0,pp: double)

fun log2p1std/num/double/log2p1: (x : double) -> double( xx: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  dbl-log2estd/num/double/dbl-log2e: double *std/core/(*).1: (double, double) -> double log1pstd/num/double/log1p: (d : double) -> double(xx: double)

fun exp2m1std/num/double/exp2m1: (x : double) -> double( xx: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  expm1std/num/double/expm1: (d : double) -> double(dbl-log2std/num/double/dbl-log2: double *std/core/(*).1: (double, double) -> double xx: double)

// Returns `log(exp(x) + exp(y))`.
// Avoids overlow/underflow errors.
pub fun logaddexpstd/num/double/logaddexp: (x : double, y : double) -> double( xx: double : doublestd/core/types/double: V, yy: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  if xx: double==std/core/(==).2: (double, double) -> boolyy: double then xx: double +std/core/(+).2: (double, double) -> double dbl-log2std/num/double/dbl-log2: double else
    val zz: double = xx: double -std/core/(-).2: (double, double) -> double yy: double
    if zz: double >std/core/(>).2: (double, double) -> bool 0.0 then xx: double +std/core/(+).2: (double, double) -> double log1pstd/num/double/log1p: (d : double) -> double(expstd/num/double/exp: (p : double) -> double(~std/core/(~).1: (d : double) -> doublezz: double))
               else yy: double +std/core/(+).2: (double, double) -> double log1pstd/num/double/log1p: (d : double) -> double(expstd/num/double/exp: (p : double) -> double(zz: double))


// Returns `log2( exp2(x) + exp2(y) )`.
// Avoids overlow/underflow errors.
pub fun logaddexp2std/num/double/logaddexp2: (x : double, y : double) -> double( xx: double : doublestd/core/types/double: V, yy: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  if xx: double==std/core/(==).2: (double, double) -> boolyy: double then xx: double +std/core/(+).2: (double, double) -> double 1.0 else
    val zz: double = xx: double -std/core/(-).2: (double, double) -> double yy: double
    if zz: double >std/core/(>).2: (double, double) -> bool 0.0 then xx: double +std/core/(+).2: (double, double) -> double log2p1std/num/double/log2p1: (x : double) -> double(exp2std/num/double/exp2: (p : double) -> double(~std/core/(~).1: (d : double) -> doublezz: double))
               else yy: double +std/core/(+).2: (double, double) -> double log2p1std/num/double/log2p1: (x : double) -> double(exp2std/num/double/exp2: (p : double) -> double(zz: double))


// Return if two doubles are nearly equal with respect to some `epsilon` (=`8*dbl-epsilon`).
// The epsilon is the nearest difference for numbers around 1.0. The routine automatically
// scales the epsilon for larger and smaller numbers, and for subnormal numbers.
pub fun nearly-eqstd/num/double/nearly-eq: (x : double, y : double, epsilon : ?double) -> bool( xx: double : doublestd/core/types/double: V, yy: double : doublestd/core/types/double: V, epsilonepsilon: ?double : doublestd/core/types/optional: V -> V = 8.0*std/core/(*).1: (double, double) -> doubledbl-epsilonstd/num/double/dbl-epsilon: double ) : boolstd/core/types/bool: V 
  if xx: double ==std/core/(==).2: (double, double) -> bool yy: double returnreturn: bool Truestd/core/types/True: bool
  val diffdiff: double = (xx: double -std/core/(-).2: (double, double) -> double yy: double).absstd/core/abs.1: (d : double) -> double
  if xx: double==std/core/(==).2: (double, double) -> bool0.0 ||std/core/types/(||): (x : bool, y : bool) -> bool yy: double==std/core/(==).2: (double, double) -> bool0.0 ||std/core/types/(||): (x : bool, y : bool) -> bool diffdiff: double <std/core/(<).4: (double, double) -> bool dbl-minstd/num/double/dbl-min: double then
    // very close to zero, scale the epsilon for denormalized numbers
    (2.0*std/core/(*).1: (double, double) -> doublediffdiff: double <std/core/(<).4: (double, double) -> bool (epsilonepsilon: double *std/core/(*).1: (double, double) -> double dbl-minstd/num/double/dbl-min: double))
  else 
    val sumsum: double = xx: double.absstd/core/abs.1: (d : double) -> double +std/core/(+).2: (double, double) -> double yy: double.absstd/core/abs.1: (d : double) -> double
    ((2.0*std/core/(*).1: (double, double) -> doublediffdiff: double /std/core/(/).1: (double, double) -> double (if sumsum: double >std/core/(>).2: (double, double) -> bool dbl-maxstd/num/double/dbl-max: double then dbl-maxstd/num/double/dbl-max: double else sumsum: double)) <std/core/(<).4: (double, double) -> bool epsilonepsilon: double)
  

// Return if two doubles are nearly equal with respect to an `epsilon` of `8*dbl-epsilon`.
// See also `nearly-eq` which takes an explicit `epsilon`.
pub fun (~=)std/num/double/(~=): (x : double, y : double) -> bool(xx: double : doublestd/core/types/double: V, yy: double : doublestd/core/types/double: V ) : boolstd/core/types/bool: V
  nearly-eqstd/num/double/nearly-eq: (x : double, y : double, epsilon : ?double) -> bool(xx: double,yy: double)

//-----------------------------------------
// Ldexp/Frexp
//-----------------------------------------

// Low-level: return the `(lo,hi)` bits of a 64-bit double.
pub extern double-to-bitsstd/num/double/double-to-bits: (d : double) -> (int32, int32)( dd: double : doublestd/core/types/double: V ) : (std/core/types/(,): (V, V) -> Vint32std/core/types/int32: V,int32std/core/types/int32: V)std/core/types/(,): (V, V) -> V
  c  "kk_double_to_bits"
  cs "Primitive.DoubleToBits"
  js "_double_to_bits"


// Low-level: create a `:double` from the given `lo` and `hi` 32-bits.
pub extern double-from-bitsstd/num/double/double-from-bits: (lo : int32, hi : int32) -> double( lolo: int32 : int32std/core/types/int32: V, hihi: int32 : int32std/core/types/int32: V ) : doublestd/core/types/double: V
  c  "kk_double_from_bits"
  cs "Primitive.DoubleFromBits"
  js "_double_from_bits"

/*
// Low-level: create a `:double` from the given `lo` and `hi` 32-bits.
pub fun double-from-bits( lo : int32, hi : int32 ) : double
  val hifrac = hi & 0xFFFFF.int32
  val exp  = (hi.shr(20.int32) & 0x07FF.int32)
  val is-pos = hi.shr(31.int32).zero?
  val p = if lo.is-zero && hifrac.is-zero then
            if exp.is-zero then 0.0
            elif exp==0x7FF.int32 then posinf
            else exp2( (exp - 1023.int32).double )

          else
            val frac = lo.double + 0x1.0p32 * hifrac.double
            if exp.is-zero then frac * exp2( -1074.0 )  // denormalized number
            elif exp==0x7FF.int32 then nan
            else (1.0 + frac * 0x1.0p-52) * exp2( (exp - 1023.int32).double )

  if is-pos then p else ~p

*/

// Calculate 2&middot;^`e`^ for an integer `e`.
// Uses efficient bit conversion for exponents between  -1022 and 1023 and
// otherwise falls back to the regular `exp2` function converting `e` to a double.
pub fun exp2std/num/double/exp2.1: (e : int) -> double( ee: int : intstd/core/types/int: V ) : doublestd/core/types/double: V
  if ee: int >=std/core/(>=).1: (x : int, y : int) -> bool -1022 &&std/core/types/(&&): (x : bool, y : bool) -> bool ee: int <=std/core/(<=).1: (x : int, y : int) -> bool 1023
    then double-from-bitsstd/num/double/double-from-bits: (lo : int32, hi : int32) -> double(zerostd/num/int32/zero: int32, (1023 +std/core/(+).4: (x : int, y : int) -> int ee: int).int32std/core/int32: (i : int) -> int32.shlstd/num/int32/shl: (int32, int32) -> int32(20.int32std/core/int32: (i : int) -> int32))
    else exp2std/num/double/exp2: (p : double) -> double( ee: int.doublestd/core/double: (i : int) -> double )


// Create a double `d` given a mantissa `man` and exponent `exp`
// such that `man`&middot;2^`exp`^ =  `d` exactly (if it is representable
// by a `:double`). See also `ldexp`.
pub fun encodestd/num/double/encode: (man : int, exp : int) -> double( manman: int : intstd/core/types/int: V, expexp: int : intstd/core/types/int: V ) : doublestd/core/types/double: V
  ldexpstd/num/double/ldexp: (x : double, e : int) -> double(manman: int.doublestd/core/double: (i : int) -> double,expexp: int)

val one-p1023std/num/double/one-p1023: double = 0x1.0p1023
val one-m1022std/num/double/one-m1022: double = 0x1.0p-1022   // = dbl-min

fun mul-exp2std/num/double/mul-exp2: (x : double, e : int) -> double( xx: double : doublestd/core/types/double: V, ee: int : intstd/core/types/int: V ) : doublestd/core/types/double: V 
  xx: double *std/core/(*).1: (double, double) -> double exp2std/num/double/exp2.1: (e : int) -> double(ee: int)

// 'Load exponent': returns `x`&middot;2^`e`^ for a `is-finite` `x` and
// otherwise `x` itself. See also `encode` which loads an integer mantissa.
pub fun ldexpstd/num/double/ldexp: (x : double, e : int) -> double( xx: double : doublestd/core/types/double: V, ee: int : intstd/core/types/int: V ) : doublestd/core/types/double: V
  if !std/core/types/(!).1: (b : bool) -> boolis-finitestd/num/double/is-finite: (d : double) -> bool(xx: double) then xx: double
  elif ee: int >=std/core/(>=).1: (x : int, y : int) -> bool -1022 then
    if ee: int <=std/core/(<=).1: (x : int, y : int) -> bool 1023 then mul-exp2std/num/double/mul-exp2: (x : double, e : int) -> double(xx: double,ee: int)  // usually this branch
    elif ee: int <=std/core/(<=).1: (x : int, y : int) -> bool 2046 then mul-exp2std/num/double/mul-exp2: (x : double, e : int) -> double( xx: double*std/core/(*).1: (double, double) -> doubleone-p1023std/num/double/one-p1023: double, ee: int -std/core/(-).4: (x : int, y : int) -> int 1023 )
    elif ee: int <=std/core/(<=).1: (x : int, y : int) -> bool 3069 then mul-exp2std/num/double/mul-exp2: (x : double, e : int) -> double( xx: double*std/core/(*).1: (double, double) -> doubleone-p1023std/num/double/one-p1023: double*std/core/(*).1: (double, double) -> doubleone-p1023std/num/double/one-p1023: double, ee: int -std/core/(-).4: (x : int, y : int) -> int 2046 )
    elif xx: double <std/core/(<).4: (double, double) -> bool 0.0 then neginfstd/num/double/neginf: double else posinfstd/num/double/posinf: double

  else
    if ee: int >=std/core/(>=).1: (x : int, y : int) -> bool -2044 then mul-exp2std/num/double/mul-exp2: (x : double, e : int) -> double(xx: double*std/core/(*).1: (double, double) -> doubleone-m1022std/num/double/one-m1022: double, ee: int +std/core/(+).4: (x : int, y : int) -> int 1022)
    elif ee: int >=std/core/(>=).1: (x : int, y : int) -> bool -3066 then mul-exp2std/num/double/mul-exp2: (x : double, e : int) -> double(xx: double*std/core/(*).1: (double, double) -> doubleone-m1022std/num/double/one-m1022: double*std/core/(*).1: (double, double) -> doubleone-m1022std/num/double/one-m1022: double, ee: int +std/core/(+).4: (x : int, y : int) -> int 2044 )
    elif xx: double <std/core/(<).4: (double, double) -> bool 0.0 then -0.0 else 0.0


// Decode a double `d` into a tuple `(m,e)` of a mantissa `m` and exponent `e`
// such that `m`&middot;2^`e`^ =  `d` exactly. The mantissa `m` is
// always either 0 or in the range [2^52^, 2^53^). See also `frexp`.
pub fun decodestd/num/double/decode: (d : double) -> (int, int)( dd: double : doublestd/core/types/double: V ) : (std/core/types/(,): (V, V) -> Vintstd/core/types/int: V,intstd/core/types/int: V)std/core/types/(,): (V, V) -> V
  if dd: double==std/core/(==).2: (double, double) -> bool0.0 then (std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)0,0)std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)
  elif dd: double.is-subnormalstd/num/double/is-subnormal: (d : double) -> bool then decode-normalizedstd/num/double/decode-normalized: (d : double, e-adjust : ?int) -> (int, int)(dd: double *std/core/(*).1: (double, double) -> double 0x1.0p54, -54)
  else decode-normalizedstd/num/double/decode-normalized: (d : double, e-adjust : ?int) -> (int, int)(dd: double,0)

// decode a normalized double (i.e. not subnormal)
fun decode-normalizedstd/num/double/decode-normalized: (d : double, e-adjust : ?int) -> (int, int)( dd: double : doublestd/core/types/double: V, e-adjuste-adjust: ?int : intstd/core/types/optional: V -> V = 0 ) : (std/core/types/(,): (V, V) -> Vintstd/core/types/int: V,intstd/core/types/int: V)std/core/types/(,): (V, V) -> V
  val (std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)lolo: int32,hihi: int32)std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)  = double-to-bitsstd/num/double/double-to-bits: (d : double) -> (int32, int32)(dd: double)
  val signsign: int  = if hihi: int32.is-negstd/num/int32/is-neg: (i : int32) -> bool then -1 else 1
  val expexp: int  = andstd/num/int32/and: (int32, int32) -> int32(hihi: int32.shrstd/num/int32/shr: (int32, int32) -> int32(20.int32std/core/int32: (i : int) -> int32),0x7FF.int32std/core/int32: (i : int) -> int32).intstd/core/int.1: (i : int32) -> int -std/core/(-).4: (x : int, y : int) -> int 1043
  val manman: int  = andstd/num/int32/and: (int32, int32) -> int32(hihi: int32,0xFFFFF.int32std/core/int32: (i : int) -> int32).intstd/core/int.1: (i : int32) -> int +std/core/(+).4: (x : int, y : int) -> int 0x100000
  (std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)signsign: int*std/core/(*): (int, int) -> int(manman: int *std/core/(*): (int, int) -> int 0x100000000 +std/core/(+).4: (x : int, y : int) -> int lolo: int32.uintstd/num/int32/uint: (i : int32) -> int), expexp: int -std/core/(-).4: (x : int, y : int) -> int 32 +std/core/(+).4: (x : int, y : int) -> int e-adjuste-adjust: int)std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)

// 'Fraction/exponent': return the normalized fraction `f` and exponent `exp`
// for a number `x` such that `x == f`&middot;2^`exp`^.
// The absolute value of the fraction `f` is always in the range [0.5, 1.0), or
// one of `0.0`, `-0.0`, `neginf`, `posinf`, or `nan`.
// See also `decode` which  decodes to an integer mantissa.
pub fun frexpstd/num/double/frexp: (x : double) -> (double, int)( xx: double : doublestd/core/types/double: V ) : (std/core/types/(,): (V, V) -> Vdoublestd/core/types/double: V, intstd/core/types/int: V)std/core/types/(,): (V, V) -> V
  if !std/core/types/(!).1: (b : bool) -> boolxx: double.is-finitestd/num/double/is-finite: (d : double) -> bool ||std/core/types/(||): (x : bool, y : bool) -> bool xx: double.is-negzerostd/num/double/is-negzero: (d : double) -> bool returnreturn: (double, int) (std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)xx: double,0)std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)
  val (std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)mm: int,ee: int)std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b) = xx: double.decodestd/num/double/decode: (d : double) -> (int, int)
  (std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)mm: int.doublestd/core/double: (i : int) -> double *std/core/(*).1: (double, double) -> double 0x1.0p-53, ee: int +std/core/(+).4: (x : int, y : int) -> int 53  )std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)


//-----------------------------------------
// Show in hexadecimal
//-----------------------------------------

/* Show a double in [hexadecimal notation](https://books.google.com/books?id=FgMsCwAAQBAJ&pg=PA41).
An advantage of this format is that it precisely represents the `:double` and can
reliably (and efficiently) be parsed back, i.e. `d.show-hex.parse-double == Just(d)`.
The significant is the _hexadecimal_ fraction while the
exponent after the `p` is the _decimal_ power of 2.
 For example, ``0xA.Fp-10`` = (10 + 15/16)&middot;2^-10^  (not 2^-16^!) = 0.01068115234375.
 Equivalently, ``0xA.Fp-10 == 0x5.78p-9 == 0x2.BCp-8 == 0x1.5Ep-7``.
```
> dbl-min.show-hex
"0x1.0p-1022"
> 0.1.show-hex
"0x1.999999999999Ap-4"
> dbl-max.show-hex
"0x1.FFFFFFFFFFFFFp+1023"
> -0.0.show-hex
"-0x0.0p+0"
> nan.show-hex
"NaN"
> 0.01068115234375.show-hex
"0x1.5Ep-7"
```
.
*/
pub fun show-hexstd/num/double/show-hex: (d : double, width : ?int, use-capitals : ?bool, pre : ?string) -> string( dd: double : doublestd/core/types/double: V, widthwidth: ?int : intstd/core/types/optional: V -> V = 1, use-capitalsuse-capitals: ?bool : boolstd/core/types/optional: V -> V = Truestd/core/types/True: bool, prepre: ?string : stringstd/core/types/optional: V -> V = "0x" ) : stringstd/core/types/string: V
  if !std/core/types/(!).1: (b : bool) -> booldd: double.is-finitestd/num/double/is-finite: (d : double) -> bool then dd: double.showstd/core/show.1: (d : double, precision : ?int) -> string else
    val (std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)mm: int,ee: int)std/core/types/(,): forall<a,b> (fst : a, snd : b) -> (a, b)  = dd: double.decodestd/num/double/decode: (d : double) -> (int, int)
    val manman: string  = mm: int.absstd/core/abs: (i : int) -> int.show-hexstd/core/show-hex: (i : int, width : ?int, use-capitals : ?bool, pre : ?string) -> string(1,use-capitalsuse-capitals: bool,"")
    val exp0exp0: int  = ee: int +std/core/(+).4: (x : int, y : int) -> int 4*std/core/(*): (int, int) -> int(manman: string.countstd/core/count.1: (s : string) -> int -std/core/(-).4: (x : int, y : int) -> int 1)
    val expexp: string  = (if exp0exp0: int >=std/core/(>=).1: (x : int, y : int) -> bool 0 then "+" else "") ++std/core/(++).1: (x : string, y : string) -> string exp0exp0: int.showstd/core/show: (i : int) -> string
    val fracfrac: string  = manman: string.tailstd/core/tail.2: (s : string) -> string.trim-rightstd/core/trim-right.1: (s : string, sub : string) -> string("0").pad-rightstd/core/pad-right: (s : string, width : int, fill : ?char) -> string(maxstd/core/max: (i : int, j : int) -> int(1,widthwidth: int),'0')
    val signsign: string  = if dd: double.is-negstd/core/is-neg.3: (d : double) -> bool ||std/core/types/(||): (x : bool, y : bool) -> bool dd: double.is-negzerostd/num/double/is-negzero: (d : double) -> bool then "-" else ""
    signsign: string ++std/core/(++).1: (x : string, y : string) -> string prepre: string ++std/core/(++).1: (x : string, y : string) -> string manman: string.headstd/core/head.3: (s : string) -> string ++std/core/(++).1: (x : string, y : string) -> string "." ++std/core/(++).1: (x : string, y : string) -> string fracfrac: string ++std/core/(++).1: (x : string, y : string) -> string "p" ++std/core/(++).1: (x : string, y : string) -> string expexp: string


//-----------------------------------------
// Parse a double
//-----------------------------------------

// Parse a double number. Can be "NaN", "Inf(inity)" (case-insensitive),
// a fix-point number (`1.2`) or in scientific notation (`-2.3e-5`).
// Also allows floats in [hexadecimal notation](https://books.google.com/books?id=FgMsCwAAQBAJ&pg=PA41) (`0xA.Fp-10`) that can
// be represented precisely (and are the preferred _round trip_ format).
pub fun parse-doublestd/num/double/parse-double: (s : string) -> maybe<double>( ss: string : stringstd/core/types/string: V ) : maybestd/core/types/maybe: V -> V<doublestd/core/types/double: V>
  ss: string.trimstd/core/trim: (s : string) -> string.to-lowerstd/core/to-lower: (s : string) -> string.slicestd/core/slice: (s : string) -> sslice.parse-eofstd/text/parse/parse-eof: forall<a,e> (input : sslice, p : () -> <parse|e> a) -> e parse-error<a>(pdoublestd/num/double/pdouble: () -> parse double).maybestd/text/parse/maybe: forall<a> (perr : parse-error<a>) -> maybe<a>

pub fun pdoublestd/num/double/pdouble: () -> parse double() : parsestd/text/parse/parse: HX doublestd/core/types/double: V
  val negneg: bool = signstd/text/parse/sign: () -> parse bool()
  val dd: double = [std/core/Cons: forall<a> (head : a, tail : list<a>) -> list<a>{ phexdoublestd/num/double/phexdouble: () -> parse double() }, { pdecdoublestd/num/double/pdecdouble: () -> parse double() }, { pspecialstd/num/double/pspecial: () -> parse double() }, { 0.0 }]std/core/Nil: forall<a> list<a>.choosestd/text/parse/choose: forall<a,e> (ps : list<parser<e,a>>) -> <parse|e> a
  if negneg: bool then ~std/core/(~).1: (d : double) -> doubledd: double else dd: double

fun phexdoublestd/num/double/phexdouble: () -> parse double() : parsestd/text/parse/parse: HX doublestd/core/types/double: V
  charstd/text/parse/char: (c : char) -> parse char('0')
  one-ofstd/text/parse/one-of: (chars : string) -> parse char("xX")
  val manman: string  = hex-digitsstd/text/parse/hex-digits: () -> parse string()
  val fracfrac: string = optionalstd/text/parse/optional: forall<a,e> (default : a, p : parser<e,a>) -> <parse|e> a( "", { charstd/text/parse/char: (c : char) -> parse char('.'); hex-digitsstd/text/parse/hex-digits: () -> parse string() }).trim-rightstd/core/trim-right.1: (s : string, sub : string) -> string("0")
  val expexp: int : intstd/core/types/int: V = optionalstd/text/parse/optional: forall<a,e> (default : a, p : parser<e,a>) -> <parse|e> a( 0, { one-ofstd/text/parse/one-of: (chars : string) -> parse char("pP"); pintstd/text/parse/pint: () -> parse int() })
  val mm: int : intstd/core/types/int: V = (manman: string ++std/core/(++).1: (x : string, y : string) -> string fracfrac: string).parse-intstd/core/parse-int: (s : string, hex : ?bool) -> maybe<int>(hex=Truestd/core/types/True: bool).defaultstd/core/default: forall<a> (m : maybe<a>, nothing : a) -> a(0)
  val ee: int : intstd/core/types/int: V = expexp: int -std/core/(-).4: (x : int, y : int) -> int 4*std/core/(*): (int, int) -> intfracfrac: string.countstd/core/count.1: (s : string) -> int
  encodestd/num/double/encode: (man : int, exp : int) -> double(mm: int,ee: int)

fun pdecdoublestd/num/double/pdecdouble: () -> parse double() : parsestd/text/parse/parse: HX doublestd/core/types/double: V
  val curcur: sslice  = current-inputstd/text/parse/current-input: () -> parse sslice()
  val manman: string  = digitsstd/text/parse/digits: () -> parse string()
  val fracfrac: string = optionalstd/text/parse/optional: forall<a,e> (default : a, p : parser<e,a>) -> <parse|e> a("", { charstd/text/parse/char: (c : char) -> parse char('.'); digits0std/text/parse/digits0: () -> parse string() }).trim-rightstd/core/trim-right.1: (s : string, sub : string) -> string("0")
  val expexp: int : intstd/core/types/int: V = optionalstd/text/parse/optional: forall<a,e> (default : a, p : parser<e,a>) -> <parse|e> a( 0, { one-ofstd/text/parse/one-of: (chars : string) -> parse char("eE"); pintstd/text/parse/pint: () -> parse int() })
  //val m : int = (man ++ frac).parse-int-default(0)
  //val e : int = exp - frac.count
  //m.double * exp10(e.double)
  curcur: sslice.stringstd/core/string.3: (slice : sslice) -> string.prim-parse-doublestd/num/double/prim-parse-double: (s : string) -> double  // more precision than simple multiply

fun pspecialstd/num/double/pspecial: () -> parse double() : parsestd/text/parse/parse: HX doublestd/core/types/double: V
  [std/core/Cons: forall<a> (head : a, tail : list<a>) -> list<a>{ pstringstd/text/parse/pstring: (s : string) -> parse string("nan"); nanstd/num/double/nan: double },
   { pstringstd/text/parse/pstring: (s : string) -> parse string("infinity"); posinfstd/num/double/posinf: double },
   { pstringstd/text/parse/pstring: (s : string) -> parse string("inf"); posinfstd/num/double/posinf: double }
  ]std/core/Nil: forall<a> list<a>.choosestd/text/parse/choose: forall<a,e> (ps : list<parser<e,a>>) -> <parse|e> a

// Return `nan` on failure
extern prim-parse-doublestd/num/double/prim-parse-double: (s : string) -> double( ss: string : stringstd/core/types/string: V ) : doublestd/core/types/double: V
  c  "kk_prim_parse_double"
  cs "Primitive.DoubleParse"
  js "parseFloat"


//-----------------------------------------
// Various
//-----------------------------------------

// Return the sum of a list of doubles.
// Uses [Kahan-Babu&scaron;kan-Neumaier summation](https://en.wikipedia.org/wiki/Kahan_summation_algorithm#Further_enhancements)
// to minimize rounding errors. This
// is more precise as Kahan summation and about as fast.\
// `[1.0e3,1.0e97,1.0e3,-1.0e97].sum == 2000.0`\
// A. Neumaier, _Rundungsfehleranalyse einiger Verfahren zur Summation endlicher Summen_.
// Math. Mechanik, 54:39--51, 1974.
pub fun sumstd/num/double/sum: (xs : list<double>) -> double( xsxs: list<double> : liststd/core/list: V -> V<doublestd/core/types/double: V> ) : doublestd/core/types/double: V
  varstd/core/hnd/local-var: forall<a,b,e,h> (init : a, action : (l : local-var<h,a>) -> <local<h>|e> b) -> <local<h>|e> b totalstd/core/types/local-scope: forall<a,e> (action : forall<h> () -> <local<h>|e> a) -> e a := 0.0
  varstd/core/hnd/local-var: forall<a,b,e,h> (init : a, action : (l : local-var<h,a>) -> <local<h>|e> b) -> <local<h>|e> b compcomp: local-var<$338,double>  := 0.0
  xsxs: list<double>.foreachstd/core/foreach: forall<a,e> (xs : list<a>, action : (a) -> e ()) -> e () fn(xx: double)
    val tt: double = totaltotal: double +std/core/(+).2: (double, double) -> double xx: double
    val cc: double = if totaltotal: double.absstd/core/abs.1: (d : double) -> double >=std/core/(>=).3: (double, double) -> bool xx: double.absstd/core/abs.1: (d : double) -> double then (totaltotal: double -std/core/(-).2: (double, double) -> double tt: double) +std/core/(+).2: (double, double) -> double xx: double else (xx: double -std/core/(-).2: (double, double) -> double tt: double) +std/core/(+).2: (double, double) -> double totaltotal: double
    compcomp: local-var<$338,double>  :=std/core/types/local-set: forall<a,e,h> (v : local-var<h,a>, assigned : a) -> <local<h>|e> () compcomp: double +std/core/(+).2: (double, double) -> double cc: double
    totaltotal: local-var<$338,double> :=std/core/types/local-set: forall<a,e,h> (v : local-var<h,a>, assigned : a) -> <local<h>|e> () tt: double
  totaltotal: double +std/core/(+).2: (double, double) -> double compcomp: double;

// The hypotenuse of `x` and `y`: `sqrt(x*x + y*y)`.
// Prevents overflow for large numbers.
pub fun hypotstd/num/double/hypot: (x : double, y : double) -> double( xx: double : doublestd/core/types/double: V, yy: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  val xxxx: double = absstd/core/abs.1: (d : double) -> double(xx: double)
  val yyyy: double = absstd/core/abs.1: (d : double) -> double(yy: double)
  val lolo: double = minstd/core/min.1: (x : double, y : double) -> double(xxxx: double,yyyy: double)
  val hihi: double = maxstd/core/max.1: (x : double, y : double) -> double(xxxx: double,yyyy: double)
  if hihi: double==std/core/(==).2: (double, double) -> bool0.0 then 0.0 else
    val zz: double  = lolo: double /std/core/(/).1: (double, double) -> double hihi: double
    hihi: double *std/core/(*).1: (double, double) -> double sqrtstd/num/double/sqrt: (d : double) -> double( 1.0 +std/core/(+).2: (double, double) -> double zz: double*std/core/(*).1: (double, double) -> doublezz: double )


// The square root of the sum of the squares of three doubles.
// Prevents overflow for large numbers.
pub fun hypotstd/num/double/hypot.1: (x : double, y : double, z : double) -> double( xx: double : doublestd/core/types/double: V, yy: double : doublestd/core/types/double: V, zz: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  val xxxx: double = absstd/core/abs.1: (d : double) -> double(xx: double)
  val yyyy: double = absstd/core/abs.1: (d : double) -> double(yy: double)
  val zzzz: double = absstd/core/abs.1: (d : double) -> double(zz: double)
  val hihi: double = maxstd/core/max.1: (x : double, y : double) -> double(maxstd/core/max.1: (x : double, y : double) -> double(xxxx: double,yyyy: double),zzzz: double)
  if hihi: double==std/core/(==).2: (double, double) -> bool0.0 then 0.0 else
    hihi: double *std/core/(*).1: (double, double) -> double sqrtstd/num/double/sqrt: (d : double) -> double( sqrstd/num/double/sqr: (x : double) -> double(xxxx: double /std/core/(/).1: (double, double) -> double hihi: double) +std/core/(+).2: (double, double) -> double sqrstd/num/double/sqr: (x : double) -> double(yyyy: double /std/core/(/).1: (double, double) -> double hihi: double) +std/core/(+).2: (double, double) -> double sqrstd/num/double/sqr: (x : double) -> double(zzzz: double /std/core/(/).1: (double, double) -> double hihi: double) )


// The square root of the sum of squares of a list of doubles.
// Prevents overflow for large numbers and uses Kahan-Babu&scaron;kan-Neumaier summation
// for precision.
pub fun hypotstd/num/double/hypot.2: (xs : list<double>) -> double( xsxs: list<double> : liststd/core/list: V -> V<doublestd/core/types/double: V> ) : doublestd/core/types/double: V
  val hihi: double = xsxs: list<double>.abs-maxstd/num/double/abs-max.1: (xs : list<double>) -> double
  if hihi: double==std/core/(==).2: (double, double) -> bool0.0 then 0.0 else
    hihi: double *std/core/(*).1: (double, double) -> double xsxs: list<double>.mapstd/core/map.5: forall<a,b,e> (xs : list<a>, f : (a) -> e b) -> e list<b>( fn(xx: double) sqrstd/num/double/sqr: (x : double) -> double(xx: double /std/core/(/).1: (double, double) -> double hihi: double) ).sumstd/num/double/sum: (xs : list<double>) -> double.sqrtstd/num/double/sqrt: (d : double) -> double


// The square of a double
pub fun sqrstd/num/double/sqr: (x : double) -> double(xx: double : doublestd/core/types/double: V ): doublestd/core/types/double: V
  xx: double*std/core/(*).1: (double, double) -> doublexx: double

// The maximum of the absolute values.
pub fun abs-maxstd/num/double/abs-max: (x : double, y : double) -> double( xx: double : doublestd/core/types/double: V, yy: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  maxstd/core/max.1: (x : double, y : double) -> double(absstd/core/abs.1: (d : double) -> double(xx: double),absstd/core/abs.1: (d : double) -> double(yy: double))

// The maximum of a list of absolute values.
pub fun abs-maxstd/num/double/abs-max.1: (xs : list<double>) -> double( xsxs: list<double> : liststd/core/list: V -> V<doublestd/core/types/double: V> ) : doublestd/core/types/double: V
  xsxs: list<double>.foldlstd/core/foldl: forall<a,b,e> (list<a>, b, (b, a) -> e b) -> e b(0.0, fn(mm: double,xx: double) maxstd/core/max.1: (x : double, y : double) -> double(absstd/core/abs.1: (d : double) -> double(xx: double),mm: double) )

//-----------------------------------------
// Trigonometry
//-----------------------------------------

val rad2degstd/num/double/rad2deg: double : doublestd/core/types/double: V = 180.0/std/core/(/).1: (double, double) -> doublepistd/num/double/pi: double
val deg2radstd/num/double/deg2rad: double : doublestd/core/types/double: V = pistd/num/double/pi: double/std/core/(/).1: (double, double) -> double180.0

// Convert radians to degrees.
pub fun degstd/num/double/deg: (rad : double) -> double( radrad: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  radrad: double *std/core/(*).1: (double, double) -> double rad2degstd/num/double/rad2deg: double

// Convert degrees to radians.
pub fun radstd/num/double/rad: (deg : double) -> double( degdeg: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  degdeg: double *std/core/(*).1: (double, double) -> double deg2radstd/num/double/deg2rad: double


// Return `x` with the sign of `y`.
pub fun with-sign-ofstd/num/double/with-sign-of: (x : double, y : double) -> double( xx: double : doublestd/core/types/double: V, yy: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  if yy: double <std/core/(<).4: (double, double) -> bool 0.0 then ~std/core/(~).1: (d : double) -> double(xx: double.absstd/core/abs.1: (d : double) -> double) else xx: double.absstd/core/abs.1: (d : double) -> double


// Return the sine of an angle in radians `d`.
pub inline extern sinstd/num/double/sin: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "sin"
  cs "Math.Sin"
  js "Math.sin"

// Return the cosine of an angle in radians `d`.
pub inline extern cosstd/num/double/cos: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "cos"
  cs "Math.Cos"
  js "Math.cos"

// Return the tangent of an angle in radians `d`.
pub inline extern tanstd/num/double/tan: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "tan"
  cs "Math.Tan"
  js "Math.tan"

// Return the arc-tangent of `d`
pub inline extern atanstd/num/double/atan: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "atan"
  cs "Math.Atan"
  js "Math.atan"

// Return the arc-tangent of a point (`x`,`y`).
pub inline extern atan2std/num/double/atan2: (x : double, y : double) -> double( x : doublestd/core/types/double: V, y : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "atan2"
  cs "Math.Atan2"
  js "Math.atan2"

// Return the arc-cosine of `d`
pub inline extern acosstd/num/double/acos: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "acos"
  cs "Math.Acos"
  js "Math.acos"

// Return the arc-sine of `d`
pub inline extern asinstd/num/double/asin: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "asin"
  cs "Math.Asin"
  js "Math.asin"


// The hyperbolic tangent of `d`
pub inline extern tanhstd/num/double/tanh: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "tanh"
  cs "Math.Tanh"
  js "Math.tanh"

// The hyperbolic cosine of `d`
pub inline extern coshstd/num/double/cosh: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "cosh"
  cs "Math.Cosh"
  js "Math.cosh"

// The hyperbolic sine of `d`
pub inline extern sinhstd/num/double/sinh: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "sinh"
  cs "Math.Sinh"
  js "Math.sinh"

// Return `log(1.0 + x)`.
// Avoids potential imprecision for small `x` where adding `1.0` explicitly
// may lead to rounding errors.
pub inline extern log1pstd/num/double/log1p: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "log1p"
  js "Math.log1p"

// Return `exp(x - 1.0)`.
// Avoids rounding errors for values of `x` very close to `1.0`.
pub inline extern expm1std/num/double/expm1: (d : double) -> double( d : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "log1p"
  js "Math.log1p"

// The area hyperbolic tangent of `d`
pub extern atanhstd/num/double/atanh: (d : double) -> double( dd: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "atanh"
  js "Math.atanh"

// The area hyperbolic cosine of `d`
pub extern acoshstd/num/double/acosh: (d : double) -> double( dd: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c  "acosh"
  js "Math.acosh"

// The area hyperbolic sine of `d`
pub extern asinhstd/num/double/asinh: (d : double) -> double( dd: double : doublestd/core/types/double: V ) : doublestd/core/types/double: V
  c "asinh"
  js "Math.asinh"


// ------------------------
// explicit implementations
// ------------------------
  
/*  
// Return `log(1.0 + x)`.
// Avoids potential imprecision for small `x` where adding `1.0` explicitly
// may lead to rounding errors.
pub fun log1p( x : double ) : double
  if x.is-posinf then x else
    val y = 1.0 + x
    val z = y - 1.0
    if z==0.0 then x else
      log(y) * (x / z)


// Return `exp(x - 1.0)`.
// Avoids rounding errors for values of `x` very close to `1.0`.
pub fun expm1( x : double ) : double
  if x.is-posinf then x else
    val y = exp(x)
    if y==1.0 then x
    elif y - 1.0 == -1.0 then -1.0
    else (y - 1.0) * (x / log(y))


// The area hyperbolic tangent of `x`
pub fun atanh( x : double ) : double
  0.5*log( (1.0 + x) / (1.0 - x))
  /*
  // 0.5*log( (1.0 + x) / (1.0 - x)) = 0.5*log1p( (2.0*x) / (1-x) )
  if x.abs >= 0.5 then 0.5*log1p(2.0*(x/(1.0 - x))) else
    val x2 = x + x
    0.5*log1p(x2 + (x2*x)/(1.0 - x))

  */

// The area hyperbolic cosine of `x`
pub fun acosh( x : double ) : double
  // log(x + sqrt((x - 1.0)*(x + 1.0)));
  if x > 0x1.0p28 then dbl-log2 + log(x)
  elif x > 2.0 then
    log(2.0*x - 1.0/(x + sqrt(x.sqr - 1.0)))

  elif x <= 1.0 then (if x < 1.0 then nan else 0.0)
  else
    val xm1 = x - 1.0
    log1p( xm1 + sqrt(2.0*xm1 + xm1*xm1) )


// The area hyperbolic sine of `x`
pub fun asinh( x : double ) : double
  //log( x + sqrt(x.sqr + 1.0))
  val xa = x.abs
  if xa > 0x1.0p28 then (log(xa) + dbl-log2).with-sign-of(x)
  elif xa >= 2.0 then (log(2.0*xa + 1.0/(xa + sqrt(xa.sqr + 1.0)))).with-sign-of(x)
  elif xa == 0.0 then 0.0
  else
    val xa2 = xa.sqr
    log1p( xa + xa2/(1.0 + sqrt( xa2 + 1.0 ))).with-sign-of(x)


*/