Skip to content

Instantly share code, notes, and snippets.

@grondilu
Last active September 26, 2023 19:18
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save grondilu/b70fb0d8abac36221f0d763f3a7317c6 to your computer and use it in GitHub Desktop.
Save grondilu/b70fb0d8abac36221f0d763f3a7317c6 to your computer and use it in GitHub Desktop.
quadratic irrationals in raku
unit class QuadraticIrrational does Real;
#
# This role cannot meaningfully be converted
# into a core type, so we'll "burn the bridge"
method Bridge {!!!}
has UInt $.discriminant;
has Rat(Cool) ($.u, $.v);
submethod TWEAK {
unless $!v == 0 {
for ^Inf {
next unless .is-prime;
last if $_**2 > $!discriminant;
if $!discriminant %% $_**2 {
$!discriminant div= $_**2;
$!v *= $_;
}
}
fail "discriminant must be at least two" unless $!discriminant > 1;
fail "discriminant must not be the square of an integer" if $!discriminant.sqrt.Int**2 == $!discriminant;
}
}
multi method gist { "$!u + $!v*√$!discriminant" }
multi method new(Rat(Cool) $r, :$discriminant) { self.new: u => $r, v => 0, :$discriminant }
sub prefix:<√>(UInt $N where $N > 1 --> ::?CLASS) is export { ::?CLASS.new: u => 0, v => 1, discriminant => $N }
proto method add (Real $ --> ::?CLASS) {*}
proto method subtract(Real $ --> ::?CLASS) {*}
proto method multiply(Real $ --> ::?CLASS) {*}
proto method divide (Real $ --> ::?CLASS) {*}
method negate { self.new: :u(-$!u), :v(-$!v), :$!discriminant }
method inverse {
my $d = $!u**2 - $!v**2*$!discriminant;
self.new: u => $!u/$d, v => -$!v/$d, :$!discriminant;
}
method sign {
if $!u & $!v ≥ 0 { return +1 }
elsif $!u & $!v ≤ 0 { return -1 }
elsif $!u & $!v == 0 { return 0 }
elsif $!u | $!v == 0 { return sign($!v) + sign($!v) }
else { return sign($!u)*sign($!u**2 - $!v**2*$!discriminant) }
}
method abs { self.sign < 0 ?? self.negate !! self }
multi method add($A: ::?CLASS $B where $A.discriminant == $B.discriminant) { $A.new: u => $A.u + $B.u, v => $A.v + $B.v, discriminant => $A.discriminant }
multi method add(::?CLASS $) {...}
multi method add($A: $x) { $A.new: u => $A.u + $x, v => $A.v, discriminant => $A.discriminant }
multi method subtract($A: ::?CLASS $B where $A.discriminant == $B.discriminant) { $A.new: u => $A.u - $B.u, v => $A.v - $B.v, discriminant => $A.discriminant }
multi method subtract(::?CLASS $) {...}
multi method subtract($A: $x) { $A.new: u => $A.u - $x, v => $A.v, discriminant => $A.discriminant }
multi method multiply($A: ::?CLASS $B where $A.discriminant == $B.discriminant) {
self.new:
u => $A.u*$B.u + $A.v*$B.v*$A.discriminant,
v => $A.u*$B.v + $A.v*$B.u,
discriminant => $A.discriminant
}
multi method multiply($A: ::?CLASS $B where $A.u & $B.u == 0) {
self.new: u => 0,
v => $A.v*$B.v,
discriminant => $A.discriminant*$B.discriminant;
}
multi method multiply(::?CLASS $) {...}
multi method multiply(Real $x) {
self.new: u => $!u*$x, v => $!v*$x, :$!discriminant
}
multi method divide(::?CLASS $) {...}
multi method divide($A: ::?CLASS $B) { $A.multiply: $B.inverse }
multi method divide(Real $x) {
self.new: u => $!u/$x, v => $!v/$x, :$!discriminant
}
multi prefix:<+>(::?CLASS $x) is export { $x }
multi prefix:<->(::?CLASS $x) is export { $x.negate }
multi infix:<+>(::?CLASS $A, Real $B --> ::?CLASS) is export { $A.add: $B }
multi infix:<+>(Real $A, ::?CLASS $B --> ::?CLASS) is export { $B.add: $A }
multi infix:<+>(::?CLASS $A, ::?CLASS $B --> ::?CLASS) is export { $A.add: $B }
multi infix:<->(::?CLASS $A, Real $B --> ::?CLASS) is export { $A.subtract: $B }
multi infix:<->(Real $A, ::?CLASS $B --> ::?CLASS) is export { $B.negate.add: $A }
multi infix:<->(::?CLASS $A, ::?CLASS $B --> ::?CLASS) is export { $A.subtract: $B }
multi infix:<*>(::?CLASS $A, Real $B --> ::?CLASS) is export { $A.multiply: $B }
multi infix:<*>(Real $A, ::?CLASS $B --> ::?CLASS) is export { $B.multiply: $A }
multi infix:</>(::?CLASS $A, Real $B --> ::?CLASS) is export { $A.divide: $B }
multi infix:</>(1, ::?CLASS $B --> ::?CLASS) is export { $B.inverse }
multi infix:<**>(::?CLASS $A, 0 --> ::?CLASS) is export {
fail "0**0 is undeterminate" if $A.u & $A.v == 0;
$A.new: :u(1), :v(0)
}
multi infix:<**>(::?CLASS $A, 1 --> ::?CLASS) is export { $A }
multi infix:<**>(::?CLASS $A, 2 --> ::?CLASS) is export { $A*$A }
multi infix:<**>(::?CLASS $A, UInt $i --> ::?CLASS) is export {
($A**2)**($i div 2) * $A**($i mod 2)
}
multi sub infix:«<=>»(::?CLASS $x, 0 --> Order) is export {
my ($u, $v) = $x.u, $x.v;
if $u & $v > 0 { return More }
elsif $u & $v < 0 { return Less }
elsif $u & $v == 0 { return Same }
elsif $u == 0 { return $v <=> 0 }
elsif $v == 0 { return $u <=> 0 }
elsif $u < 0 { $v**2*$x.discriminant <=> $u**2 }
else { $u**2 <=> $v**2*$x.discriminant }
}
multi infix:«<=>»(::?CLASS $a, ::?CLASS $b --> Order) is export { $a - $b <=> 0 }
multi infix:«<=>»(::?CLASS $A, ::?CLASS $B where $A.discriminant == $B.discriminant) is export {
given ($A - $B).sign {
when -1 { return Less }
when 0 { return Same }
when 1 { return More }
}
}
multi infix:«<=>»(Real $r, ::?CLASS $B) is export {
samewith $B.new($r, discriminant => $B.discriminant), $B
}
multi infix:«<=>»(::?CLASS $A, Real $r) is export {
samewith $A, $A.new($r, discriminant => $A.discriminant)
}
multi infix:«>»(::?CLASS $A, Real $r) is export { $A <=> $r ~~ More }
multi infix:«<»(::?CLASS $A, Real $r) is export { $A <=> $r ~~ Less }
multi infix:«≤»(::?CLASS $A, Real $r) is export { $A <=> $r ~~ Less|Same }
multi infix:«≥»(::?CLASS $A, Real $r) is export { $A <=> $r ~~ More|Same }
multi sub infix:<==>(::?CLASS $x, 0) is export { $x.sign == 0 }
multi sub infix:<==>(::?CLASS $x, Inf) is export { False }
multi infix:<==>(::?CLASS $A, Real $x) is export {
if $x ~~ ::?CLASS { return $x <=> $A ~~ Same }
else { return False }
}
multi infix:<==>(Real $x, ::?CLASS $B) is export { samewith $B, $x }
method floor returns Int {
if self > 0 {
my $max-power = (1, 2, 4, 8 ...^ self < *).tail;
return $max-power ?? $max-power + samewith(self - $max-power.Rat) !! 0;
} else { -samewith(-self) - 1 }
}
method ceiling returns Int { self.floor + 1 }
CHECK {
.say for
√2,
√2 + √2,
√2 + 1,
-√2,
3*√2,
5*√2 + 5,
√3 / 57,
√2 - 0.2*√2,
√2 <=> 1.415,
1 + √2,
1/(1 + √2),
(√31).floor,
31.sqrt,
(-√31).floor,
√2 == 3;
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment