math_utils/math_utils.mqh

//+------------------------------------------------------------------+
//|                                                   math_utils.mqh |
//|                                        Copyright © 2026, Amr Ali |
//|                             https://www.mql5.com/en/users/amrali |
//+------------------------------------------------------------------+
#property copyright "Copyright © 2026, Amr Ali"
#property link      "https://www.mql5.com/en/users/amrali"
#property version   "4.50"
#property description "Handy functions for comparison, rounding, formatting, and debugging of doubles (prices, lots, and money)."

#ifdef __MQL4__
#property strict
#endif

#ifndef MATH_UTILS_UNIQUE_HEADER_ID_H
#define MATH_UTILS_UNIQUE_HEADER_ID_H

//+------------------------------------------------------------------+
//| IEEE-754 binary64 bit-layout constants (used in Section IV).     |
//+------------------------------------------------------------------+
#define MU_MANTISSA_BITS    52                     // fraction bits
#define MU_EXPONENT_BIAS    1023                   // exponent bias
#define MU_EXPONENT_MASK    0x7FF                  // 11-bit exponent field
#define MU_MANTISSA_MASK    0xFFFFFFFFFFFFF        // 52-bit fraction field
#define MU_SIGN_BIT         0x8000000000000000     // sign bit, as a long/ulong bit pattern


//+------------------------------------------------------------------+
//| Functions in the library.                                        |
//+------------------------------------------------------------------+
/*
I. Handy functions for comparison of doubles
--------------------------------------------
int    Compare(const double a, const double b, const int digits);
bool   EQ(const double a, const double b, const int digits = 8);
bool   NE(const double a, const double b, const int digits = 8);
bool   GT(const double a, const double b, const int digits = 8);
bool   LT(const double a, const double b, const int digits = 8);
bool   GTE(const double a, const double b, const int digits = 8);
bool   LTE(const double a, const double b, const int digits = 8);
bool   IsZero(const double value, const int digits = 8);
bool   AlmostEqual(const double a, const double b);
bool   EqualDoubles(const double a, const double b, const int significantDigits = 15);
bool   IsClose(const double a, const double b, const int maxDifferentSignificantDigits = 2);
int    GetEqualDigits(const double a, const double b);
int    GetEqualSigDigits(const double a, const double b);

II. Handy functions for rounding of doubles
-------------------------------------------
double Ceil(const double v);
double Ceil(const double value, const int digits);
double Ceil(const double value, const double step);
double Floor(const double v);
double Floor(const double value, const int digits);
double Floor(const double value, const double step);
double Round(const double v);
double Round(const double value, const int digits);
double Round(const double value, const double step);
double Trunc(const double v);
double Trunc(const double value, const int digits);
double Trunc(const double value, const double step);
double RoundToSignificantDigits(const double value, int digits);
double RoundPrice(const double pPrice, const string pSymbol = NULL);
double RoundVolume(const double pVolume, const string pSymbol = NULL);
double StripError(const double value);

III. Miscellaneous handy functions for doubles
----------------------------------------------
bool   IsInteger(const double value);
bool   IsRound(const double value, const int digits);
bool   IsRound(const double value, const double step);
double Frac(const double value);
int    GetDigits(const double value);
int    GetIntegerDigits(const double value);
int    GetSignificantDigits(double value);
double Remap(const double value, const double min, const double max, const double destMin, const double destMax);
double Lerp(const double start, const double end, const double t);
double Normalize(const double value, const double min, const double max);
double Denormalize(const double value, const double destMin, const double destMax);
T      Clamp(const T value, const T min, const T max);
T      ClampMax(const T value, const T max);
T      ClampMin(const T value, const T min);
T      Wrap(const T value, const T min, const T max);
double safeDiv(const double a, const double b);
int    FloorLog10(double value);
double GetPower10(const int power);
int    SignDbl(const double value);

IV. Handy functions for low-level binary operations on doubles
--------------------------------------------------------------
long   DoubleToLongBits(const double value);
double LongBitsToDouble(const long bits);
long   DoubleToOrderedLong(const double value);
double OrderedLongToDouble(const long bits);
int    RawExponent(const double value);
long   RawSignificand(const double value);
int    GetExponent(const double value);
double GetSignificand(const double value);
bool   IsExactDouble(const double value);
double NextAwayFromZero(const double value);
double DoubleAdvance(const double value, const long distance);
double Ulp(const double value);
long   UlpDiff(const double a, const double b);

V. Handy functions for formatting of doubles to string
------------------------------------------------------
string DoubleToHexadecimal(const double value);
string DoubleToHexFloatConstant(const double value);
string DoubleToStringExact(const double value);
string DoubleToExponential(const double value, int digits = -1);
string Repr(const double value);
string Repr(const float value);
string DoubleDebugInfo(const double value);
string FormatDouble(const double value, const int digits, const string separator=",");
*/

//+------------------------------------------------------------------+
//| I. Handy functions for comparison of doubles.                    |
//+------------------------------------------------------------------+

// GT
//---------------------> double a
// EQ (1/2 point zone)
//---------------------> double b
// LT

//+------------------------------------------------------------------+
//| Compares two doubles up to a specified decimal precision.        |
//| Returns 0 if equal, 1 if a > b, -1 if a < b.                     |
//| digits: number of decimal places of tolerance. Unlike            |
//|         NormalizeDouble(), this is NOT limited to 0-8, so it is  |
//|         safe to pass larger precisions.                          |
//+------------------------------------------------------------------+
int Compare(const double a, const double b, const int digits)
  {
// if(NormalizeDouble(a - b, digits) == 0.0)  
   if(MathAbs(a - b) < GetPower10(-digits) * 0.5)
      return 0;

   return (a > b) ? 1 : -1;
  }
//+------------------------------------------------------------------+
//| Is Equal.                                                        |
//+------------------------------------------------------------------+
bool EQ(const double a, const double b, const int digits = 8)
  {
   return Compare(a, b, digits) == 0;
  }
//+------------------------------------------------------------------+
//| Is Not Equal.                                                    |
//+------------------------------------------------------------------+
bool NE(const double a, const double b, const int digits = 8)
  {
   return Compare(a, b, digits) != 0;
  }
//+------------------------------------------------------------------+
//| Is Greater Than.                                                 |
//+------------------------------------------------------------------+
bool GT(const double a, const double b, const int digits = 8)
  {
   return Compare(a, b, digits) > 0;
  }
//+------------------------------------------------------------------+
//| Is Less Than.                                                    |
//+------------------------------------------------------------------+
bool LT(const double a, const double b, const int digits = 8)
  {
   return Compare(a, b, digits) < 0;
  }
//+------------------------------------------------------------------+
//| Greater Than or Equal.                                           |
//+------------------------------------------------------------------+
bool GTE(const double a, const double b, const int digits = 8)
  {
   return Compare(a, b, digits) >= 0;
  }
//+------------------------------------------------------------------+
//| Is Less Than or Equal.                                           |
//+------------------------------------------------------------------+
bool LTE(const double a, const double b, const int digits = 8)
  {
   return Compare(a, b, digits) <= 0;
  }
//+------------------------------------------------------------------+
//| Is (nearly) Zero, at the given decimal precision.                |
//+------------------------------------------------------------------+
bool IsZero(const double value, const int digits = 8)
  {
   return Compare(value, 0.0, digits) == 0;
  }
//+------------------------------------------------------------------+
//| Check almost equality (ignore tiny roundoff errors <= 1 ULP).    |
//| If binary representations of two doubles differ by more than     |
//| one least significant bit (ulp), the function returns false.     |
//| For example, AlmostEqual(0.1+0.2, 0.3) => true                   |
//+------------------------------------------------------------------+
bool AlmostEqual(const double a, const double b)
  {
// return (a == b) || MathAbs(a - b) <= DBL_EPSILON * MathMax(MathAbs(a), MathAbs(b));
   return (a == b) || MathAbs(UlpDiff(a, b)) <= 1;
  }
//+------------------------------------------------------------------+
//| Check whether two numbers match up to "n" significant digits.    |
//| For example, EqualDoubles(3.124, 3.122, 3) => true               |
//+------------------------------------------------------------------+
bool EqualDoubles(const double a, const double b, const int significantDigits = 15)
  {
//--- https://stackoverflow.com/a/17382806
   return (a == b) || MathAbs(a - b) <= GetPower10(-significantDigits) * MathMax(MathAbs(a), MathAbs(b));
  }
//+------------------------------------------------------------------+
//| Check whether two floating point numbers are close in value.     |
//| n: the max number of least significant digits that do not agree. |
//| Something like 2 or 3 usually works in practice.                 |
//+------------------------------------------------------------------+
bool IsClose(const double a, const double b, const int maxDifferentSignificantDigits = 2)
  {
//--- https://stackoverflow.com/a/17382806
   return (a == b) || MathAbs(a - b) <= GetPower10(maxDifferentSignificantDigits) * DBL_EPSILON * MathMax(MathAbs(a), MathAbs(b));
  }
//+------------------------------------------------------------------+
//| Get number of fractional digits that agree after decimal point.  |
//| For example, GetEqualDigits(3.124, 3.122) => 2                   |
//| Returns 0 if either input is NaN or Infinity.                    |
//+------------------------------------------------------------------+
int GetEqualDigits(const double a, const double b)
  {
   if(!MathIsValidNumber(a) || !MathIsValidNumber(b))
      return 0;

//--- https://stackoverflow.com/a/49242211
   return (a == b) ? 30 : MathMax(-FloorLog10(MathAbs(a - b)) - 1, 0);
  }
