Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with HTTPS or Subversion.

Download ZIP

Loading…

std.complex update: added expi() and moved some functions to module level #559

Merged
merged 5 commits into from

3 participants

@kyllingstad
Collaborator

I have just copied the implementation of std.math.expi() more or less verbatim. I don't know any assembler, so someone should check that the implementation is still valid when the function returns a struct as opposed to a built-in creal. I can't see why it shouldn't be, but you never know... Since there are several functions in std.complex now that are basically ripped from std.math, I have added Don as module author.

Furthermore, I have moved abs(), arg() and conj() down to module level. (They were previously properties of Complex.) This is mainly for consistency with std.math. Backwards compatibility is preserved thanks to the magic of UFCS. (Not that there is much compatibility to worry about; it doesn't seem like too many people are using this module yet.)

kyllingstad added some commits
@kyllingstad kyllingstad std.complex: Moved abs, arg & conj to module level
Backwards compatibility is preserved through the magic of UFCS.
5b4351d
@kyllingstad kyllingstad std.complex: Cleanup
Moved some stuff around, removed some extraneous whitespace, made
function documentation more consistent.
34bc329
@kyllingstad kyllingstad Added std.complex.expi()
Also added Don as author, since there are several functions in
std.complex, including expi(), that are copied more or less verbatim
from std.math.
03295d1
@kyllingstad kyllingstad std.complex: Minor doc improvement c80bb85
@donc donc was assigned
@kyllingstad
Collaborator

Don, I just noticed there is an "assign" feature on GitHub pull requests now. You are probably the one who should review this request, that is why I assigned it to you. (Please let me know if you'd rather I'd not do that in the future.)

@donc
Collaborator

I don't know any assembler, so someone should check that the implementation is still valid when the function returns a struct as opposed to a built-in creal. I can't see why it shouldn't be, but you never know...

Unfortunately, it isn't valid. The ABI says that creals are returned on the top two elements of the x87 stack. There is no other type which does that. In fact this instruction is the only example of built-in complex type on any architecture that I know of...
They will have to be returned in a struct. Not worth using the asm at all, just use a separate sin() and cos() function. We can maybe try to get the backend to recognize and create an fsincos at some point, but I doubt this function is really performance critical anyway.

@kyllingstad
Collaborator

Too bad. That means the only reason to have this function in std.complex is to preserve backwards compatibility when the built-in types are deprecated. I'll fix it later today.

@kyllingstad
Collaborator

Done.

@donc donc merged commit 944987d into from
@klickverbot
Collaborator

It probably really doesn't matter much at this point, but I just wanted to note that this is a breaking change, even with UFCS in place – D doesn't have ADL.

@kyllingstad
Collaborator

Please explain? (Or perhaps give an example of where it will fail.)

@klickverbot
Collaborator
import std.complex : complex;

auto c = complex(2.0);
c.abs();

Alternatively, consider that a module doesn't import std.complex at all, for example because Complex is merely passed as a template parameter (this is arguably not so likely with Complex, but it could well happen in similar other cases).

@kyllingstad
Collaborator

Ok, I see. Good point. I don't think there will be too many cases of this happening, but maybe I should reinstate the struct methods (only deprecated and written in terms of the module-level functions), just in case.

@ghost Unknown referenced this pull request from a commit
Commit has since been removed from the repository and is no longer available.
@ghost Unknown referenced this pull request from a commit
Commit has since been removed from the repository and is no longer available.
@ghost Unknown referenced this pull request from a commit
Commit has since been removed from the repository and is no longer available.
@ghost Unknown referenced this pull request from a commit
Commit has since been removed from the repository and is no longer available.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Commits on May 1, 2012
  1. @kyllingstad

    std.complex: Moved abs, arg & conj to module level

    kyllingstad authored
    Backwards compatibility is preserved through the magic of UFCS.
  2. @kyllingstad

    std.complex: Cleanup

    kyllingstad authored
    Moved some stuff around, removed some extraneous whitespace, made
    function documentation more consistent.
  3. @kyllingstad

    Added std.complex.expi()

    kyllingstad authored
    Also added Don as author, since there are several functions in
    std.complex, including expi(), that are copied more or less verbatim
    from std.math.
  4. @kyllingstad
Commits on May 11, 2012
  1. @kyllingstad
This page is out of date. Refresh to see the latest.
Showing with 119 additions and 135 deletions.
  1. +119 −135 std/complex.d
View
254 std/complex.d
@@ -1,9 +1,13 @@
// Written in the D programming language.
-/** Module that will replace the built-in types $(D cfloat), $(D cdouble),
- $(D creal), $(D ifloat), $(D idouble), and $(D ireal).
+/** This module contains the $(LREF Complex) type, which is used to represent
+ _complex numbers, along with related mathematical operations and functions.
- Authors: Lars Tandle Kyllingstad
+ $(LREF Complex) will eventually $(LINK2 ../deprecate.html, replace)
+ the built-in types $(D cfloat), $(D cdouble), $(D creal), $(D ifloat),
+ $(D idouble), and $(D ireal).
+
+ Authors: Lars Tandle Kyllingstad, Don Clugston
Copyright: Copyright (c) 2010, Lars T. Kyllingstad.
License: $(WEB boost.org/LICENSE_1_0.txt, Boost License 1.0)
Source: $(PHOBOSSRC std/_complex.d)
@@ -17,8 +21,6 @@ import std.numeric;
import std.traits;
-
-
/** Helper function that returns a _complex number with the specified
real and imaginary parts.
@@ -61,7 +63,6 @@ auto complex(R, I)(R re, I im) @safe pure nothrow
return Complex!double(re, im);
}
-
unittest
{
auto a = complex(1.0);
@@ -101,8 +102,6 @@ unittest
}
-
-
/** A complex number parametrised by a type $(D T), which must be either
$(D float), $(D double) or $(D real).
*/
@@ -114,38 +113,45 @@ struct Complex(T) if (isFloatingPoint!T)
/** The imaginary part of the number. */
T im;
+ /** Converts the complex number to a string representation.
-@safe pure nothrow // The following functions depend only on std.math.
-{
-
- /** Calculate the absolute value (or modulus) of the number. */
- @property T abs() const
- {
- return hypot(re, im);
- }
+ If a $(D sink) delegate is specified, the string is passed to it
+ and this function returns $(D null). Otherwise, this function
+ returns the string representation directly.
+ The output format is controlled via $(D formatSpec), which should consist
+ of a single POSIX format specifier, including the percent (%) character.
+ Note that complex numbers are floating point numbers, so the only
+ valid format characters are 'e', 'f', 'g', 'a', and 's', where 's'
+ gives the default behaviour. Positional parameters are not valid
+ in this context.
- /** Calculate the argument (or phase) of the number. */
- @property T arg() const
+ See the $(LINK2 std_format.html, std.format documentation) for
+ more information.
+ */
+ string toString(scope void delegate(const(char)[]) sink = null,
+ string formatSpec = "%s")
+ const
{
- return atan2(im, re);
- }
-
+ if (sink == null)
+ {
+ char[] buf;
+ buf.reserve(100);
+ toString((const(char)[] s) { buf ~= s; }, formatSpec);
+ return cast(string) buf;
+ }
- /** Return the complex conjugate of the number. */
- @property Complex conj() const
- {
- return Complex(re, -im);
+ formattedWrite(sink, formatSpec, re);
+ if (signbit(im) == 0) sink("+");
+ formattedWrite(sink, formatSpec, im);
+ sink("i");
+ return null;
}
-
-
+@safe pure nothrow:
// ASSIGNMENT OPERATORS
- // TODO: Make operators return by ref when DMD bug 2460 is fixed.
-
-
// this = complex
ref Complex opAssign(R : T)(Complex!R z)
{
@@ -154,7 +160,6 @@ struct Complex(T) if (isFloatingPoint!T)
return this;
}
-
// this = numeric
ref Complex opAssign(R : T)(R r)
{
@@ -163,30 +168,22 @@ struct Complex(T) if (isFloatingPoint!T)
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)(R r) const
{
return re == r && im == 0;
}
-
-
// UNARY OPERATORS
-
// +complex
Complex opUnary(string op)() const
if (op == "+")
@@ -194,7 +191,6 @@ struct Complex(T) if (isFloatingPoint!T)
return this;
}
-
// -complex
Complex opUnary(string op)() const
if (op == "-")
@@ -202,11 +198,8 @@ struct Complex(T) if (isFloatingPoint!T)
return Complex(-re, -im);
}
-
-
// BINARY OPERATORS
-
// complex op complex
Complex!(CommonType!(T,R)) opBinary(string op, R)(Complex!R z) const
{
@@ -215,7 +208,6 @@ struct Complex(T) if (isFloatingPoint!T)
return w.opOpAssign!(op)(z);
}
-
// complex op numeric
Complex!(CommonType!(T,R)) opBinary(string op, R)(R r) const
if (isNumeric!R)
@@ -225,7 +217,6 @@ struct Complex(T) if (isFloatingPoint!T)
return w.opOpAssign!(op)(r);
}
-
// numeric + complex, numeric * complex
Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(R r) const
if ((op == "+" || op == "*") && (isNumeric!R))
@@ -233,7 +224,6 @@ struct Complex(T) if (isFloatingPoint!T)
return opBinary!(op)(r);
}
-
// numeric - complex
Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(R r) const
if (op == "-" && isNumeric!R)
@@ -241,7 +231,6 @@ struct Complex(T) if (isFloatingPoint!T)
return Complex(r - re, -im);
}
-
// numeric / complex
Complex!(CommonType!(T, R)) opBinaryRight(string op, R)(R r) const
if (op == "/" && isNumeric!R)
@@ -269,10 +258,7 @@ struct Complex(T) if (isFloatingPoint!T)
return w;
}
-
-
- // OPASSIGN OPERATORS
-
+ // OP-ASSIGN OPERATORS
// complex += complex, complex -= complex
ref Complex opOpAssign(string op, C)(C z)
@@ -283,7 +269,6 @@ struct Complex(T) if (isFloatingPoint!T)
return this;
}
-
// complex *= complex
ref Complex opOpAssign(string op, C)(C z)
if (op == "*" && is(C R == Complex!R))
@@ -294,7 +279,6 @@ struct Complex(T) if (isFloatingPoint!T)
return this;
}
-
// complex /= complex
ref Complex opOpAssign(string op, C)(C z)
if (op == "/" && is(C R == Complex!R))
@@ -320,13 +304,12 @@ struct Complex(T) if (isFloatingPoint!T)
return this;
}
-
// complex ^^= complex
ref Complex opOpAssign(string op, C)(C z)
if (op == "^^" && is(C R == Complex!R))
{
- FPTemporary!T r = abs;
- FPTemporary!T t = arg;
+ FPTemporary!T r = abs(this);
+ FPTemporary!T t = arg(this);
FPTemporary!T ab = r^^z.re * exp(-t*z.im);
FPTemporary!T ar = t*z.re + log(r)*z.im;
@@ -335,7 +318,6 @@ struct Complex(T) if (isFloatingPoint!T)
return this;
}
-
// complex += numeric, complex -= numeric
ref Complex opOpAssign(string op, U : T)(U a)
if (op == "+" || op == "-")
@@ -344,7 +326,6 @@ struct Complex(T) if (isFloatingPoint!T)
return this;
}
-
// complex *= numeric, complex /= numeric
ref Complex opOpAssign(string op, U : T)(U a)
if (op == "*" || op == "/")
@@ -354,19 +335,17 @@ struct Complex(T) if (isFloatingPoint!T)
return this;
}
-
// complex ^^= real
ref Complex opOpAssign(string op, R)(R r)
if (op == "^^" && isFloatingPoint!R)
{
- FPTemporary!T ab = abs^^r;
- FPTemporary!T ar = arg*r;
+ FPTemporary!T ab = abs(this)^^r;
+ FPTemporary!T ar = arg(this)*r;
re = ab*std.math.cos(ar);
im = ab*std.math.sin(ar);
return this;
}
-
// complex ^^= int
ref Complex opOpAssign(string op, U)(U i)
if (op == "^^" && isIntegral!U)
@@ -393,60 +372,12 @@ struct Complex(T) if (isFloatingPoint!T)
}
return this;
}
-
-} // @safe pure nothrow
-
-
-
- /** Convert the complex number to a string representation.
-
- If a $(D sink) delegate is specified, the string is passed to it
- and this function returns $(D null). Otherwise, this function
- returns the string representation directly.
-
- The output format is controlled via $(D formatSpec), which should consist
- of a single POSIX format specifier, including the percent (%) character.
- Note that complex numbers are floating point numbers, so the only
- valid format characters are 'e', 'f', 'g', 'a', and 's', where 's'
- gives the default behaviour. Positional parameters are not valid
- in this context.
-
- See the $(LINK2 std_format.html, std.format documentation) for
- more information.
- */
- string toString
- (scope void delegate(const(char)[]) sink = null, string formatSpec = "%s")
- const
- {
- if (sink == null)
- {
- char[] buf;
- buf.reserve(100);
- toString((const(char)[] s) { buf ~= s; }, formatSpec);
- return cast(string) buf;
- }
-
- formattedWrite(sink, formatSpec, re);
- if (signbit(im) == 0) sink("+");
- formattedWrite(sink, formatSpec, im);
- sink("i");
- return null;
- }
}
-
unittest
{
enum EPS = double.epsilon;
-
- // Check abs() and arg()
- auto c1 = Complex!double(1.0, 1.0);
- assert (approxEqual(c1.abs, std.math.sqrt(2.0), EPS));
- assert (approxEqual(c1.arg, PI_4, EPS));
-
- auto c1c = c1.conj;
- assert (c1c.re == 1.0 && c1c.im == -1.0);
-
+ auto c1 = complex(1.0, 1.0);
// Check unary operations.
auto c2 = Complex!double(0.5, 2.0);
@@ -457,7 +388,6 @@ unittest
assert ((-c2).im == -(c2.im));
assert (c2 == -(-c2));
-
// Check complex-complex operations.
auto cpc = c1 + c2;
assert (cpc.re == c1.re + c2.re);
@@ -468,18 +398,17 @@ unittest
assert (cmc.im == c1.im - c2.im);
auto ctc = c1 * c2;
- assert (approxEqual(ctc.abs, c1.abs*c2.abs, EPS));
- assert (approxEqual(ctc.arg, c1.arg+c2.arg, EPS));
+ assert (approxEqual(abs(ctc), abs(c1)*abs(c2), EPS));
+ assert (approxEqual(arg(ctc), arg(c1)+arg(c2), EPS));
auto cdc = c1 / c2;
- assert (approxEqual(cdc.abs, c1.abs/c2.abs, EPS));
- assert (approxEqual(cdc.arg, c1.arg-c2.arg, EPS));
+ assert (approxEqual(abs(cdc), abs(c1)/abs(c2), EPS));
+ assert (approxEqual(arg(cdc), arg(c1)-arg(c2), EPS));
auto cec = c1^^c2;
assert (approxEqual(cec.re, 0.11524131979943839881, EPS));
assert (approxEqual(cec.im, 0.21870790452746026696, EPS));
-
// Check complex-real operations.
double a = 123.456;
@@ -496,8 +425,8 @@ unittest
assert (ctr.im == c1.im*a);
auto cdr = c1 / a;
- assert (approxEqual(cdr.abs, c1.abs/a, EPS));
- assert (approxEqual(cdr.arg, c1.arg, EPS));
+ assert (approxEqual(abs(cdr), abs(c1)/a, EPS));
+ assert (approxEqual(arg(cdr), arg(c1), EPS));
auto rpc = a + c1;
assert (rpc == cpr);
@@ -510,25 +439,23 @@ unittest
assert (rtc == ctr);
auto rdc = a / c1;
- assert (approxEqual(rdc.abs, a/c1.abs, EPS));
- assert (approxEqual(rdc.arg, -c1.arg, EPS));
+ assert (approxEqual(abs(rdc), a/abs(c1), EPS));
+ assert (approxEqual(arg(rdc), -arg(c1), EPS));
auto cer = c1^^3.0;
- assert (approxEqual(cer.abs, c1.abs^^3, EPS));
- assert (approxEqual(cer.arg, c1.arg*3, EPS));
-
+ assert (approxEqual(abs(cer), abs(c1)^^3, EPS));
+ assert (approxEqual(arg(cer), arg(c1)*3, EPS));
// Check Complex-int operations.
foreach (i; 0..6)
{
auto cei = c1^^i;
- assert (approxEqual(cei.abs, c1.abs^^i, EPS));
+ assert (approxEqual(abs(cei), abs(c1)^^i, EPS));
// Use cos() here to deal with arguments that go outside
// the (-pi,pi] interval (only an issue for i>3).
- assert (approxEqual(std.math.cos(cei.arg), std.math.cos(c1.arg*i), EPS));
+ assert (approxEqual(std.math.cos(arg(cei)), std.math.cos(arg(c1)*i), EPS));
}
-
// Check operations between different complex types.
auto cf = Complex!float(1.0, 1.0);
auto cr = Complex!real(1.0, 1.0);
@@ -557,7 +484,6 @@ unittest
assert (z == 1.0L);
assert (z.re == 1.0 && z.im == 0.0);
-
auto w = Complex!real(1.0, 1.0);
z = w;
assert (z == w);
@@ -581,7 +507,7 @@ unittest
assert (s1 == z1.toString());
// Using custom format specifier
- auto z2 = z1.conj;
+ auto z2 = conj(z1);
char[] s2;
z2.toString((const(char)[] c) { s2 ~= c; }, "%.8e");
assert (s2 == "1.23456789e-01-1.23456789e-01i");
@@ -589,9 +515,7 @@ unittest
}
-
-
-/* Fold Complex!(Complex!T) to Complex!T.
+/* Makes Complex!(Complex!T) fold to Complex!T.
The rationale for this is that just like the real line is a
subspace of the complex plane, the complex plane is a subspace
@@ -627,9 +551,48 @@ unittest
}
+/** Calculates the absolute value (or modulus) of a complex number. */
+T abs(T)(Complex!T z) @safe pure nothrow
+{
+ return hypot(z.re, z.im);
+}
+
+unittest
+{
+ assert (abs(complex(1.0)) == 1.0);
+ assert (abs(complex(0.0, 1.0)) == 1.0);
+ assert (abs(complex(1.0L, -2.0L)) == std.math.sqrt(5.0L));
+}
-/** Construct a complex number given its absolute value and argument. */
+/** Calculates the argument (or phase) of a complex number. */
+T arg(T)(Complex!T z) @safe pure nothrow
+{
+ return atan2(z.im, z.re);
+}
+
+unittest
+{
+ assert (arg(complex(1.0)) == 0.0);
+ assert (arg(complex(0.0L, 1.0L)) == PI_2);
+ assert (arg(complex(1.0L, 1.0L)) == PI_4);
+}
+
+
+/** Returns the complex conjugate of a complex number. */
+Complex!T conj(T)(Complex!T z) @safe pure nothrow
+{
+ return Complex!T(z.re, -z.im);
+}
+
+unittest
+{
+ assert (conj(complex(1.0)) == complex(1.0));
+ assert (conj(complex(1.0, 2.0)) == complex(1.0, -2.0));
+}
+
+
+/** Constructs a complex number given its absolute value and argument. */
Complex!(CommonType!(T, U)) fromPolar(T, U)(T modulus, U argument)
@safe pure nothrow
{
@@ -645,8 +608,6 @@ unittest
}
-
-
/** Trigonometric functions. */
Complex!T sin(T)(Complex!T z) @safe pure nothrow
{
@@ -677,6 +638,29 @@ unittest{
}
+/** Calculates cos(y) + i sin(y).
+
+ Note:
+ $(D expi) is included here for convenience and for easy migration of code
+ that uses $(XREF math,_expi). Unlike $(XREF math,_expi), which uses the
+ x87 $(I fsincos) instruction when possible, this function is no faster
+ than calculating cos(y) and sin(y) separately.
+*/
+Complex!real expi(real y) @trusted pure nothrow
+{
+ return Complex!real(std.math.cos(y), std.math.sin(y));
+}
+
+unittest
+{
+ assert(expi(1.3e5L) == complex(std.math.cos(1.3e5L), std.math.sin(1.3e5L)));
+ assert(expi(0.0L) == 1.0L);
+ auto z1 = expi(1.234);
+ auto z2 = std.math.expi(1.234);
+ assert(z1.re == z2.re && z1.im == z2.im);
+}
+
+
/** Square root. */
Complex!T sqrt(T)(Complex!T z) @safe pure nothrow
{
Something went wrong with that request. Please try again.