Some accumulated math changes

5 files changed, 54 insertions, 15 deletions

diff --git a/std/Generic/Number.bruijn b/std/Generic/Number.bruijn
index e39fc91..08effa6 100644
--- a/std/Generic/Number.bruijn
+++ b/std/Generic/Number.bruijn
@@ -21,17 +21,17 @@ not-eq? not! ∘∘ eq? ⧗ Generic → Generic → Boolean
 # prefix for comparing functions
 ?‣ &eq?
 
-# returns true if number is lts than other number
+# returns true if number is less than other number
 lt? \gt? ⧗ Generic → Generic → Boolean
 
 …<?… lt?
 
-# returns true if number is lts than or equal to other number
+# returns true if number is less than or equal to other number
 le? not! ∘∘ gt? ⧗ Generic → Generic → Boolean
 
 …≤?… le?
 
-# returns true if number is gtater than or equal to other number
+# returns true if number is greater than or equal to other number
 ge? \le? ⧗ Generic → Generic → Boolean
 
 …≥?… ge?
@@ -47,10 +47,10 @@ compare compare-case (+0) (+1) (-1) ⧗ Generic → Generic → Generic
 # returns true if comparison result is equal (EQ)
 c-eq? eq? (+0) ⧗ Generic → Generic
 
-# returns true if comparison result is lts than (LT)
+# returns true if comparison result is less than (LT)
 c-lt? eq? (-1) ⧗ Generic → Generic
 
-# returns true if comparison result is gtater than (GT)
+# returns true if comparison result is greater than (GT)
 c-gt? eq? (+1) ⧗ Generic → Generic
 
 # returns max number of two
diff --git a/std/Math.bruijn b/std/Math.bruijn
index 3e32bef..84cfae2 100644
--- a/std/Math.bruijn
+++ b/std/Math.bruijn
@@ -1,5 +1,6 @@
 # MIT License, Copyright (c) 2022 Marvin Borner
 # experimental functions; sometimes list-based; could work on any base
+# TODO: some functions should be moved to respective bases
 
 :input std/Number
 
@@ -60,19 +61,24 @@ product L.foldl mul (+1) ⧗ (List Number) → Number
 
 :test (∏ (+1) → (+3) | ++‣) ((+24))
 
-# greatest common divisor using repeated subtraction
-gcd* z [[[1 =? 0 case-eq (1 >? 0 case-gre case-les)]]] ⧗ Number → Number → Number
-	case-eq 1
-	case-gre 2 (1 - 0) 0
-	case-les 2 1 (0 - 1)
-
-# greatest common divisor using modulo (mostly faster than gcd*)
-gcd z [[[=?0 1 (2 0 (1 % 0))]]] ⧗ Number → Number → Number
+# greatest common divisor using repeated subtraction enhanced by ternary shifts
+# (temporary)
+gcd z [[[(1 =? 0) 1 (=?1 0 (=?0 1 else))]]] ⧗ Number → Number → Number
+	else [[(1 ⋀? 0) ((+3) ⋅ (4 /³3 /³2)) (1 (4 /³3 2) (0 (4 3 /³2) else))]] (t⁰? (lst 1)) (t⁰? (lst 0))
+		else 3 >? 2 (4 (3 - 2) 2) (4 3 (2 - 3))
 
 :test ((gcd (+2) (+4)) =? (+2)) (true)
 :test ((gcd (+10) (+5)) =? (+5)) (true)
 :test ((gcd (+3) (+8)) =? (+1)) (true)
 
+# greatest common divisor using modulo (mostly slower than gcd)
+# TODO: would be faster if ternary quot-rem was efficient!
+gcd* z [[[=?0 1 (2 0 (1 % 0))]]] ⧗ Number → Number → Number
+
+:test ((gcd* (+2) (+4)) =? (+2)) (true)
+:test ((gcd* (+10) (+5)) =? (+5)) (true)
+:test ((gcd* (+3) (+8)) =? (+1)) (true)
+
 # least common multiple using gcd
 lcm [[=?1 1 (=?0 0 |(1 / (gcd 1 0) ⋅ 0))]] ⧗ Number → Number → Number
 
@@ -101,7 +107,9 @@ pow* z [[[rec]]] ⧗ Number → Number → Number
 		case-end (+1)
 
 # factorial function
-fac [∏ (+1) → 0 | i] ⧗ Number → Number
+# fac [∏ (+1) → 0 | i] ⧗ Number → Number
+
+fac [0 0] [[=?0 (+1) (0 ⋅ (1 1 --0))]] ⧗ Number → Number
 
 :test ((fac (+3)) =? (+6)) (true)
 