//+------------------------------------------------------------------+
//| Get number of significant digits that agree between two numbers. |
//| For example, GetEqualSigDigits(3.124, 3.122) => 3                |
//| Returns 0 if either input is NaN or Infinity.                    |
//+------------------------------------------------------------------+
int GetEqualSigDigits(const double a, const double b)
  {
   if(!MathIsValidNumber(a) || !MathIsValidNumber(b))
      return 0;

   return (a == b) ? 17 : FloorLog10(MathMax(MathAbs(a), MathAbs(b))) - FloorLog10(MathAbs(a - b));
  }


//+------------------------------------------------------------------+
//| II. Handy functions for rounding of doubles.                     |
//+------------------------------------------------------------------+

//+------------------------------------------------------------------+
//| Largest magnitude for which not every double is already an       |
//| integer. Beyond 2^52, every representable double has no          |
//| fractional part, so Ceil/Floor/Round/Trunc are no-ops there --   |
//| and casting a value that large to "long" would overflow anyway  |
//| (long only spans up to ~9.2e18, and the cast is undefined        |
//| behavior once the value falls outside that range).               |
//+------------------------------------------------------------------+
#define MU_MAX_FRACTIONAL_DOUBLE   9007199254740992.0   // 2^52

//+------------------------------------------------------------------+
//| Fast ceil, floor, and round using arithmetic operators.          |
//| NaN, Infinity, and values with |v| >= 2^52 are returned          |
//| unchanged (see MU_MAX_FRACTIONAL_DOUBLE above).                  |
//+------------------------------------------------------------------+
double Ceil(const double v)
  {
   if(!MathIsValidNumber(v) || MathAbs(v) >= MU_MAX_FRACTIONAL_DOUBLE)
      return v;
   const double k = (double)(long)v;
   return k + (v > k);
  }
double Floor(const double v)
  {
   if(!MathIsValidNumber(v) || MathAbs(v) >= MU_MAX_FRACTIONAL_DOUBLE)
      return v;
   const double k = (double)(long)v;
   return k - (v < k);
  }
double Round(const double v)
  {
   if(!MathIsValidNumber(v) || MathAbs(v) >= MU_MAX_FRACTIONAL_DOUBLE)
      return v;
   return (double)(long)((v < 0) ? v - 0.5 : v + 0.5);
  }
double Trunc(const double v)
  {
   if(!MathIsValidNumber(v) || MathAbs(v) >= MU_MAX_FRACTIONAL_DOUBLE)
      return v;
   return (double)(long)v;
  }

//+------------------------------------------------------------------+
//| Round to a specified number of decimal digits.                   |
//+------------------------------------------------------------------+
/*
 * https://stackoverflow.com/a/48764436/4208440
 *
 * A simple drop in solution that provides accurate decimal rounding.
 * It treats doubles more like decimals by fixing the binary rounding
 * issues to avoid unexpected results.
 */
double Ceil (const double value, const int digits)  { double p = GetPower10(digits); return Ceil ((value * p) * (1.0 - DBL_EPSILON * SignDbl(value))) / p; }
double Floor(const double value, const int digits)  { double p = GetPower10(digits); return Floor((value * p) * (1.0 + DBL_EPSILON * SignDbl(value))) / p; }
double Round(const double value, const int digits)  { double p = GetPower10(digits); return Round((value * p) * (1.0 + DBL_EPSILON)) / p; }
double Trunc(const double value, const int digits)  { double p = GetPower10(digits); return Trunc((value * p) * (1.0 + DBL_EPSILON)) / p; }

//+------------------------------------------------------------------+
//| Round to a whole multiple of some increment step.                |
//+------------------------------------------------------------------+
double Ceil (const double value, const double step)  { return (step <= 0) ? value : StripError(Ceil (StripError(value / step)) * step); }
double Floor(const double value, const double step)  { return (step <= 0) ? value : StripError(Floor(StripError(value / step)) * step); }
double Round(const double value, const double step)  { return (step <= 0) ? value : StripError(Round(StripError(value / step)) * step); }
double Trunc(const double value, const double step)  { return (step <= 0) ? value : StripError(Trunc(StripError(value / step)) * step); }

//+------------------------------------------------------------------+
//| Round to a specified number of significant figures.              |
//+------------------------------------------------------------------+
/**
* https://en.wikipedia.org/wiki/Significant_figures
* https://github.com/php/php-src/blob/PHP-7.4.10/ext/standard/math.c#L127-L206
* https://stackoverflow.com/a/374470/4208440
*/
double RoundToSignificantDigits(const double value, int digits)
  {
   if(!value || !MathIsValidNumber(value))
      return value;

//--- Round to significant digits
   digits -= FloorLog10(value) + 1;
   const double p = GetPower10(MathAbs(digits));
   return digits > 0 ? Round(value * p) / p : Round(value / p) * p;
  }
