// Written in the D programming language. /** This module contains the $(LREF Complex) type, which is used to represent complex numbers, along with related mathematical operations and functions. $(LREF Complex) will eventually $(DDLINK deprecate, Deprecated Features, replace) the built-in types `cfloat`, `cdouble`, `creal`, `ifloat`, `idouble`, and `ireal`. Macros: TABLE_SV = <table border="1" cellpadding="4" cellspacing="0"> <caption>Special Values</caption> $0</table> PLUSMN = ± NAN = $(RED NAN) INFIN = ∞ PI = π Authors: Lars Tandle Kyllingstad, Don Clugston Copyright: Copyright (c) 2010, Lars T. Kyllingstad. License: $(HTTP boost.org/LICENSE_1_0.txt, Boost License 1.0) Source: $(PHOBOSSRC std/complex.d)*/module std.complex; import std.traits; /** Helper function that returns a complex number with the specified real and imaginary parts. Params: R = (template parameter) type of real part of complex number I = (template parameter) type of imaginary part of complex number re = real part of complex number to be constructed im = (optional) imaginary part of complex number, 0 if omitted. Returns: `Complex` instance with real and imaginary parts set to the values provided as input. If neither `re` nor `im` are floating-point numbers, the return type will be `Complex!double`. Otherwise, the return type is deduced using $(D std.traits.CommonType!(R, I)).*/auto complex(R)(const R re) @safe pure nothrow @nogcif (is(R : double)){ static if (isFloatingPoint!R) return Complex!R(re, 0); else return Complex!double(re, 0);} /// dittoauto complex(R, I)(const R re, const I im) @safe pure nothrow @nogcif (is(R : double) && is(I : double)){ static if (isFloatingPoint!R || isFloatingPoint!I) return Complex!(CommonType!(R, I))(re, im); else return Complex!double(re, im);} ///@safe pure nothrow unittest{ auto a = complex(1.0); static assert(is(typeof(a) == Complex!double)); assert(a.re == 1.0); assert(a.im == 0.0); auto b = complex(2.0L); static assert(is(typeof(b) == Complex!real)); assert(b.re == 2.0L); assert(b.im == 0.0L); auto c = complex(1.0, 2.0); static assert(is(typeof(c) == Complex!double)); assert(c.re == 1.0); assert(c.im == 2.0); auto d = complex(3.0, 4.0L); static assert(is(typeof(d) == Complex!real)); assert(d.re == 3.0); assert(d.im == 4.0L); auto e = complex(1); static assert(is(typeof(e) == Complex!double)); assert(e.re == 1); assert(e.im == 0); auto f = complex(1L, 2); static assert(is(typeof(f) == Complex!double)); assert(f.re == 1L); assert(f.im == 2); auto g = complex(3, 4.0L); static assert(is(typeof(g) == Complex!real)); assert(g.re == 3); assert(g.im == 4.0L);} /** A complex number parametrised by a type `T`, which must be either `float`, `double` or `real`.*/struct Complex(T)if (isFloatingPoint!T){ import std.format.spec : FormatSpec; import std.range.primitives : isOutputRange; /** The real part of the number. */ T re; /** The imaginary part of the number. */ T im; /** Converts the complex number to a string representation. The second form of this function is usually not called directly; instead, it is used via $(REF format, std,string), as shown in the examples below. Supported format characters are 'e', 'f', 'g', 'a', and 's'. See the $(MREF std, format) and $(REF format, std,string) documentation for more information. */ string toString() const @safe /* TODO: pure nothrow */ { import std.exception : assumeUnique; char[] buf; buf.reserve(100); auto fmt = FormatSpec!char("%s"); toString((const(char)[] s) { buf ~= s; }, fmt); static trustedAssumeUnique(T)(T t) @trusted { return assumeUnique(t); } return trustedAssumeUnique(buf); } static if (is(T == double)) /// @safe unittest { auto c = complex(1.2, 3.4); // Vanilla toString formatting: assert(c.toString() == "1.2+3.4i"); // Formatting with std.string.format specs: the precision and width // specifiers apply to both the real and imaginary parts of the // complex number. import std.format : format; assert(format("%.2f", c) == "1.20+3.40i"); assert(format("%4.1f", c) == " 1.2+ 3.4i"); } /// ditto void toString(Writer, Char)(scope Writer w, scope const ref FormatSpec!Char formatSpec) const if (isOutputRange!(Writer, const(Char)[])) { import std.format.write : formatValue; import std.math.traits : signbit; import std.range.primitives : put; formatValue(w, re, formatSpec); if (signbit(im) == 0) put(w, "+"); formatValue(w, im, formatSpec); put(w, "i"); } @safe pure nothrow @nogc: /** Construct a complex number with the specified real and imaginary parts. In the case where a single argument is passed that is not complex, the imaginary part of the result will be zero. */ this(R : T)(Complex!R z) { re = z.re; im = z.im; } /// ditto this(Rx : T, Ry : T)(const Rx x, const Ry y) { re = x; im = y; } /// ditto this(R : T)(const R r) { re = r; im = 0; } // ASSIGNMENT OPERATORS // this = complex ref Complex opAssign(R : T)(Complex!R z) { re = z.re; im = z.im; return this; } // this = numeric ref Complex opAssign(R : T)(const R r) { re = r; im = 0; return this; } // COMPARISON OPERATORS // this == complex bool opEquals(R : T)(Complex!R z) const { return re == z.re && im == z.im; } // this == numeric bool opEquals(R : T)(const R r) const { return re == r && im == 0; } // UNARY OPERATORS // +complex Complex opUnary(string op)() const if (op == "+") { return this; } // -complex Complex opUnary(string op)() const if (op == "-") { return Complex(-re, -im); } // BINARY OPERATORS // complex op complex Complex!(CommonType!(T,R)) opBinary(string op, R)(Complex!R z) const { alias C = typeof(return); auto w = C(this.re, this.im); return w.opOpAssign!(op)(z); } // complex op numeric Complex!(CommonType!(T,R)) opBinary(string op, R)(const R r) const if (isNumeric!R) { alias C = typeof(return); auto w = C(this.re, this.im); return w.opOpAssign!(op)(r); } // numeric + complex, numeric * complex Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(const R r) const if ((op == "+" || op == "*") && (isNumeric!R)) { return opBinary!(op)(r); } // numeric - complex Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(const R r) const if (op == "-" && isNumeric!R) { return Complex(r - re, -im); } // numeric / complex Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(const R r) const if (op == "/" && isNumeric!R) { version (FastMath) { // Compute norm(this) immutable norm = re * re + im * im; // Compute r * conj(this) immutable prod_re = r * re; immutable prod_im = r * -im; // Divide the product by the norm typeof(return) w = void; w.re = prod_re / norm; w.im = prod_im / norm; return w; } else { import core.math : fabs; typeof(return) w = void; if (fabs(re) < fabs(im)) { immutable ratio = re/im; immutable rdivd = r/(re*ratio + im); w.re = rdivd*ratio; w.im = -rdivd; } else { immutable ratio = im/re; immutable rdivd = r/(re + im*ratio); w.re = rdivd; w.im = -rdivd*ratio; } return w; } } // numeric ^^ complex Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(const R lhs) const if (op == "^^" && isNumeric!R) { import core.math : cos, sin; import std.math.exponential : exp, log; import std.math.constants : PI; Unqual!(CommonType!(T, R)) ab = void, ar = void; if (lhs >= 0) { // r = lhs // theta = 0 ab = lhs ^^ this.re; ar = log(lhs) * this.im; } else { // r = -lhs // theta = PI ab = (-lhs) ^^ this.re * exp(-PI * this.im); ar = PI * this.re + log(-lhs) * this.im; } return typeof(return)(ab * cos(ar), ab * sin(ar)); } // OP-ASSIGN OPERATORS // complex += complex, complex -= complex ref Complex opOpAssign(string op, C)(const C z) if ((op == "+" || op == "-") && is(C R == Complex!R)) { mixin ("re "~op~"= z.re;"); mixin ("im "~op~"= z.im;"); return this; } // complex *= complex ref Complex opOpAssign(string op, C)(const C z) if (op == "*" && is(C R == Complex!R)) { auto temp = re*z.re - im*z.im; im = im*z.re + re*z.im; re = temp; return this; } // complex /= complex ref Complex opOpAssign(string op, C)(const C z) if (op == "/" && is(C R == Complex!R)) { version (FastMath) { // Compute norm(z) immutable norm = z.re * z.re + z.im * z.im; // Compute this * conj(z) immutable prod_re = re * z.re - im * -z.im; immutable prod_im = im * z.re + re * -z.im; // Divide the product by the norm re = prod_re / norm; im = prod_im / norm; return this; } else { import core.math : fabs; if (fabs(z.re) < fabs(z.im)) { immutable ratio = z.re/z.im; immutable denom = z.re*ratio + z.im; immutable temp = (re*ratio + im)/denom; im = (im*ratio - re)/denom; re = temp; } else { immutable ratio = z.im/z.re; immutable denom = z.re + z.im*ratio; immutable temp = (re + im*ratio)/denom; im = (im - re*ratio)/denom; re = temp; } return this; } } // complex ^^= complex ref Complex opOpAssign(string op, C)(const C z) if (op == "^^" && is(C R == Complex!R)) { import core.math : cos, sin; import std.math.exponential : exp, log; immutable r = abs(this); immutable t = arg(this); immutable ab = r^^z.re * exp(-t*z.im); immutable ar = t*z.re + log(r)*z.im; re = ab*cos(ar); im = ab*sin(ar); return this; } // complex += numeric, complex -= numeric ref Complex opOpAssign(string op, U : T)(const U a) if (op == "+" || op == "-") { mixin ("re "~op~"= a;"); return this; } // complex *= numeric, complex /= numeric ref Complex opOpAssign(string op, U : T)(const U a) if (op == "*" || op == "/") { mixin ("re "~op~"= a;"); mixin ("im "~op~"= a;"); return this; } // complex ^^= real ref Complex opOpAssign(string op, R)(const R r) if (op == "^^" && isFloatingPoint!R) { import core.math : cos, sin; immutable ab = abs(this)^^r; immutable ar = arg(this)*r; re = ab*cos(ar); im = ab*sin(ar); return this; } // complex ^^= int ref Complex opOpAssign(string op, U)(const U i) if (op == "^^" && isIntegral!U) { switch (i) { case 0: re = 1.0; im = 0.0; break; case 1: // identity; do nothing break; case 2: this *= this; break; case 3: auto z = this; this *= z; this *= z; break; default: this ^^= cast(real) i; } return this; } /** Returns a complex number instance that correponds in size and in ABI to the associated C compiler's `_Complex` type. */ auto toNative() { import core.stdc.config : c_complex_float, c_complex_double, c_complex_real; static if (is(T == float)) return c_complex_float(re, im); else static if (is(T == double)) return c_complex_double(re, im); else return c_complex_real(re, im); }} @safe pure nothrow unittest{ import std.complex; static import core.math;Lines 1–500 of 1972. Load more
// Written in the D programming language. /** This module contains the $(LREF Complex) type, which is used to represent complex numbers, along with related mathematical operations and functions. $(LREF Complex) will eventually $(DDLINK deprecate, Deprecated Features, replace) the built-in types `cfloat`, `cdouble`, `creal`, `ifloat`, `idouble`, and `ireal`. Macros: TABLE_SV = <table border="1" cellpadding="4" cellspacing="0"> <caption>Special Values</caption> $0</table> PLUSMN = ± NAN = $(RED NAN) INFIN = ∞ PI = π Authors: Lars Tandle Kyllingstad, Don Clugston Copyright: Copyright (c) 2010, Lars T. Kyllingstad. License: $(HTTP boost.org/LICENSE_1_0.txt, Boost License 1.0) Source: $(PHOBOSSRC std/complex.d) */ module std.complex; import std.traits; /** Helper function that returns a complex number with the specified real and imaginary parts. Params: R = (template parameter) type of real part of complex number I = (template parameter) type of imaginary part of complex number re = real part of complex number to be constructed im = (optional) imaginary part of complex number, 0 if omitted. Returns: `Complex` instance with real and imaginary parts set to the values provided as input. If neither `re` nor `im` are floating-point numbers, the return type will be `Complex!double`. Otherwise, the return type is deduced using $(D std.traits.CommonType!(R, I)). */ auto complex(R)(const R re) @safe pure nothrow @nogc if (is(R : double)) { static if (isFloatingPoint!R) return Complex!R(re, 0); else return Complex!double(re, 0); } /// ditto auto complex(R, I)(const R re, const I im) @safe pure nothrow @nogc if (is(R : double) && is(I : double)) { static if (isFloatingPoint!R || isFloatingPoint!I) return Complex!(CommonType!(R, I))(re, im); else return Complex!double(re, im); } /// @safe pure nothrow unittest { auto a = complex(1.0); static assert(is(typeof(a) == Complex!double)); assert(a.re == 1.0); assert(a.im == 0.0); auto b = complex(2.0L); static assert(is(typeof(b) == Complex!real)); assert(b.re == 2.0L); assert(b.im == 0.0L); auto c = complex(1.0, 2.0); static assert(is(typeof(c) == Complex!double)); assert(c.re == 1.0); assert(c.im == 2.0); auto d = complex(3.0, 4.0L); static assert(is(typeof(d) == Complex!real)); assert(d.re == 3.0); assert(d.im == 4.0L); auto e = complex(1); static assert(is(typeof(e) == Complex!double)); assert(e.re == 1); assert(e.im == 0); auto f = complex(1L, 2); static assert(is(typeof(f) == Complex!double)); assert(f.re == 1L); assert(f.im == 2); auto g = complex(3, 4.0L); static assert(is(typeof(g) == Complex!real)); assert(g.re == 3); assert(g.im == 4.0L); } /** A complex number parametrised by a type `T`, which must be either `float`, `double` or `real`. */ struct Complex(T) if (isFloatingPoint!T) { import std.format.spec : FormatSpec; import std.range.primitives : isOutputRange; /** The real part of the number. */ T re; /** The imaginary part of the number. */ T im; /** Converts the complex number to a string representation. The second form of this function is usually not called directly; instead, it is used via $(REF format, std,string), as shown in the examples below. Supported format characters are 'e', 'f', 'g', 'a', and 's'. See the $(MREF std, format) and $(REF format, std,string) documentation for more information. */ string toString() const @safe /* TODO: pure nothrow */ { import std.exception : assumeUnique; char[] buf; buf.reserve(100); auto fmt = FormatSpec!char("%s"); toString((const(char)[] s) { buf ~= s; }, fmt); static trustedAssumeUnique(T)(T t) @trusted { return assumeUnique(t); } return trustedAssumeUnique(buf); } static if (is(T == double)) /// @safe unittest { auto c = complex(1.2, 3.4); // Vanilla toString formatting: assert(c.toString() == "1.2+3.4i"); // Formatting with std.string.format specs: the precision and width // specifiers apply to both the real and imaginary parts of the // complex number. import std.format : format; assert(format("%.2f", c) == "1.20+3.40i"); assert(format("%4.1f", c) == " 1.2+ 3.4i"); } /// ditto void toString(Writer, Char)(scope Writer w, scope const ref FormatSpec!Char formatSpec) const if (isOutputRange!(Writer, const(Char)[])) { import std.format.write : formatValue; import std.math.traits : signbit; import std.range.primitives : put; formatValue(w, re, formatSpec); if (signbit(im) == 0) put(w, "+"); formatValue(w, im, formatSpec); put(w, "i"); } @safe pure nothrow @nogc: /** Construct a complex number with the specified real and imaginary parts. In the case where a single argument is passed that is not complex, the imaginary part of the result will be zero. */ this(R : T)(Complex!R z) { re = z.re; im = z.im; } /// ditto this(Rx : T, Ry : T)(const Rx x, const Ry y) { re = x; im = y; } /// ditto this(R : T)(const R r) { re = r; im = 0; } // ASSIGNMENT OPERATORS // this = complex ref Complex opAssign(R : T)(Complex!R z) { re = z.re; im = z.im; return this; } // this = numeric ref Complex opAssign(R : T)(const R r) { re = r; im = 0; return this; } // COMPARISON OPERATORS // this == complex bool opEquals(R : T)(Complex!R z) const { return re == z.re && im == z.im; } // this == numeric bool opEquals(R : T)(const R r) const { return re == r && im == 0; } // UNARY OPERATORS // +complex Complex opUnary(string op)() const if (op == "+") { return this; } // -complex Complex opUnary(string op)() const if (op == "-") { return Complex(-re, -im); } // BINARY OPERATORS // complex op complex Complex!(CommonType!(T,R)) opBinary(string op, R)(Complex!R z) const { alias C = typeof(return); auto w = C(this.re, this.im); return w.opOpAssign!(op)(z); } // complex op numeric Complex!(CommonType!(T,R)) opBinary(string op, R)(const R r) const if (isNumeric!R) { alias C = typeof(return); auto w = C(this.re, this.im); return w.opOpAssign!(op)(r); } // numeric + complex, numeric * complex Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(const R r) const if ((op == "+" || op == "*") && (isNumeric!R)) { return opBinary!(op)(r); } // numeric - complex Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(const R r) const if (op == "-" && isNumeric!R) { return Complex(r - re, -im); } // numeric / complex Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(const R r) const if (op == "/" && isNumeric!R) { version (FastMath) { // Compute norm(this) immutable norm = re * re + im * im; // Compute r * conj(this) immutable prod_re = r * re; immutable prod_im = r * -im; // Divide the product by the norm typeof(return) w = void; w.re = prod_re / norm; w.im = prod_im / norm; return w; } else { import core.math : fabs; typeof(return) w = void; if (fabs(re) < fabs(im)) { immutable ratio = re/im; immutable rdivd = r/(re*ratio + im); w.re = rdivd*ratio; w.im = -rdivd; } else { immutable ratio = im/re; immutable rdivd = r/(re + im*ratio); w.re = rdivd; w.im = -rdivd*ratio; } return w; } } // numeric ^^ complex Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(const R lhs) const if (op == "^^" && isNumeric!R) { import core.math : cos, sin; import std.math.exponential : exp, log; import std.math.constants : PI; Unqual!(CommonType!(T, R)) ab = void, ar = void; if (lhs >= 0) { // r = lhs // theta = 0 ab = lhs ^^ this.re; ar = log(lhs) * this.im; } else { // r = -lhs // theta = PI ab = (-lhs) ^^ this.re * exp(-PI * this.im); ar = PI * this.re + log(-lhs) * this.im; } return typeof(return)(ab * cos(ar), ab * sin(ar)); } // OP-ASSIGN OPERATORS // complex += complex, complex -= complex ref Complex opOpAssign(string op, C)(const C z) if ((op == "+" || op == "-") && is(C R == Complex!R)) { mixin ("re "~op~"= z.re;"); mixin ("im "~op~"= z.im;"); return this; } // complex *= complex ref Complex opOpAssign(string op, C)(const C z) if (op == "*" && is(C R == Complex!R)) { auto temp = re*z.re - im*z.im; im = im*z.re + re*z.im; re = temp; return this; } // complex /= complex ref Complex opOpAssign(string op, C)(const C z) if (op == "/" && is(C R == Complex!R)) { version (FastMath) { // Compute norm(z) immutable norm = z.re * z.re + z.im * z.im; // Compute this * conj(z) immutable prod_re = re * z.re - im * -z.im; immutable prod_im = im * z.re + re * -z.im; // Divide the product by the norm re = prod_re / norm; im = prod_im / norm; return this; } else { import core.math : fabs; if (fabs(z.re) < fabs(z.im)) { immutable ratio = z.re/z.im; immutable denom = z.re*ratio + z.im; immutable temp = (re*ratio + im)/denom; im = (im*ratio - re)/denom; re = temp; } else { immutable ratio = z.im/z.re; immutable denom = z.re + z.im*ratio; immutable temp = (re + im*ratio)/denom; im = (im - re*ratio)/denom; re = temp; } return this; } } // complex ^^= complex ref Complex opOpAssign(string op, C)(const C z) if (op == "^^" && is(C R == Complex!R)) { import core.math : cos, sin; import std.math.exponential : exp, log; immutable r = abs(this); immutable t = arg(this); immutable ab = r^^z.re * exp(-t*z.im); immutable ar = t*z.re + log(r)*z.im; re = ab*cos(ar); im = ab*sin(ar); return this; } // complex += numeric, complex -= numeric ref Complex opOpAssign(string op, U : T)(const U a) if (op == "+" || op == "-") { mixin ("re "~op~"= a;"); return this; } // complex *= numeric, complex /= numeric ref Complex opOpAssign(string op, U : T)(const U a) if (op == "*" || op == "/") { mixin ("re "~op~"= a;"); mixin ("im "~op~"= a;"); return this; } // complex ^^= real ref Complex opOpAssign(string op, R)(const R r) if (op == "^^" && isFloatingPoint!R) { import core.math : cos, sin; immutable ab = abs(this)^^r; immutable ar = arg(this)*r; re = ab*cos(ar); im = ab*sin(ar); return this; } // complex ^^= int ref Complex opOpAssign(string op, U)(const U i) if (op == "^^" && isIntegral!U) { switch (i) { case 0: re = 1.0; im = 0.0; break; case 1: // identity; do nothing break; case 2: this *= this; break; case 3: auto z = this; this *= z; this *= z; break; default: this ^^= cast(real) i; } return this; } /** Returns a complex number instance that correponds in size and in ABI to the associated C compiler's `_Complex` type. */ auto toNative() { import core.stdc.config : c_complex_float, c_complex_double, c_complex_real; static if (is(T == float)) return c_complex_float(re, im); else static if (is(T == double)) return c_complex_double(re, im); else return c_complex_real(re, im); } } @safe pure nothrow unittest { import std.complex; static import core.math;