@@ -226,6 +234,9 @@ primes L.map fst (L.filter snd (enumerate characteristic-primes)) ⧗ (List Numb
 # slower but cooler prime number sequence
 primes* L.nub ((…≠?… (+1)) ∘∘ gcd) (L.iterate ++‣ (+2)) ⧗ (List Number)
 
+# primality test
+prime? L.index characteristic-primes ⧗ Number → Boolean
+
 # prime factors
 factors \divs primes ⧗ Number → (List Number)
 	divs y [[&[[&[[3 ⋅ 3 >? 4 case-1 (=?0 case-2 case-3)]] (quot-rem 2 1)]]]]
diff --git a/std/Math/Rational.bruijn b/std/Math/Rational.bruijn
index f1509b0..e72b6a9 100644
--- a/std/Math/Rational.bruijn
+++ b/std/Math/Rational.bruijn
@@ -43,6 +43,12 @@ add &[[&[[p : q]]]] ⧗ Rational → Rational → Rational
 :test ((+1.8) + (+1.2) =? (+3.0)) (true)
 :test ((-1.8) + (+1.2) =? (-0.6)) (true)
 
+reduce [[[(N.div 2 0) : N.--(N.div 1 0)] (N.gcd 1 0)]]
+
+add' &[[&[[reduce p q]]]]
+	p N.add (N.mul 3 N.++0) (N.mul 1 N.++2)
+	q N.mul N.++0 N.++2
+
 # subtracts two rational numbers
 sub &[[&[[p : q]]]] ⧗ Rational → Rational → Rational
 	p N.sub (N.mul 3 N.++0) (N.mul 1 N.++2)
diff --git a/std/Math/Real.bruijn b/std/Math/Real.bruijn
index 4b3e7d0..9a4b1bd 100644
--- a/std/Math/Real.bruijn
+++ b/std/Math/Real.bruijn
@@ -139,6 +139,12 @@ sqrt [[y [[[N.=?0 guess go]]] (1 0) 0]] ⧗ Real → Real
 	guess (+1.0q)
 	go [Q.div (Q.add 0 (Q.div 2 0)) (+2.0q)] (2 1 N.--0)
 
+# Newton's/Heron's method, quadratic convergence
+# tex: x_{n+1}=\frac{3x_n+a/(x_n^2)}{4}
+cbrt [[y [[[N.=?0 guess go]]] (1 0) 0]] ⧗ Real → Real
+	guess (+1.0q)
+	go [Q.div (Q.add (Q.mul (+3.0q) 0) (Q.div 2 (Q.mul 0 0))) (+4.0q)] (2 1 N.--0)
+
 # hypotenuse
 hypot [[sqrt ((0 ⋅ 0) + (1 ⋅ 1))]] ⧗ Real → Real → Real
 
@@ -164,6 +170,18 @@ hypot [[sqrt ((0 ⋅ 0) + (1 ⋅ 1))]] ⧗ Real → Real → Real
 		p Q.mul (+2.0q) 2
 	final [[[[Q.div (Q.pow-n (Q.add 3 2) (+2)) (Q.mul (+4.0q) 1)]]]]
 
+# golden ratio from direct formula
+φ ++(sqrt (+5.0r)) / (+2.0r)
+
+# conjugate golden ratio
+ψ -(~φ)
+
+# golden ratio from fibonacci convergence
+φ* [(L.index (L.iterate &[[0 : (N.add 1 0)]] ((+0) : (+1))) 0) [[1 : N.--0]]]
+
+# real fibonacci
+fib [((pow φ 0) - (pow ψ 0)) / (sqrt (+5.0r))]
+
 # arctan by Taylor expansion, only for |x|<=1
 # tex: \sum_{n=0}^\infty(-1)^n \frac{x^{2n+1}}{2n+1}
 arctan* [[[L.nth-iterate &[[[[op]]]] start 1] (1 0) [[[[3]]]]]] ⧗ Real → Real
diff --git a/std/Number/Unary.bruijn b/std/Number/Unary.bruijn
index b764e87..214bfb9 100644
--- a/std/Number/Unary.bruijn
+++ b/std/Number/Unary.bruijn
@@ -34,7 +34,7 @@ mod [[[[3 &k (3 [3 [[[0 (2 (5 1)) 1]]] [1] 1] [1]) ki]]]] ⧗ Unary → Unary �
 :test ((+3u) % (+5u)) ((+3u))
 
 # returns true if the number is even (remainder mod 2 == 0)
-even? [=?(0 % (+2u))] ⧗ Unary → Boolean
+even? [0 not! true] ⧗ Unary → Boolean
 
 =²?‣ even?
 
@@ -161,3 +161,7 @@ fac [[1 [[0 (1 [[2 1 (1 0)]])]] [1] i]] ⧗ Unary → Unary
 hyperfac [[1 [[(0 0) (1 [[2 1 (1 0)]])]] [1] i]] ⧗ Unary → Unary
 
 :test (hyperfac (+3u)) ((+108u))
+
+# Wilson's theorem
+# very inefficient but very golfable
+prime? [=?(++(fac --0) % 0)] ⧗ Unary → Boolean