//+------------------------------------------------------------------+
//| Strip intermediate binary floating-point round-off errors.       |
//+------------------------------------------------------------------+
double StripError(const double value)
  {
//--- if the number is whole integer
   if(Floor(value) == value)
     {
      return value;
     }
//--- DBL_DIG: Number of significant decimal digits for double type = 15
   return RoundToSignificantDigits(value, DBL_DIG);
  }
//+------------------------------------------------------------------+
//| Round a calculated price to the nearest tick price.              |
//+------------------------------------------------------------------+
double RoundPrice(const double pPrice, const string pSymbol = NULL)
  {
   const string sym = (pSymbol == NULL) ? _Symbol : pSymbol;
   const double ticksize = SymbolInfoDouble(sym, SYMBOL_TRADE_TICK_SIZE);

   if(ticksize <= 0.0)
     {
      Print(__FUNCTION__, ": Warning, cannot retrieve ticksize for " + sym);
      return pPrice;
     }

   return Round(pPrice, ticksize);
  }
//+------------------------------------------------------------------+
//| Round a volume according to broker volume min/max/step limits.   |
//| Note: Returns 0.0 if volume is strictly below SYMBOL_VOLUME_MIN. |
//+------------------------------------------------------------------+
double RoundVolume(const double pVolume, const string pSymbol = NULL)
  {
   const string sym = (pSymbol == NULL) ? _Symbol : pSymbol;

   const double minlot  = SymbolInfoDouble(sym, SYMBOL_VOLUME_MIN);
   const double maxlot  = SymbolInfoDouble(sym, SYMBOL_VOLUME_MAX);
   const double lotstep = SymbolInfoDouble(sym, SYMBOL_VOLUME_STEP);

   if(minlot == 0.0 || maxlot == 0.0 || lotstep == 0.0)
     {
      Print(__FUNCTION__, ": Warning, cannot retrieve volume info for " + sym);
      return pVolume;
     }

//--- tolerance scales with lotstep instead of a fixed constant, which
//--- would be too coarse for instruments whose volume step is smaller
//--- (e.g. some crypto CFDs quote SYMBOL_VOLUME_STEP as 0.00000001)
   const double Epsilon = lotstep * 0.5;
   double volume = pVolume - minlot;

   if(volume < -Epsilon)
      return 0.0;

   if(volume < Epsilon)
      return minlot;

   volume = minlot + Round(volume, lotstep);

   if(volume > maxlot)
      volume = maxlot;

   return volume;
  }


//+------------------------------------------------------------------+
//| III. Miscellaneous handy functions for doubles.                  |
//+------------------------------------------------------------------+

//+------------------------------------------------------------------+
//| Determines whether the passed value is an integer.               |
//+------------------------------------------------------------------+
bool IsInteger(const double value)
  {
// return (MathMod(value, 1.0) == 0)
   return Round(value) == value;
  }
//+------------------------------------------------------------------+
//| Checks if a number has a specified number of decimal digits.     |
//+------------------------------------------------------------------+
bool IsRound(const double value, const int digits)
  {
// return value == StringToDouble(DoubleToString(value, digits));
   const double p = GetPower10(digits);
   return MathRound(value * p) / p == value;
  }
//+------------------------------------------------------------------+
//| Checks if a number is a whole multiple of some increment.        |
//+------------------------------------------------------------------+
bool IsRound(const double value, const double step)
  {
   return Round(value, step) == value;
  }
//+------------------------------------------------------------------+
//| Get the decimal fractional part (retains sign).                  |
//+------------------------------------------------------------------+
double Frac(const double value)
  {
   return MathMod(value, 1.0);
  }
//+------------------------------------------------------------------+
//| Get number of decimal digits after the decimal point.            |
//+------------------------------------------------------------------+
int GetDigits(const double value)
  {
   int d = 0;
// for(double p = 1; value != MathRound(value * p) / p && (d < 308); p *= 10)
   while(!IsRound(value, d) && d < DBL_MAX_10_EXP)
      d++;
   return d;
  }
//+------------------------------------------------------------------+
//| Get number of integer digits to the left of decimal point.       |
//| Returns 0 for NaN or Infinity.                                   |
//+------------------------------------------------------------------+
int GetIntegerDigits(const double value)
  {
   if(!MathIsValidNumber(value) || value == 0.0 || MathAbs(value) < 1.0)
      return 0;

   return FloorLog10(value) + 1;
  }
//+------------------------------------------------------------------+
//| Get number of significant digits. Significant digits or figures  |
//| is the sum of integer and decimal digits (left and right of the  |
//| decimal point), excluding leading and trailing zeros.            |
//| For example, the number 1.23 has 2 decimal places and 3 s.f.     |
//| Hint: Change number to scientific notation. It is easier to see. |
//+------------------------------------------------------------------+
int GetSignificantDigits(double value)
  {
   if(!value || !MathIsValidNumber(value))
      return 0;

//--- sum of decimal and integral digits
   int digits = GetDigits(value) + FloorLog10(value) + 1;

//--- excluding trailing zeros
   while(MathMod(value, 10.0) == 0.0)
     {
      value /= 10.0;
      digits--;
     }

   return digits;
  }
