1. my class X::Numeric::Real { ... };
  2. my class Complex is Cool does Numeric {
  3. has num $.re;
  4. has num $.im;
  5. method !SET-SELF(Num() \re, Num() \im) {
  6. $!re = re;
  7. $!im = im;
  8. self
  9. }
  10. proto method new(|) {*}
  11. multi method new() { self.new: 0, 0 }
  12. multi method new(Real \re, Real \im) { nqp::create(self)!SET-SELF(re, im) }
  13. multi method WHICH(Complex:D:) {
  14. nqp::box_s(
  15. nqp::concat(
  16. nqp::if(
  17. nqp::eqaddr(self.WHAT,Complex),
  18. 'Complex|',
  19. nqp::concat(nqp::unbox_s(self.^name), '|')
  20. ),
  21. nqp::concat($!re, nqp::concat('|', $!im))
  22. ),
  23. ValueObjAt
  24. )
  25. }
  26. method reals(Complex:D:) {
  27. (self.re, self.im);
  28. }
  29. method isNaN(Complex:D:) {
  30. self.re.isNaN || self.im.isNaN;
  31. }
  32. method coerce-to-real(Complex:D: $exception-target) {
  33. $!im ≅ 0e0
  34. ?? $!re
  35. !! Failure.new(X::Numeric::Real.new(target => $exception-target, reason => "imaginary part not zero", source => self))
  36. }
  37. multi method Real(Complex:D:) { self.coerce-to-real(Real); }
  38. # should probably be eventually supplied by role Numeric
  39. method Num(Complex:D:) { self.coerce-to-real(Num).Num; }
  40. method Int(Complex:D:) { self.coerce-to-real(Int).Int; }
  41. method Rat(Complex:D: $epsilon?) {
  42. self.coerce-to-real(Rat).Rat( |($epsilon // Empty) );
  43. }
  44. method FatRat(Complex:D: $epsilon?) {
  45. self.coerce-to-real(FatRat).FatRat( |($epsilon // Empty) );
  46. }
  47. multi method Bool(Complex:D:) {
  48. $!re != 0e0 || $!im != 0e0;
  49. }
  50. method Complex() { self }
  51. multi method Str(Complex:D:) {
  52. nqp::concat(
  53. $!re,
  54. nqp::concat(
  55. nqp::if(nqp::iseq_i(nqp::ord($!im),45),'','+'),
  56. nqp::concat(
  57. $!im,
  58. nqp::if(nqp::isnanorinf($!im),'\\i','i')
  59. )
  60. )
  61. )
  62. }
  63. multi method perl(Complex:D:) {
  64. '<' ~ self.Str ~ '>';
  65. }
  66. method conj(Complex:D:) {
  67. Complex.new($.re, -$.im);
  68. }
  69. method abs(Complex $x:) {
  70. nqp::p6box_n(nqp::sqrt_n(
  71. nqp::add_n(
  72. nqp::mul_n($!re, $!re),
  73. nqp::mul_n($!im, $!im),
  74. )
  75. ))
  76. }
  77. method polar() {
  78. $.abs, $!im.atan2($!re);
  79. }
  80. multi method log(Complex:D:) {
  81. my Num ($mag, $angle) = self.polar;
  82. Complex.new($mag.log, $angle);
  83. }
  84. method cis(Complex:D:) {
  85. self.cos + self.sin*Complex.new(0,1)
  86. }
  87. method sqrt(Complex:D:) {
  88. my Num $abs = self.abs;
  89. my Num $re = (($abs + self.re)/2).sqrt;
  90. my Num $im = (($abs - self.re)/2).sqrt;
  91. Complex.new($re, self.im < 0 ?? -$im !! $im);
  92. }
  93. multi method exp(Complex:D:) {
  94. my Num $mag = $!re.exp;
  95. Complex.new($mag * $!im.cos, $mag * $!im.sin);
  96. }
  97. method roots(Complex:D: Int() $n) {
  98. return NaN if $n < 1;
  99. return self if $n == 1;
  100. for $!re, $!im {
  101. return NaN if $_ eq 'Inf' || $_ eq '-Inf' || $_ eq 'NaN';
  102. }
  103. my ($mag, $angle) = self.polar;
  104. $mag **= 1e0 / $n;
  105. (^$n).map: { $mag.unpolar( ($angle + $_ * 2e0 * pi) / $n) };
  106. }
  107. method sin(Complex:D:) {
  108. $!re.sin * $!im.cosh + ($!re.cos * $!im.sinh)i;
  109. }
  110. method asin(Complex:D:) {
  111. (Complex.new(0e0, -1e0) * log((self)i + sqrt(1e0 - self * self)));
  112. }
  113. method cos(Complex:D:) {
  114. $!re.cos * $!im.cosh - ($!re.sin * $!im.sinh)i;
  115. }
  116. method acos(Complex:D:) {
  117. (pi / 2e0) - self.asin;
  118. }
  119. method tan(Complex:D:) {
  120. self.sin / self.cos;
  121. }
  122. method atan(Complex:D:) {
  123. ((log(1e0 - (self)i) - log(1e0 + (self)i))i / 2e0);
  124. }
  125. method sec(Complex:D:) {
  126. 1e0 / self.cos;
  127. }
  128. method asec(Complex:D:) {
  129. (1e0 / self).acos;
  130. }
  131. method cosec(Complex:D:) {
  132. 1e0 / self.sin;
  133. }
  134. method acosec(Complex:D:) {
  135. (1e0 / self).asin;
  136. }
  137. method cotan(Complex:D:) {
  138. self.cos / self.sin;
  139. }
  140. method acotan(Complex:D:) {
  141. (1e0 / self).atan;
  142. }
  143. method sinh(Complex:D:) {
  144. -((Complex.new(0e0, 1e0) * self).sin)i;
  145. }
  146. method asinh(Complex:D:) {
  147. (self + sqrt(1e0 + self * self)).log;
  148. }
  149. method cosh(Complex:D:) {
  150. (Complex.new(0e0, 1e0) * self).cos;
  151. }
  152. method acosh(Complex:D:) {
  153. (self + sqrt(self * self - 1e0)).log;
  154. }
  155. method tanh(Complex:D:) {
  156. -((Complex.new(0e0, 1e0) * self).tan)i;
  157. }
  158. method atanh(Complex:D:) {
  159. (((1e0 + self) / (1e0 - self)).log / 2e0);
  160. }
  161. method sech(Complex:D:) {
  162. 1e0 / self.cosh;
  163. }
  164. method asech(Complex:D:) {
  165. (1e0 / self).acosh;
  166. }
  167. method cosech(Complex:D:) {
  168. 1e0 / self.sinh;
  169. }
  170. method acosech(Complex:D:) {
  171. (1e0 / self).asinh;
  172. }
  173. method cotanh(Complex:D:) {
  174. 1e0 / self.tanh;
  175. }
  176. method acotanh(Complex:D:) {
  177. (1e0 / self).atanh;
  178. }
  179. method floor(Complex:D:) {
  180. Complex.new( self.re.floor, self.im.floor );
  181. }
  182. method ceiling(Complex:D:) {
  183. Complex.new( self.re.ceiling, self.im.ceiling );
  184. }
  185. proto method round(|) {*}
  186. multi method round(Complex:D:) {
  187. Complex.new( self.re.round, self.im.round );
  188. }
  189. multi method round(Complex:D: Real() $scale) {
  190. Complex.new( self.re.round($scale), self.im.round($scale) );
  191. }
  192. method truncate(Complex:D:) {
  193. Complex.new( self.re.truncate, self.im.truncate );
  194. }
  195. method narrow(Complex:D:) {
  196. self == 0e0 ?? 0 !!
  197. $!re == 0e0 ?? self !!
  198. $!im / $!re ≅ 0e0
  199. ?? $!re.narrow
  200. !! self;
  201. }
  202. }
  203. multi sub prefix:<->(Complex:D \a --> Complex:D) {
  204. my $new := nqp::create(Complex);
  205. nqp::bindattr_n( $new, Complex, '$!re',
  206. nqp::neg_n(
  207. nqp::getattr_n(nqp::decont(a), Complex, '$!re')
  208. )
  209. );
  210. nqp::bindattr_n( $new, Complex, '$!im',
  211. nqp::neg_n(
  212. nqp::getattr_n(nqp::decont(a), Complex, '$!im')
  213. )
  214. );
  215. $new;
  216. }
  217. multi sub abs(Complex:D \a --> Num:D) {
  218. my num $re = nqp::getattr_n(nqp::decont(a), Complex, '$!re');
  219. my num $im = nqp::getattr_n(nqp::decont(a), Complex, '$!im');
  220. nqp::p6box_n(nqp::sqrt_n(nqp::add_n(nqp::mul_n($re, $re), nqp::mul_n($im, $im))));
  221. }
  222. multi sub infix:<+>(Complex:D \a, Complex:D \b --> Complex:D) {
  223. my $new := nqp::create(Complex);
  224. nqp::bindattr_n( $new, Complex, '$!re',
  225. nqp::add_n(
  226. nqp::getattr_n(nqp::decont(a), Complex, '$!re'),
  227. nqp::getattr_n(nqp::decont(b), Complex, '$!re'),
  228. )
  229. );
  230. nqp::bindattr_n( $new, Complex, '$!im',
  231. nqp::add_n(
  232. nqp::getattr_n(nqp::decont(a), Complex, '$!im'),
  233. nqp::getattr_n(nqp::decont(b), Complex, '$!im'),
  234. )
  235. );
  236. $new;
  237. }
  238. multi sub infix:<+>(Complex:D \a, Num(Real) \b --> Complex:D) {
  239. my $new := nqp::create(Complex);
  240. nqp::bindattr_n( $new, Complex, '$!re',
  241. nqp::add_n(
  242. nqp::getattr_n(nqp::decont(a), Complex, '$!re'),
  243. nqp::unbox_n(b)
  244. )
  245. );
  246. nqp::bindattr_n($new, Complex, '$!im',
  247. nqp::getattr_n(nqp::decont(a), Complex, '$!im'),
  248. );
  249. $new
  250. }
  251. multi sub infix:<+>(Num(Real) \a, Complex:D \b --> Complex:D) {
  252. my $new := nqp::create(Complex);
  253. nqp::bindattr_n($new, Complex, '$!re',
  254. nqp::add_n(
  255. nqp::unbox_n(a),
  256. nqp::getattr_n(nqp::decont(b), Complex, '$!re'),
  257. )
  258. );
  259. nqp::bindattr_n($new, Complex, '$!im',
  260. nqp::getattr_n(nqp::decont(b), Complex, '$!im'),
  261. );
  262. $new;
  263. }
  264. multi sub infix:<->(Complex:D \a, Complex:D \b --> Complex:D) {
  265. my $new := nqp::create(Complex);
  266. nqp::bindattr_n( $new, Complex, '$!re',
  267. nqp::sub_n(
  268. nqp::getattr_n(nqp::decont(a), Complex, '$!re'),
  269. nqp::getattr_n(nqp::decont(b), Complex, '$!re'),
  270. )
  271. );
  272. nqp::bindattr_n($new, Complex, '$!im',
  273. nqp::sub_n(
  274. nqp::getattr_n(nqp::decont(a), Complex, '$!im'),
  275. nqp::getattr_n(nqp::decont(b), Complex, '$!im'),
  276. )
  277. );
  278. $new
  279. }
  280. multi sub infix:<->(Complex:D \a, Num(Real) \b --> Complex:D) {
  281. my $new := nqp::create(Complex);
  282. nqp::bindattr_n( $new, Complex, '$!re',
  283. nqp::sub_n(
  284. nqp::getattr_n(nqp::decont(a), Complex, '$!re'),
  285. b,
  286. )
  287. );
  288. nqp::bindattr_n($new, Complex, '$!im',
  289. nqp::getattr_n(nqp::decont(a), Complex, '$!im')
  290. );
  291. $new
  292. }
  293. multi sub infix:<->(Num(Real) \a, Complex:D \b --> Complex:D) {
  294. my $new := nqp::create(Complex);
  295. nqp::bindattr_n( $new, Complex, '$!re',
  296. nqp::sub_n(
  297. a,
  298. nqp::getattr_n(nqp::decont(b), Complex, '$!re'),
  299. )
  300. );
  301. nqp::bindattr_n($new, Complex, '$!im',
  302. nqp::neg_n(
  303. nqp::getattr_n(nqp::decont(b), Complex, '$!im')
  304. )
  305. );
  306. $new
  307. }
  308. multi sub infix:<*>(Complex:D \a, Complex:D \b --> Complex:D) {
  309. my num $a_re = nqp::getattr_n(nqp::decont(a), Complex, '$!re');
  310. my num $a_im = nqp::getattr_n(nqp::decont(a), Complex, '$!im');
  311. my num $b_re = nqp::getattr_n(nqp::decont(b), Complex, '$!re');
  312. my num $b_im = nqp::getattr_n(nqp::decont(b), Complex, '$!im');
  313. my $new := nqp::create(Complex);
  314. nqp::bindattr_n($new, Complex, '$!re',
  315. nqp::sub_n(nqp::mul_n($a_re, $b_re), nqp::mul_n($a_im, $b_im)),
  316. );
  317. nqp::bindattr_n($new, Complex, '$!im',
  318. nqp::add_n(nqp::mul_n($a_re, $b_im), nqp::mul_n($a_im, $b_re)),
  319. );
  320. $new;
  321. }
  322. multi sub infix:<*>(Complex:D \a, Num(Real) \b --> Complex:D) {
  323. my $new := nqp::create(Complex);
  324. my num $b_num = b;
  325. nqp::bindattr_n($new, Complex, '$!re',
  326. nqp::mul_n(
  327. nqp::getattr_n(nqp::decont(a), Complex, '$!re'),
  328. $b_num,
  329. )
  330. );
  331. nqp::bindattr_n($new, Complex, '$!im',
  332. nqp::mul_n(
  333. nqp::getattr_n(nqp::decont(a), Complex, '$!im'),
  334. $b_num,
  335. )
  336. );
  337. $new
  338. }
  339. multi sub infix:<*>(Num(Real) \a, Complex:D \b --> Complex:D) {
  340. my $new := nqp::create(Complex);
  341. my num $a_num = a;
  342. nqp::bindattr_n($new, Complex, '$!re',
  343. nqp::mul_n(
  344. $a_num,
  345. nqp::getattr_n(nqp::decont(b), Complex, '$!re'),
  346. )
  347. );
  348. nqp::bindattr_n($new, Complex, '$!im',
  349. nqp::mul_n(
  350. $a_num,
  351. nqp::getattr_n(nqp::decont(b), Complex, '$!im'),
  352. )
  353. );
  354. $new
  355. }
  356. multi sub infix:</>(Complex:D \a, Complex:D \b --> Complex:D) {
  357. my num $a_re = nqp::getattr_n(nqp::decont(a), Complex, '$!re');
  358. my num $a_im = nqp::getattr_n(nqp::decont(a), Complex, '$!im');
  359. my num $b_re = nqp::getattr_n(nqp::decont(b), Complex, '$!re');
  360. my num $b_im = nqp::getattr_n(nqp::decont(b), Complex, '$!im');
  361. my num $d = nqp::add_n(nqp::mul_n($b_re, $b_re), nqp::mul_n($b_im, $b_im));
  362. my $new := nqp::create(Complex);
  363. nqp::bindattr_n($new, Complex, '$!re',
  364. nqp::div_n(
  365. nqp::add_n(nqp::mul_n($a_re, $b_re), nqp::mul_n($a_im, $b_im)),
  366. $d,
  367. )
  368. );
  369. nqp::bindattr_n($new, Complex, '$!im',
  370. nqp::div_n(
  371. nqp::sub_n(nqp::mul_n($a_im, $b_re), nqp::mul_n($a_re, $b_im)),
  372. $d,
  373. )
  374. );
  375. $new;
  376. }
  377. multi sub infix:</>(Complex:D \a, Real \b --> Complex:D) {
  378. Complex.new(a.re / b, a.im / b);
  379. }
  380. multi sub infix:</>(Real \a, Complex:D \b --> Complex:D) {
  381. Complex.new(a, 0e0) / b;
  382. }
  383. multi sub infix:<**>(Complex:D \a, Complex:D \b --> Complex:D) {
  384. (a.re == 0e0 && a.im == 0e0)
  385. ?? ( b.re == 0e0 && b.im == 0e0
  386. ?? Complex.new(1e0, 0e0)
  387. !! Complex.new(0e0, 0e0)
  388. )
  389. !! (b * a.log).exp
  390. }
  391. multi sub infix:<**>(Num(Real) \a, Complex:D \b --> Complex:D) {
  392. a == 0e0
  393. ?? ( b.re == 0e0 && b.im == 0e0
  394. ?? Complex.new(1e0, 0e0)
  395. !! Complex.new(0e0, 0e0)
  396. )
  397. !! (b * a.log).exp
  398. }
  399. multi sub infix:<**>(Complex:D \a, Num(Real) \b --> Complex:D) {
  400. b == 0e0 ?? Complex.new(1e0, 0e0) !! (b * a.log).exp
  401. }
  402. multi sub infix:<==>(Complex:D \a, Complex:D \b --> Bool:D) { a.re == b.re && a.im == b.im }
  403. multi sub infix:<==>(Complex:D \a, Num(Real) \b --> Bool:D) { a.re == b && a.im == 0e0 }
  404. multi sub infix:<==>(Num(Real) \a, Complex:D \b --> Bool:D) { a == b.re && 0e0 == b.im }
  405. multi sub infix:<===>(Complex:D \a, Complex:D \b --> Bool:D) {
  406. a.WHAT =:= b.WHAT && a.re === b.re && a.im === b.im
  407. }
  408. multi sub infix:<≅>(Complex:D \a, Complex:D \b --> Bool:D) { a.re ≅ b.re && a.im ≅ b.im || a <=> b =:= Same }
  409. multi sub infix:<≅>(Complex:D \a, Num(Real) \b --> Bool:D) { a ≅ b.Complex }
  410. multi sub infix:<≅>(Num(Real) \a, Complex:D \b --> Bool:D) { a.Complex ≅ b }
  411. # Meaningful only for sorting purposes, of course.
  412. # We delegate to Real::cmp rather than <=> because parts might be NaN.
  413. multi sub infix:<cmp>(Complex:D \a, Complex:D \b --> Order:D) { a.re cmp b.re || a.im cmp b.im }
  414. multi sub infix:<cmp>(Num(Real) \a, Complex:D \b --> Order:D) { a cmp b.re || 0 cmp b.im }
  415. multi sub infix:<cmp>(Complex:D \a, Num(Real) \b --> Order:D) { a.re cmp b || a.im cmp 0 }
  416. multi sub infix:«<=>»(Complex:D \a, Complex:D \b --> Order:D) {
  417. my $tolerance = a && b
  418. ?? (a.re.abs + b.re.abs) / 2 * $*TOLERANCE # Scale slop to average real parts.
  419. !! $*TOLERANCE; # Don't want tolerance 0 if either arg is 0.
  420. # Fail unless imaginary parts are relatively negligible, compared to real parts.
  421. infix:<≅>(a.im, 0e0, :$tolerance) && infix:<≅>(b.im, 0e0, :$tolerance)
  422. ?? a.re <=> b.re
  423. !! Failure.new(X::Numeric::Real.new(target => Real, reason => "Complex is not numerically orderable", source => "Complex"))
  424. }
  425. multi sub infix:«<=>»(Num(Real) \a, Complex:D \b --> Order:D) { a.Complex <=> b }
  426. multi sub infix:«<=>»(Complex:D \a, Num(Real) \b --> Order:D) { a <=> b.Complex }
  427. proto sub postfix:<i>(\a --> Complex:D) is pure {*}
  428. multi sub postfix:<i>(Real \a --> Complex:D) { Complex.new(0e0, a); }
  429. multi sub postfix:<i>(Complex:D \a --> Complex:D) { Complex.new(-a.im, a.re) }
  430. multi sub postfix:<i>(Numeric \a --> Complex:D) { a * Complex.new(0e0, 1e0) }
  431. multi sub postfix:<i>(Cool \a --> Complex:D) { a.Numeric * Complex.new(0e0, 1e0) }
  432. constant i = Complex.new(0e0, 1e0);