-
Notifications
You must be signed in to change notification settings - Fork 57
/
Copy pathcomplex_numbers.complex.m
189 lines (148 loc) · 5.21 KB
/
complex_numbers.complex.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
%-----------------------------------------------------------------------------%
% vim: ft=mercury ts=4 sw=4 et
%-----------------------------------------------------------------------------%
% Copyright (C) 1997-1998, 2001, 2005-2006 The University of Melbourne.
% Copyright (C) 2015, 2018 The Mercury team.
% This file is distributed under the terms specified in COPYING.LIB.
%-----------------------------------------------------------------------------%
%
% File: complex.m.
% Main author: fjh.
% Stability: medium.
%
% Complex numbers.
%
% Note that the overloaded versions of the binary operators that
% provide mixed-type arithmetic are defined in other modules.
%
% See also:
% complex_float.m, float_complex.m
% imag.m, complex_imag.m, imag_complex.m
%
%-----------------------------------------------------------------------------%
:- module complex_numbers.complex.
:- interface.
%-----------------------------------------------------------------------------%
% The constructor cmplx/2 is made public, but generally it is most convenient
% to use the syntax `X + Y*i' for complex numbers, where `i' is declared in
% module `imag'. Due to the wonders of logic programming, this works fine for
% both constructing and pattern matching; with intermodule optimization
% enabled, the compiler should generate equally good code for it.
:- type complex
---> cmplx(
float, % real part
float % imag part
).
%-----------------------------------------------------------------------------%
% Convert float to complex.
%
:- func complex(float) = complex.
% Extract real part.
%
:- func real(complex) = float.
% Extract imaginary part.
%
:- func imag(complex) = float.
% Square of absolute value.
%
:- func abs2(complex) = float.
% Absolute value (a.k.a. modulus).
%
:- func abs(complex) = float.
% Argument (a.k.a. phase, or amplitude, or angle).
% This function returns the principle value:
%
% for all Z, -pi < arg(Z) and arg(Z) =< pi.
%
:- func arg(complex) = float.
% Complex conjugate.
%
:- func conj(complex) = complex.
% Addition.
%
:- func complex + complex = complex.
:- mode in + in = uo is det.
% Subtraction.
%
:- func complex - complex = complex.
:- mode in - in = uo is det.
% Multiplication.
%
:- func complex * complex = complex.
:- mode in * in = uo is det.
% Division.
%
:- func complex / complex = complex.
:- mode in / in = uo is det.
% Unary plus.
%
:- func + complex = complex.
:- mode + in = uo is det.
% Unary minus.
%
:- func - complex = complex.
:- mode - in = uo is det.
% sqr(X) = X * X.
%
:- func sqr(complex) = complex.
% Square root.
%
:- func sqrt(complex) = complex.
% cis(Theta) = cos(Theta) + i * sin(Theta)
%
:- func cis(float) = complex.
% polar_to_complex(R, Theta).
% Conversion from polar coordinates.
%
:- func polar_to_complex(float, float) = complex.
% polar_to_complex(Z, R, Theta).
% Conversion to polar coordinates.
%
:- pred complex_to_polar(complex::in, float::out, float::out) is det.
%---------------------------------------------------------------------------%
%---------------------------------------------------------------------------%
:- implementation.
:- import_module float.
:- import_module math.
%---------------------------------------------------------------------------%
complex(Real) = cmplx(Real, 0.0).
real(cmplx(Real, _Imag)) = Real.
imag(cmplx(_Real, Imag)) = Imag.
cmplx(Xr, Xi) + cmplx(Yr, Yi) = cmplx(Xr + Yr, Xi + Yi).
cmplx(Xr, Xi) - cmplx(Yr, Yi) = cmplx(Xr - Yr, Xi - Yi).
cmplx(Xr, Xi) * cmplx(Yr, Yi) =
cmplx(Xr * Yr - Xi * Yi, Xr * Yi + Xi * Yr).
cmplx(Xr, Xi) / cmplx(Yr, Yi) =
cmplx((Xr * Yr + Xi * Yi) / Div, (Xi * Yr - Xr * Yi) / Div) :-
Div = (Yr * Yr + Yi * Yi).
% Here's the derivation of the formula for complex division:
% cmplx(Xr, Xi) / cmplx(Yr, Yi) =
% (cmplx(Xr, Xi) / cmplx(Yr, Yi)) * 1.0 =
% (cmplx(Xr, Xi) / cmplx(Yr, Yi)) * (cmplx(Yr, -Yi) / cmplx(Yr, -Yi)) =
% (cmplx(Xr, Xi) * (cmplx(Yr, -Yi)) / (cmplx(Yr, Yi) * cmplx(Yr, -Yi)) =
% (cmplx(Xr, Xi) * (cmplx(Yr, -Yi)) / (Yr * Yr + Yi * Yi) =
% cmplx(Xr * Yr + Xi * Yi, Xi * Yr - Xr * Yi) / (Yr * Yr + Yi * Yi) =
% cmplx((Xr * Yr + Xi * Yi) / Div, (Xi * Yr - Xr * Yi) / Div) :-
% Div = (Yr * Yr + Yi * Yi).
+ cmplx(R, I) = cmplx(+ R, + I).
- cmplx(R, I) = cmplx(- R, - I).
abs2(cmplx(R, I)) = R*R + I*I.
abs(Z) = sqrt(abs2(Z)).
arg(cmplx(R, I)) = atan2(I, R).
conj(cmplx(R, I)) = cmplx(R, -I).
sqr(cmplx(Re0, Im0)) = cmplx(Re, Im) :-
Re = Re0 * Re0 - Im0 * Im0,
Im = 2.0 * Re0 * Im0.
sqrt(Z0) = Z :-
complex_to_polar(Z0, Magnitude0, Theta0),
Magnitude = sqrt(Magnitude0),
Theta = Theta0 / 2.0,
Z = polar_to_complex(Magnitude, Theta).
complex_to_polar(Z, abs(Z), arg(Z)).
polar_to_complex(Magnitude, Theta) = cmplx(Real, Imag) :-
Real = Magnitude * cos(Theta),
Imag = Magnitude * sin(Theta).
cis(Theta) = cmplx(cos(Theta), sin(Theta)).
%-----------------------------------------------------------------------------%
:- end_module complex_numbers.complex.
%-----------------------------------------------------------------------------%