//+------------------------------------------------------------------+
//| Linear Transformation (RESCALE)                                  |
//| Convert x in range [min, max] to y in another range [destMin,    |
//| destMax]. This can be used to rescale indicator values or        |
//| graphic objects. Common scales are [0,1], [0,100] or [-1,1].     |
//+------------------------------------------------------------------+
double Remap(const double value, const double min, const double max, const double destMin, const double destMax)
  {
//--- to maintain ratios
   double y = safeDiv((value - min), (max - min)) * (destMax - destMin) + destMin;
   return Clamp(y, destMin, destMax);
  }
//+------------------------------------------------------------------+
//| Linear Interpolation (Lerp) between start and end by factor t.   |
//+------------------------------------------------------------------+
double Lerp(const double start, const double end, const double t)
  {
   return start + t * (end - start);
  }

//+------------------------------------------------------------------+
//| Remap value from scale [min, max] to normalized scale [0, 1].    |
//+------------------------------------------------------------------+
double Normalize(const double value, const double min, const double max)
  {
   return Remap(value, min, max, 0.0, 1.0);
  }
//+------------------------------------------------------------------+
//| Remap normalized value [0, 1] to a scale [destMin, destMax].     |
//+------------------------------------------------------------------+
double Denormalize(const double value, const double destMin, const double destMax)
  {
   return Remap(value, 0, 1, destMin, destMax);
  }
//+------------------------------------------------------------------+
//| Clamp a value into closed interval [min, max].                   |
//+------------------------------------------------------------------+
template<typename T>
T Clamp(const T value, const T min, const T max)
  {
   //if(value < min) return min;
   //if(value > max) return max;
   //return value;
   T lo = MathMin(min, max);
   T hi = MathMax(min, max);
   return MathMin(MathMax(value, lo), hi);
  }
//+------------------------------------------------------------------+
//| Clamp a value to a maximum upper threshold.                      |
//+------------------------------------------------------------------+
template<typename T>
T ClampMax(const T value, const T max)
  {
   return MathMin(value, max);
  }
//+------------------------------------------------------------------+
//| Clamp a value to a minimum lower threshold.                      |
//+------------------------------------------------------------------+
template<typename T>
T ClampMin(const T value, const T min)
  {
   return MathMax(value, min);
  }
//+------------------------------------------------------------------+
//| Wrap a value into the inclusive interval [min, max].             |
//| Requires max >= min; returns min unchanged otherwise.            |
//+------------------------------------------------------------------+
// https://stackoverflow.com/a/14416133
template<typename T>
T Wrap(const T value, const T min, const T max)
  {
   if(max < min)
      return min;

   const T range = max - min + 1;
   //while(value < min) value += range;
   //while(value > max) value -= range;
   //return value;
   if(value < min)
      return max - (T)MathMod((max - value), range);
   if(value > max)
      return min + (T)MathMod((value - min), range);
   return value;
  }

// Overload for double
double Wrap(const double value, const double min, const double max)
  {
   const double range = max - min;
   if(range <= 0.0)
      return min;

   return min + MathMod(MathMod((value - min), range) + range, range);
  }
//+------------------------------------------------------------------+
//| Division helper preventing zero-division errors.                 |
//| Returns zero if denominator is zero.                             |
//+------------------------------------------------------------------+
double safeDiv(const double a, const double b)
  {
   if(b == 0.0)
      return 0.0;

   return a / b;
  }
//+------------------------------------------------------------------+
//| Fast lookup table for exact powers of 10.                        |
//+------------------------------------------------------------------+
const double PowersOf10[] =
  {
   1e00, 1e01, 1e02, 1e03, 1e04, 1e05, 1e06, 1e07,
   1e08, 1e09, 1e10, 1e11, 1e12, 1e13, 1e14, 1e15,
   1e16, 1e17, 1e18, 1e19, 1e20, 1e21, 1e22
  };
//+------------------------------------------------------------------+
//| Fast base-10 exponent solver (floor(log10(value))).              |
//|                                                                  |
//| Returns INT_MIN as a sentinel for value == 0.0, NaN, or Infinity |
//| (log10 is undefined for these).                                  |
//|                                                                  |
//| Uses the relationship log10(x) = log2(x) * log10(2)              |
//| http://graphics.stanford.edu/~seander/bithacks.html#IntegerLog10 |
//+------------------------------------------------------------------+
int FloorLog10(double value)
  {
   if(value == 0.0 || !MathIsValidNumber(value))
      return INT_MIN;

   int result = 0;
   value = MathAbs(value);

   /* Not in lookup table */
   if(value < 1e-22 || value > 1e22)
     {
      result = (int)MathFloor(MathLog10(value));
     }
   else
     {
      const int e = GetExponent(value);           /* floor(log2(x)) */
      const int t = ((e + 1) * 20201781) >> 26;   /* log10(2) ~= 20201781/2^26 */
      if(t >= 0)
         result = t - (value < PowersOf10[t]);    /* t may be off by one */
      else
         result = t - (value < 1.0 / PowersOf10[-t]);
     }

   return result;
  }
