// Written in the D programming language. /**This is a submodule of $(MREF std, math). It contains classical algebraic functions like `abs`, `sqrt`, and `poly`. Copyright: Copyright The D Language Foundation 2000 - 2011.License:   $(HTTP www.boost.org/LICENSE_1_0.txt, Boost License 1.0).Authors:   $(HTTP digitalmars.com, Walter Bright), Don Clugston,           Conversion of CEPHES math library to D by Iain Buclaw and David NadlingerSource: $(PHOBOSSRC std/math/algebraic.d) Macros:    TABLE_SV = <table border="1" cellpadding="4" cellspacing="0">               <caption>Special Values</caption>               $0</table>    NAN = $(RED NAN)    POWER = $1<sup>$2</sup>    SUB = $1<sub>$2</sub>    PLUSMN = &plusmn;    INFIN = &infin;    PLUSMNINF = &plusmn;&infin;    LT = &lt;  */ module std.math.algebraic; static import core.math;static import core.stdc.math;import std.traits : CommonType, isFloatingPoint, isIntegral, isSigned, Unqual; /*********************************** * Calculates the absolute value of a number. * * Params: *     Num = (template parameter) type of number *       x = real number value * * Returns: *     The absolute value of the number. If floating-point or integral, *     the return type will be the same as the input. * * Limitations: *     When x is a signed integral equal to `Num.min` the value of x will be returned instead. *     Note for 2's complement; `-Num.min` (= `Num.max + 1`) is not representable due to overflow. */auto abs(Num)(Num x) @nogc nothrow pureif (isIntegral!Num || (is(typeof(Num.init >= 0)) && is(typeof(-Num.init)))){    static if (isFloatingPoint!(Num))        return fabs(x);    else    {        static if (isIntegral!Num)            return x >= 0 ? x : cast(Num) -x;        else            return x >= 0 ? x : -x;    }} ///@safe pure nothrow @nogc unittest{    import std.math.traits : isIdentical, isNaN;     assert(isIdentical(abs(-0.0L), 0.0L));    assert(isNaN(abs(real.nan)));    assert(abs(-real.infinity) == real.infinity);    assert(abs(-56) == 56);    assert(abs(2321312L)  == 2321312L);    assert(abs(23u) == 23u);} @safe pure nothrow @nogc unittest{    assert(abs(byte(-8)) == 8);    assert(abs(ubyte(8u)) == 8);    assert(abs(short(-8)) == 8);    assert(abs(ushort(8u)) == 8);    assert(abs(int(-8)) == 8);    assert(abs(uint(8u)) == 8);    assert(abs(long(-8)) == 8);    assert(abs(ulong(8u)) == 8);    assert(is(typeof(abs(byte(-8))) == byte));    assert(is(typeof(abs(ubyte(8u))) == ubyte));    assert(is(typeof(abs(short(-8))) == short));    assert(is(typeof(abs(ushort(8u))) == ushort));    assert(is(typeof(abs(int(-8))) == int));    assert(is(typeof(abs(uint(8u))) == uint));    assert(is(typeof(abs(long(-8))) == long));    assert(is(typeof(abs(ulong(8u))) == ulong));} @safe pure nothrow @nogc unittest{    import std.meta : AliasSeq;    static foreach (T; AliasSeq!(float, double, real))    {{        T f = 3;        assert(abs(f) == f);        assert(abs(-f) == f);    }}} // see https://issues.dlang.org/show_bug.cgi?id=20205// to avoid falling into the trap again@safe pure nothrow @nogc unittest{    assert(50 - abs(-100) == -50);} // https://issues.dlang.org/show_bug.cgi?id=19162@safe unittest{    struct Vector(T, int size)    {        T x, y, z;    }     static auto abs(T, int size)(auto ref const Vector!(T, size) v)    {        return v;    }    Vector!(int, 3) v;    assert(abs(v) == v);} /******************************* * Returns |x| * *      $(TABLE_SV *      $(TR $(TH x)                 $(TH fabs(x))) *      $(TR $(TD $(PLUSMN)0.0)      $(TD +0.0) ) *      $(TR $(TD $(PLUSMN)$(INFIN)) $(TD +$(INFIN)) ) *      ) */pragma(inline, true)real fabs(real x) @safe pure nothrow @nogc { return core.math.fabs(x); } ///dittopragma(inline, true)double fabs(double x) @safe pure nothrow @nogc { return core.math.fabs(x); } ///dittopragma(inline, true)float fabs(float x) @safe pure nothrow @nogc { return core.math.fabs(x); } ///@safe unittest{    import std.math.traits : isIdentical;     assert(isIdentical(fabs(0.0f), 0.0f));    assert(isIdentical(fabs(-0.0f), 0.0f));    assert(fabs(-10.0f) == 10.0f);     assert(isIdentical(fabs(0.0), 0.0));    assert(isIdentical(fabs(-0.0), 0.0));    assert(fabs(-10.0) == 10.0);     assert(isIdentical(fabs(0.0L), 0.0L));    assert(isIdentical(fabs(-0.0L), 0.0L));    assert(fabs(-10.0L) == 10.0L);} @safe unittest{    real function(real) pfabs = &fabs;    assert(pfabs != null);} @safe pure nothrow @nogc unittest{    float f = fabs(-2.0f);    assert(f == 2);     double d = fabs(-2.0);    assert(d == 2);     real r = fabs(-2.0L);    assert(r == 2);} /*************************************** * Compute square root of x. * *      $(TABLE_SV *      $(TR $(TH x)         $(TH sqrt(x))   $(TH invalid?)) *      $(TR $(TD -0.0)      $(TD -0.0)      $(TD no)) *      $(TR $(TD $(LT)0.0)  $(TD $(NAN))    $(TD yes)) *      $(TR $(TD +$(INFIN)) $(TD +$(INFIN)) $(TD no)) *      ) */pragma(inline, true)float sqrt(float x) @nogc @safe pure nothrow { return core.math.sqrt(x); } /// dittopragma(inline, true)double sqrt(double x) @nogc @safe pure nothrow { return core.math.sqrt(x); } /// dittopragma(inline, true)real sqrt(real x) @nogc @safe pure nothrow { return core.math.sqrt(x); } ///@safe pure nothrow @nogc unittest{    import std.math.operations : feqrel;    import std.math.traits : isNaN;     assert(sqrt(2.0).feqrel(1.4142) > 16);    assert(sqrt(9.0).feqrel(3.0) > 16);     assert(isNaN(sqrt(-1.0f)));    assert(isNaN(sqrt(-1.0)));    assert(isNaN(sqrt(-1.0L)));} @safe unittest{    // https://issues.dlang.org/show_bug.cgi?id=5305    float function(float) psqrtf = &sqrt;    assert(psqrtf != null);    double function(double) psqrtd = &sqrt;    assert(psqrtd != null);    real function(real) psqrtr = &sqrt;    assert(psqrtr != null);     //ctfe    enum ZX80 = sqrt(7.0f);    enum ZX81 = sqrt(7.0);    enum ZX82 = sqrt(7.0L);} @safe pure nothrow @nogc unittest{    float f = sqrt(2.0f);    assert(fabs(f * f - 2.0f) < .00001);     double d = sqrt(2.0);    assert(fabs(d * d - 2.0) < .00001);     real r = sqrt(2.0L);    assert(fabs(r * r - 2.0) < .00001);} /*************** * Calculates the cube root of x. * *      $(TABLE_SV *      $(TR $(TH $(I x))            $(TH cbrt(x))           $(TH invalid?)) *      $(TR $(TD $(PLUSMN)0.0)      $(TD $(PLUSMN)0.0)      $(TD no) ) *      $(TR $(TD $(NAN))            $(TD $(NAN))            $(TD yes) ) *      $(TR $(TD $(PLUSMN)$(INFIN)) $(TD $(PLUSMN)$(INFIN)) $(TD no) ) *      ) */real cbrt(real x) @trusted pure nothrow @nogc{    version (CRuntime_Microsoft)    {        import std.math.traits : copysign;        import std.math.exponential : exp2;         version (INLINE_YL2X)            return copysign(exp2(core.math.yl2x(fabs(x), 1.0L/3.0L)), x);        else            return core.stdc.math.cbrtl(x);    }    else        return core.stdc.math.cbrtl(x);} ///@safe pure unittest{    import std.math.operations : feqrel;     assert(cbrt(1.0).feqrel(1.0) > 16);    assert(cbrt(27.0).feqrel(3.0) > 16);    assert(cbrt(15.625).feqrel(2.5) > 16);} /*********************************************************************** * Calculates the length of the * hypotenuse of a right-angled triangle with sides of length x and y. * The hypotenuse is the value of the square root of * the sums of the squares of x and y: * *      sqrt($(POWER x, 2) + $(POWER y, 2)) * * Note that hypot(x, y), hypot(y, x) and * hypot(x, -y) are equivalent. * *  $(TABLE_SV *  $(TR $(TH x)            $(TH y)                 $(TH hypot(x, y)) $(TH invalid?)) *  $(TR $(TD x)            $(TD $(PLUSMN)0.0)      $(TD |x|)         $(TD no)) *  $(TR $(TD $(PLUSMNINF)) $(TD y)                 $(TD +$(INFIN))   $(TD no)) *  $(TR $(TD $(PLUSMNINF)) $(TD $(NAN))            $(TD +$(INFIN))   $(TD no)) *  $(TR $(TD $(NAN))       $(TD y != $(PLUSMNINF)) $(TD $(NAN))      $(TD no)) *  ) */T hypot(T)(const T x, const T y) @safe pure nothrow @nogcif (isFloatingPoint!T){    // Scale x and y to avoid underflow and overflow.    // If one is huge and the other tiny, return the larger.    // If both are huge, avoid overflow by scaling by 2^^-N.    // If both are tiny, avoid underflow by scaling by 2^^N.    import core.math : fabs, sqrt;    import std.math.traits : floatTraits, RealFormat, isNaN;     alias F = floatTraits!T;     T u = fabs(x);    T v = fabs(y);    if (!(u >= v))    {        v = u;        u = fabs(y);        if (u == T.infinity) return u; // hypot(inf, nan) == inf        if (v == T.infinity) return v; // hypot(nan, inf) == inf        if (u.isNaN || v.isNaN)            return T.nan;    }    assert(!(u.isNaN || v.isNaN), "Comparison to NaN always fails, thus is is always handled in the branch above");     static if (F.realFormat == RealFormat.ieeeSingle)    {        enum SQRTMIN = 0x1p-60f;        enum SQRTMAX = 0x1p+60f;        enum SCALE_UNDERFLOW = 0x1p+90f;        enum SCALE_OVERFLOW = 0x1p-90f;    }    else static if (F.realFormat == RealFormat.ieeeDouble ||                    F.realFormat == RealFormat.ieeeExtended53 ||                    F.realFormat == RealFormat.ibmExtended)    {        enum SQRTMIN = 0x1p-450L;        enum SQRTMAX = 0x1p+500L;        enum SCALE_UNDERFLOW = 0x1p+600L;        enum SCALE_OVERFLOW = 0x1p-600L;    }    else static if (F.realFormat == RealFormat.ieeeExtended ||                    F.realFormat == RealFormat.ieeeQuadruple)    {        enum SQRTMIN = 0x1p-8000L;        enum SQRTMAX = 0x1p+8000L;        enum SCALE_UNDERFLOW = 0x1p+10000L;        enum SCALE_OVERFLOW = 0x1p-10000L;    }    else        assert(0, "hypot not implemented");     if (u * T.epsilon > v)    {        // hypot (huge, tiny) = huge        // also: hypot(x, 0) = x        return u;    }     // Now u >= v, or else one is NaN.    T ratio = 1.0;    if (v >= SQRTMAX)    {        // hypot(huge, huge) -- avoid overflow        ratio = SCALE_UNDERFLOW;        u *= SCALE_OVERFLOW;        v *= SCALE_OVERFLOW;    }    else if (u <= SQRTMIN)    {        // hypot (tiny, tiny) -- avoid underflow        // This is only necessary to avoid setting the underflow        // flag.        ratio = SCALE_OVERFLOW;        u *= SCALE_UNDERFLOW;        v *= SCALE_UNDERFLOW;    }     // both are in the normal range    return ratio * sqrt(u*u + v*v);} ///@safe unittest{    import std.math.operations : feqrel;    import std.math.traits : isNaN;     assert(hypot(1.0, 1.0).feqrel(1.4142) > 16);    assert(hypot(3.0, 4.0).feqrel(5.0) > 16);    assert(hypot(real.infinity, 1.0L) == real.infinity);    assert(hypot(1.0L, real.infinity) == real.infinity);    assert(hypot(real.infinity, real.nan) == real.infinity);    assert(hypot(real.nan, real.infinity) == real.infinity);    assert(hypot(real.nan, 1.0L).isNaN);    assert(hypot(1.0L, real.nan).isNaN);} @safe unittest{    import std.math.operations : feqrel;     assert(hypot(1.0f, 1.0f).feqrel(1.4142f) > 16);    assert(hypot(3.0f, 4.0f).feqrel(5.0f) > 16);    assert(hypot(float.infinity, 1.0f) == float.infinity);    assert(hypot(float.infinity, float.nan) == float.infinity);     assert(hypot(1.0L, 1.0L).feqrel(1.4142L) > 16);    assert(hypot(3.0L, 4.0L).feqrel(5.0L) > 16);    assert(hypot(double.infinity, 1.0) == double.infinity);    assert(hypot(double.infinity, double.nan) == double.infinity);} // https://github.com/dlang/phobos/issues/10491@safe pure nothrow unittest{    import std.math : isClose;     enum small = 5.016556e-20f;    assert(hypot(small, 0).isClose(small));    assert(hypot(small, float.min_normal).isClose(small));} @safe unittest{    import std.math.operations : feqrel;    import std.math.traits : isIdentical;    import std.meta : AliasSeq;     static foreach (T; AliasSeq!(float, double, real))    {{        static T[3][] vals =     // x,y,hypot        [            [ 0.0,     0.0,   0.0],            [ 0.0,    -0.0,   0.0],            [ -0.0,   -0.0,   0.0],            [ 3.0,     4.0,   5.0],            [ -300,   -400,   500],            [0.0,      7.0,   7.0],            [9.0,   9*T.epsilon,   9.0],            [88/(64*sqrt(T.min_normal)), 105/(64*sqrt(T.min_normal)), 137/(64*sqrt(T.min_normal))],            [88/(128*sqrt(T.min_normal)), 105/(128*sqrt(T.min_normal)), 137/(128*sqrt(T.min_normal))],            [3*T.min_normal*T.epsilon, 4*T.min_normal*T.epsilon, 5*T.min_normal*T.epsilon],            [ T.min_normal, T.min_normal, sqrt(2.0L)*T.min_normal],            [ T.max/sqrt(2.0L), T.max/sqrt(2.0L), T.max],            [ T.infinity, T.nan, T.infinity],            [ T.nan, T.infinity, T.infinity],            [ T.nan, T.nan, T.nan],            [ T.nan, T.max, T.nan],            [ T.max, T.nan, T.nan],        ];        for (int i = 0; i < vals.length; i++)        {            T x = vals[i][0];            T y = vals[i][1];            T z = vals[i][2];            T h = hypot(x, y);            assert(isIdentical(z,h) || feqrel(z, h) >= T.mant_dig - 1);        }     }}} /*********************************************************************** * Calculates the distance of the point (x, y, z) from the origin (0, 0, 0) * in three-dimensional space. * The distance is the value of the square root of the sums of the squares * of x, y, and z: * *      sqrt($(POWER x, 2) + $(POWER y, 2) + $(POWER z, 2)) * * Note that the distance between two points (x1, y1, z1) and (x2, y2, z2) * in three-dimensional space can be calculated as hypot(x2-x1, y2-y1, z2-z1). * * Params: *     x = floating point value *     y = floating point value *     z = floating point value * * Returns: *     The square root of the sum of the squares of the given arguments. */T hypot(T)(const T x, const T y, const T z) @safe pure nothrow @nogcif (isFloatingPoint!T){    import core.math : fabs, sqrt;    import std.math.operations : fmax;    const absx = fabs(x);    const absy = fabs(y);    const absz = fabs(z);     // Scale all parameters to avoid overflow.    const ratio = fmax(absx, fmax(absy, absz));    if (ratio == 0.0)        return ratio;     return ratio * sqrt((absx / ratio) * (absx / ratio)                        + (absy / ratio) * (absy / ratio)

Lines 1–500 of 1028. Load more