//+------------------------------------------------------------------+
//| Fast integer powers of 10 lookup routine.                        |
//| https://github.com/php/php-src/blob/master/ext/standard/math     |
//+------------------------------------------------------------------+
double GetPower10(const int power)
  {
   /* Not in lookup table */
   if(power < -22 || power > 22)
      return MathPow(10.0, (double)power);

   if(power >= 0)
      return PowersOf10[power];

   return 1.0 / PowersOf10[-power];
  }
//+------------------------------------------------------------------+
//| Computes the sign of a value as 1, 0, or -1.                     |
//+------------------------------------------------------------------+
int SignDbl(const double value)
  {
   return (value > 0.0) - (value < 0.0);
  }


//+------------------------------------------------------------------+
//| IV. Handy functions for low-level binary operations on doubles.  |
//+------------------------------------------------------------------+

//+------------------------------------------------------------------+
//| Returns the bit representation corresponding to a double value . |
//+------------------------------------------------------------------+
long DoubleToLongBits(const double value)
  {
   union _d {double value; long bits;} dbl;
   dbl.value = value;
   return dbl.bits;
  }
//+------------------------------------------------------------------+
//| Returns the double value corresponding to a bit representation.  |
//+------------------------------------------------------------------+
double LongBitsToDouble(const long bits)
  {
   union _d {double value; long bits;} dbl;
   dbl.bits = bits;
   return dbl.value;
  }
//+------------------------------------------------------------------+
//| Converts the IEEE-754 bit representation of a double into an     |
//| integer ordering that preserves the numeric ordering of doubles. |
//+------------------------------------------------------------------+
long DoubleToOrderedLong(const double value)
  {
   long bits = DoubleToLongBits(value);
   if(bits < 0)
      bits = MU_SIGN_BIT - bits;
   return bits;
  }
//+------------------------------------------------------------------+
//| Converts an ordered integer representation produced by           |
//| DoubleToOrderedLong() back into the original IEEE-754 double.    |
//+------------------------------------------------------------------+
double OrderedLongToDouble(long bits)
  {
   if(bits < 0)
      bits = MU_SIGN_BIT - bits;
   return LongBitsToDouble(bits);
  }
//+------------------------------------------------------------------+
//| Returns the raw encoding of exponent in the bit representation.  |
//|                                                                  |
//| Valid only for finite values.                                    |
//| Results for NaN and ±Infinity are unspecified.                   |
//+------------------------------------------------------------------+
int RawExponent(const double value)
  {
   return (int)((DoubleToLongBits(value) >> MU_MANTISSA_BITS) & MU_EXPONENT_MASK);
  }
//+------------------------------------------------------------------+
//| Returns raw encoding of significand in the bit representation.   |
//|                                                                  |
//| Valid only for finite values.                                    |
//| Results for NaN and ±Infinity are unspecified.                   |
//+------------------------------------------------------------------+
long RawSignificand(const double value)
  {
   return (long)(DoubleToLongBits(value) & MU_MANTISSA_MASK);
  }
//+------------------------------------------------------------------+
//| Returns the unbiased (adjusted) exponent of a double value.      |
//|                                                                  |
//| Valid only for finite, normal values.                            |
//| Results for NaN, ±Infinity, and subnormal values are unspecified.|
//+------------------------------------------------------------------+
int GetExponent(const double value)
  {
// return (int) MathFloor(MathLog(MathAbs(value)) / M_LN2);
   return RawExponent(value) - MU_EXPONENT_BIAS;
  }
//+------------------------------------------------------------------+
//| Returns the significand of a double value. For finite nonzero    |
//| numbers, the significand is in the range 1.0 ..< 2.0.            |
//| An IEEE-754 double number = (Sign) 2 ^ Exponent * Significand.   |
//|                                                                  |
//| Valid only for finite, normal values.                            |
//| Results for NaN, ±Infinity, and subnormal values are unspecified.|
//+------------------------------------------------------------------+
double GetSignificand(const double value)
  {
// return 1.0 + RawSignificand(value) / MathPow(2, 52);
   return LongBitsToDouble(RawSignificand(value) | 0x3FF0000000000000);  // 0x3FF0... = exponent bits for 2^0
  }
//+------------------------------------------------------------------+
//| Determine whether number is exactly representable in double.     |
//| i.e., No rounding to an approximation during the conversion.     |
//| Results are valid for numbers in the range [2^-24, 2^52].        |
//| Returns false for NaN and Infinity.                              |
//+------------------------------------------------------------------+
//https://stackoverflow.com/a/67952907/4208440
bool IsExactDouble(const double value)
  {
   if(!MathIsValidNumber(value))
      return false;

   if(value == 0.0)
      return true;

   const int exp2 = GetExponent(value);
   const int exp10 = FloorLog10(value);

//--- check for any mismatch between the exact decimal and
//--- the round-trip representation.
   const int rightmost_bits = (52 - exp2) - (16 - exp10);
//--- guard against undefined shift behavior at the boundaries
   if(rightmost_bits <= 0)
      return true;                                 // integer or exact power-of-two
   if(rightmost_bits >= 63)
      return false;                                // outside [2^-24, 2^52] range

//--- create bitmask for rightmost bits
   const ulong mask = ((ulong)1 << rightmost_bits) - 1;

//--- test if all rightmost bits are 0's (i.e., no rounding)
   return (DoubleToLongBits(value) & mask) == 0;
  }
//+------------------------------------------------------------------+
//| Returns the next representable value after x away from zero.     |
//|                                                                  |
//| Valid for all finite values, including zero and subnormals.      |
//| NaN is returned unchanged.                                       |
//| Advancing beyond ±DBL_MAX returns ±Infinity.                     |
//+------------------------------------------------------------------+
double NextAwayFromZero(const double value)
  {
   return DoubleAdvance(value, value >= 0.0 ? 1 : -1);
  }
//+------------------------------------------------------------------+
//| Advances a floating-point number by a specified number of ULPs.  |
//| distance may be negative.                                        |
//|                                                                  |
//| Valid for all finite values, including zero and subnormals.      |
//| NaN is returned unchanged.                                       |
//| Advancing beyond ±DBL_MAX returns ±Infinity.                     |
//+------------------------------------------------------------------+
double DoubleAdvance(const double value, const long distance)
  {
   if(!MathIsValidNumber(value))
      return value;

   const long bits = DoubleToOrderedLong(value);

   return OrderedLongToDouble(bits + distance);
  }
//+------------------------------------------------------------------+
//| Returns the size of a unit in the last place (machine epsilon)   |
//| for a specified double value. This is the unit of the least      |
//| significant digit in this value’s significand. For most numbers, |
//| this is the positive distance between this double value and the  |
//| representable value next larger in magnitude.                    |
//| Note that, Ulp(1.0) == DBL_EPSILON.  (i.e., epsilon at 1.0)      |
//+------------------------------------------------------------------+
double Ulp(const double value)
  {
   if(!MathIsValidNumber(value))
      return value;

   return MathAbs(NextAwayFromZero(value) - value);
  }
//+------------------------------------------------------------------+
//| Returns the distance between two doubles a and b expressed as    |
//| the number of gaps/bits/ULP (Units in the Last Place) between a  |
//| and b. Note that the function returns a signed value indicating  |
//| whether a > b or not.                                            |
//+------------------------------------------------------------------+
long UlpDiff(const double a, const double b)
  {
   return DoubleToOrderedLong(a) -
          DoubleToOrderedLong(b);
  }

//+------------------------------------------------------------------+
//| V. Handy functions for formatting of doubles to string.          |
//+------------------------------------------------------------------+

//+------------------------------------------------------------------+
//| Converting numeric value into the raw hexadecimal text string.   |
//| https://baseconvert.com/ieee-754-floating-point                  |
//+------------------------------------------------------------------+
string DoubleToHexadecimal(const double value)
  {
   return StringFormat("0x%.16I64X", DoubleToLongBits(value));
  }
//+------------------------------------------------------------------+
//| Converting numeric value into the hex float constant string.     |
//| A real number in format [−]0xh.hhhh P±dd, where h.hhhh – mantissa|
//| in the form of hexadecimal digits, using "ABCDEF", dd - One or   |
//| more digits of exponent.                                         |
//+------------------------------------------------------------------+
string DoubleToHexFloatConstant(const double value)
  {
   return StringFormat("%.13a", value);
  }
//+------------------------------------------------------------------+
//| Converting numeric value into the exact decimal text string of   |
//| the binary approximation stored by the machine.                  |
//| https://baseconvert.com/ieee-754-floating-point                  |
//+------------------------------------------------------------------+
string DoubleToStringExact(const double value)
  {
   if(!value || !MathIsValidNumber(value))
      return (string)value;

//--- maximum number of decimal places = number of fraction bits
   const int exponent = (int)MathFloor(MathLog(MathAbs(value)) / M_LN2);
   const int digits   = MathMax(0, 52 - exponent);

//--- https://www.exploringbinary.com/number-of-decimal-digits-in-a-binary-fraction/
#ifdef __MQL4__
   return StringFormat("%." + (string)digits + "f", value);
#else
   return StringFormat("%.*f", digits, value);
#endif
  }
//+------------------------------------------------------------------+
//| Converting numeric value into a string in scientific notation    |
//| with one digit before the decimal point (e.g., 6.22e-23).        |
//| Digits : Optional. The number of digits after the decimal point. |
//| Defaults to as many digits as necessary to represent the value.  |
//+------------------------------------------------------------------+
string DoubleToExponential(const double value, int digits = -1)
  {
   if(!value || !MathIsValidNumber(value) || digits > 20)
      return (string)value;

//--- as many digits as necessary to represent the value
//--- https://www.w3schools.com/jsref/jsref_toexponential.asp
   if(digits < 0)
      digits = GetSignificantDigits(value) - 1;
//---
#ifdef __MQL4__
   return StringFormat("%." + (string)digits + "e", value);
#else
   return StringFormat("%.*e", digits, value);
#endif
  }
//+------------------------------------------------------------------+
//| Converting numeric value into the shortest string representation |
//| that round-trips into the same numeric value. The result will    |
//| contain at most 17 significant digits, discarding trailing zeros.|
//| The round-trip ("%.17g") format specifier ensures that a numeric |
//| value converted to a string is always parsed back into the same  |
//| numeric value, StringToDouble(Repr(f)) == f.                     |
//| Note: results are consistent with David M. Gay's dtoa.c library. |
//+------------------------------------------------------------------+
// toString()
string Repr(const double value)
  {
//--- https://stackoverflow.com/a/35708911/4208440
//--- https://www.exploringbinary.com/number-of-digits-required-for-round-trip-conversions/
//--- Try format with 15, 16, 17 significant digits to return the shortest
//--- decimal numeric string which round-trips to the specified value.
   string str = NULL;
   for(int sig = 15; sig <= 17; sig++)
     {
#ifdef __MQL4__
      if(value == StringToDouble(str = StringFormat("%." + (string)sig + "g", value)))
         break;
#else
      if(value == StringToDouble(str = StringFormat("%.*g", sig, value)))
         break;
#endif
     }
   return str;
  }
//+------------------------------------------------------------------+
//| Converting float value into the shortest string representation   |
//| that round-trips into the same numeric value. This ensures that  |
//| a float value converted to a string is always parsed back into   |
//| the same numeric value. i.e., (float)(Repr(f)) == f.             |
//+------------------------------------------------------------------+
// toString()
string Repr(const float value)
  {
//--- Try format with 6, 7, 8, 9 significant digits to return the shortest
//--- decimal numeric string which round-trips to the specified value.
   string str = NULL;
   for(int sig = 6; sig <= 9; sig++)
     {
#ifdef __MQL4__
      if(value == float(str = StringFormat("%." + (string)sig + "g", value)))
         break;
#else
      if(value == float(str = StringFormat("%.*g", sig, value)))
         break;
#endif
     }
   return str;
  }
//+------------------------------------------------------------------+
//| Make a diagnostic dump of a double for use in Print()/debugging. |
//| Combines Repr, hex bit pattern, exponent, and significand into   |
//| a single readable string, e.g.:                                  |
//|   DoubleDebugInfo(0.1)                                           |
//|   => "0.1 [hex=0x3FB999999999999A, exp=-4, significand=1.6]"     |
//+------------------------------------------------------------------+
string DoubleDebugInfo(const double value)
  {
   if(!MathIsValidNumber(value))
      return Repr(value) + " [non-finite]";

   if(value == 0.0)
      return Repr(value) + " [hex=" + DoubleToHexadecimal(value) + "]";

   return StringFormat("%s [hex=%s, exp=%d, significand=%s]",
                        Repr(value),
                        DoubleToHexadecimal(value),
                        GetExponent(value),
                        Repr(GetSignificand(value)));
  }
//+------------------------------------------------------------------+
//| Formats double with thousands separator and specified decimals.  |
//| Note: digits can be only 0-8 decimal places.                     |
//+------------------------------------------------------------------+
/**
 * Alert( FormatDouble(balance, 2) );
 */
string FormatDouble(const double value, const int digits, const string separator=",")
  {
   const double numb = MathAbs(value);
   const string sign = (value < 0.0) ? "-" : "";
   string num_str = DoubleToString(NormalizeDouble(numb, digits), digits);  // FormatDouble(1.005, 2) => "1.01"

//--- Find the length of the integer part
   int pos = StringFind(num_str, ".");
   int length = (pos > -1) ? pos : StringLen(num_str);

//--- Format the integer part with thousand separators
   for(int i = length - 3; i > 0; i -= 3)
     {
      // insert a separator at position i
      string tmp = StringSubstr(num_str, 0, i);
      tmp += separator;
      num_str = tmp + StringSubstr(num_str, i);
     }

   return sign + num_str;
  }


#endif // #ifndef MATH_UTILS_UNIQUE_HEADER_ID_H