RationalComputability.lean

1 import RealizationBoundaries.PrimeComputability
2 import Mathlib.FieldTheory.Finite.Basic
3 import Mathlib.Data.Rat.Lemmas
4 import Mathlib.Tactic.Ring
5 import Mathlib.Tactic.FieldSimp
6 
7 /-! # Native rational-code arithmetic and total prime reduction -/
8 set_option autoImplicit false
9 noncomputable section
10 namespace RealizationBoundaries
11 
12 /-- Positive-minus-negative numerator and positive denominator predecessor.
13 The carrier is an ordinary native primitive-codable product of naturals. -/
14 abbrev RationalCode := (ℕ × ℕ) × ℕ
15 namespace RationalCode
16 
17 def positive (c : RationalCode) : ℕ := c.1.1
18 def negative (c : RationalCode) : ℕ := c.1.2
19 def denominator (c : RationalCode) : ℕ := c.2 + 1
20 def value (c : RationalCode) : ℚ :=
21   mkRat ((positive c : ℤ) - (negative c : ℤ)) (denominator c)
22 
23 theorem value_eq (c : RationalCode) :
24     value c = ((positive c : ℚ) - (negative c : ℚ)) / (denominator c : ℚ) := by
25   simp only [value, Rat.mkRat_eq_div, Int.cast_sub, Int.cast_natCast]
26 
27 theorem denominator_pos (c : RationalCode) : 0 < denominator c := Nat.succ_pos _
28 
29 theorem value_surjective : Function.Surjective value := by
30   intro q
31   refine ⟨((q.num.toNat, (-q.num).toNat), q.den - 1), ?_⟩
32   unfold value positive negative denominator
33   rw [Int.toNat_sub_toNat_neg, Nat.sub_add_cancel q.den_pos, Rat.mkRat_self]
34 
35 def add (c e : RationalCode) : RationalCode :=
36   ((positive c * denominator e + positive e * denominator c,
37     negative c * denominator e + negative e * denominator c),
38     denominator c * denominator e - 1)
39 
40 def neg (c : RationalCode) : RationalCode := ((negative c, positive c), c.2)
41 
42 @[simp] theorem value_add (c e : RationalCode) : value (add c e) = value c + value e := by
43   rw [value_eq, value_eq, value_eq]
44   have hd : 0 < denominator c * denominator e := Nat.mul_pos (denominator_pos c) (denominator_pos e)
45   simp only [add, positive, negative, denominator]
46   simp only [denominator] at hd
47   rw [Nat.sub_add_cancel hd]
48   push_cast
49   have hc : (c.2 + 1 : ℚ) ≠ 0 := by positivity
50   have he : (e.2 + 1 : ℚ) ≠ 0 := by positivity
51   field_simp
52   ring
53 
54 @[simp] theorem value_neg (c : RationalCode) : value (neg c) = -value c := by
55   simp only [value_eq, neg, positive, negative, denominator]
56   ring
57 
58 /-- Cross multiplication decides equality of the decoded rational values. -/
59 def equal (c e : RationalCode) : Bool := decide
60   (positive c * denominator e + negative e * denominator c =
61     positive e * denominator c + negative c * denominator e)
62 
63 theorem equal_iff (c e : RationalCode) : equal c e = true ↔ value c = value e := by
64   simp only [equal, decide_eq_true_eq, value_eq]
65   have hc : (denominator c : ℚ) ≠ 0 := by exact_mod_cast (denominator_pos c).ne'
66   have he : (denominator e : ℚ) ≠ 0 := by exact_mod_cast (denominator_pos e).ne'
67   rw [div_eq_div_iff hc he]
68   constructor
69   · intro h
70     have h' : (positive c : ℚ) * denominator e + negative e * denominator c =
71         positive e * denominator c + negative c * denominator e := by exact_mod_cast h
72     nlinarith
73   · intro h
74     have h' : (positive c : ℚ) * denominator e + negative e * denominator c =
75         positive e * denominator c + negative c * denominator e := by nlinarith
76     exact_mod_cast h'
77 
78 theorem primrec_positive : Primrec positive := Primrec.fst.comp Primrec.fst
79 theorem primrec_negative : Primrec negative := Primrec.snd.comp Primrec.fst
80 theorem primrec_denominator : Primrec denominator := Primrec.succ.comp Primrec.snd
81 
82 theorem primrec_add : Primrec₂ add := by
83   let a : Primrec (fun p : RationalCode × RationalCode => positive p.1) := primrec_positive.comp Primrec.fst
84   let b : Primrec (fun p : RationalCode × RationalCode => negative p.1) := primrec_negative.comp Primrec.fst
85   let d : Primrec (fun p : RationalCode × RationalCode => denominator p.1) := primrec_denominator.comp Primrec.fst
86   let e : Primrec (fun p : RationalCode × RationalCode => positive p.2) := primrec_positive.comp Primrec.snd
87   let f : Primrec (fun p : RationalCode × RationalCode => negative p.2) := primrec_negative.comp Primrec.snd
88   let g : Primrec (fun p : RationalCode × RationalCode => denominator p.2) := primrec_denominator.comp Primrec.snd
89   exact ((Primrec.nat_add.comp (Primrec.nat_mul.comp a g) (Primrec.nat_mul.comp e d)).pair
90     (Primrec.nat_add.comp (Primrec.nat_mul.comp b g) (Primrec.nat_mul.comp f d))).pair
91     (Primrec.nat_sub.comp (Primrec.nat_mul.comp d g) (Primrec.const 1))
92 
93 theorem primrec_neg : Primrec neg :=
94   (primrec_negative.pair primrec_positive).pair Primrec.snd
95 
96 theorem primrec_equal : Primrec₂ equal := by
97   let a : Primrec (fun p : RationalCode × RationalCode => positive p.1) := primrec_positive.comp Primrec.fst
98   let b : Primrec (fun p : RationalCode × RationalCode => negative p.1) := primrec_negative.comp Primrec.fst
99   let d : Primrec (fun p : RationalCode × RationalCode => denominator p.1) := primrec_denominator.comp Primrec.fst
100   let e : Primrec (fun p : RationalCode × RationalCode => positive p.2) := primrec_positive.comp Primrec.snd
101   let f : Primrec (fun p : RationalCode × RationalCode => negative p.2) := primrec_negative.comp Primrec.snd
102   let g : Primrec (fun p : RationalCode × RationalCode => denominator p.2) := primrec_denominator.comp Primrec.snd
103   exact (Primrec.eq.comp
104     (Primrec.nat_add.comp (Primrec.nat_mul.comp a g) (Primrec.nat_mul.comp f d))
105     (Primrec.nat_add.comp (Primrec.nat_mul.comp e d) (Primrec.nat_mul.comp b g))).decide
106 
107 end RationalCode
108 
109 /-- Euclid's gcd is native primitive recursive: bounded common-divisor search
110 returns the same greatest common divisor. The both-zero case is retained. -/
111 theorem primrec_nat_gcd : Primrec₂ Nat.gcd := by
112   have hpred : PrimrecRel (fun p : ℕ × ℕ => fun d : ℕ => d ∣ p.1 ∧ d ∣ p.2) :=
113     (primrec_nat_dvd.comp Primrec.snd (Primrec.fst.comp Primrec.fst)).and
114       (primrec_nat_dvd.comp Primrec.snd (Primrec.snd.comp Primrec.fst))
115   refine (Primrec.nat_findGreatest Primrec.nat_add hpred).of_eq ?_
116   intro p
117   apply Nat.findGreatest_eq_iff.mpr
118   refine ⟨?_, fun _ => ⟨Nat.gcd_dvd_left _ _, Nat.gcd_dvd_right _ _⟩, ?_⟩
119   · by_cases hb : p.2 = 0
120     · simp [hb]
121     · exact (Nat.gcd_le_right p.1 (Nat.pos_of_ne_zero hb)).trans (Nat.le_add_left _ _)
122   · intro n hlt hle hn
123     have hdn : n ∣ Nat.gcd p.1 p.2 := Nat.dvd_gcd hn.1 hn.2
124     by_cases hg : Nat.gcd p.1 p.2 = 0
125     · have hz := Nat.gcd_eq_zero_iff.mp hg
126       have hsum : p.1 + p.2 = 0 := by simp [hz.1, hz.2]
127       change n ≤ p.1 + p.2 at hle
128       omega
129     · exact (not_lt_of_ge (Nat.le_of_dvd (Nat.pos_of_ne_zero hg) hdn)) hlt
130 
131 namespace RationalCode
132 
133 def magnitude (c : RationalCode) : ℕ :=
134   (positive c - negative c) + (negative c - positive c)
135 def divisor (c : RationalCode) : ℕ := Nat.gcd (denominator c) (magnitude c)
136 def normalPositive (c : RationalCode) : ℕ := (positive c - negative c) / divisor c
137 def normalNegative (c : RationalCode) : ℕ := (negative c - positive c) / divisor c
138 def normalDenominator (c : RationalCode) : ℕ := denominator c / divisor c
139 
140 theorem primrec_magnitude : Primrec magnitude :=
141   Primrec.nat_add.comp (Primrec.nat_sub.comp primrec_positive primrec_negative)
142     (Primrec.nat_sub.comp primrec_negative primrec_positive)
143 theorem primrec_divisor : Primrec divisor :=
144   primrec_nat_gcd.comp primrec_denominator primrec_magnitude
145 theorem primrec_normalPositive : Primrec normalPositive :=
146   Primrec.nat_div.comp (Primrec.nat_sub.comp primrec_positive primrec_negative) primrec_divisor
147 theorem primrec_normalNegative : Primrec normalNegative :=
148   Primrec.nat_div.comp (Primrec.nat_sub.comp primrec_negative primrec_positive) primrec_divisor
149 theorem primrec_normalDenominator : Primrec normalDenominator :=
150   Primrec.nat_div.comp primrec_denominator primrec_divisor
151 
152 /-- The native normalized fields are the literal reduced Rat numerator and
153 denominator, including zero and negative numerators. -/
154 theorem normal_fields (c : RationalCode) :
155     (value c).num = (normalPositive c : ℤ) - (normalNegative c : ℤ) ∧
156     (value c).den = normalDenominator c := by
157   by_cases h : negative c ≤ positive c
158   · have hsub : negative c - positive c = 0 := Nat.sub_eq_zero_of_le h
159     have he : (positive c : ℤ) - (negative c : ℤ) = ((positive c - negative c : ℕ) : ℤ) :=
160       (Int.ofNat_sub h).symm
161     simp only [value, he, Rat.num_mkRat, Rat.den_mkRat,
162       (denominator_pos c).ne', ↓reduceIte, Int.natAbs_natCast,
163       normalPositive, normalNegative, normalDenominator, divisor, magnitude,
164       hsub, Nat.add_zero, Nat.zero_div, Int.natCast_ediv, Nat.cast_zero, sub_zero]
165     exact ⟨trivial, trivial⟩
166   · have hle : positive c ≤ negative c := Nat.le_of_lt (Nat.lt_of_not_ge h)
167     have hsub : positive c - negative c = 0 := Nat.sub_eq_zero_of_le hle
168     have he : (positive c : ℤ) - (negative c : ℤ) = -((negative c - positive c : ℕ) : ℤ) := by
169       rw [Int.ofNat_sub hle]
170       ring
171     have hdvd : (divisor c : ℤ) ∣ ((negative c - positive c : ℕ) : ℤ) := by
172       apply Int.ofNat_dvd.mpr
173       simpa [divisor, magnitude, hsub] using Nat.gcd_dvd_right (denominator c) (negative c - positive c)
174     simp only [value, he, Rat.num_mkRat, Rat.den_mkRat,
175       (denominator_pos c).ne', ↓reduceIte, Int.natAbs_neg, Int.natAbs_natCast,
176       normalPositive, normalNegative, normalDenominator, divisor, magnitude,
177       hsub, Nat.zero_add, Nat.zero_div, Nat.cast_zero, zero_sub]
178     constructor
179     · have hdvd' : ((denominator c).gcd (negative c - positive c) : ℤ) ∣
180           ((negative c - positive c : ℕ) : ℤ) := by
181         simpa [divisor, magnitude, hsub] using hdvd
182       rw [Int.neg_ediv_of_dvd hdvd', Int.natCast_ediv]
183     · trivial
184 
185 /-- A nonnegative representative of the reduced signed numerator modulo `p`. -/
186 def numeratorResidue (c : RationalCode) (p : ℕ) : ℕ :=
187   (normalPositive c % p + p - normalNegative c % p) % p
188 
189 /-- Fermat inversion is a total natural-number operation. The zero-denominator
190 branch is explicit and agrees with the Rat-to-ZMod cast at denominator primes. -/
191 def primeResidue (c : RationalCode) (p : ℕ) : ℕ :=
192   if normalDenominator c % p = 0 then 0
193   else (numeratorResidue c p * normalDenominator c ^ (p - 2)) % p
194 
195 theorem primrec_numeratorResidue : Primrec₂ numeratorResidue := by
196   let a : Primrec (fun z : RationalCode × ℕ => normalPositive z.1) :=
197     primrec_normalPositive.comp Primrec.fst
198   let b : Primrec (fun z : RationalCode × ℕ => normalNegative z.1) :=
199     primrec_normalNegative.comp Primrec.fst
200   exact Primrec.nat_mod.comp
201     (Primrec.nat_sub.comp (Primrec.nat_add.comp (Primrec.nat_mod.comp a Primrec.snd) Primrec.snd)
202       (Primrec.nat_mod.comp b Primrec.snd)) Primrec.snd
203 
204 theorem primrec_primeResidue : Primrec₂ primeResidue := by
205   have hpow : Primrec₂ (fun a b : ℕ => a ^ b) := Primrec₂.unpaired'.mp Nat.Primrec.pow
206   let d : Primrec (fun z : RationalCode × ℕ => normalDenominator z.1) :=
207     primrec_normalDenominator.comp Primrec.fst
208   exact Primrec.ite
209     (Primrec.eq.comp (Primrec.nat_mod.comp d Primrec.snd) (Primrec.const 0))
210     (Primrec.const 0)
211     (Primrec.nat_mod.comp
212       (Primrec.nat_mul.comp primrec_numeratorResidue
213         (hpow.comp d (Primrec.nat_sub.comp Primrec.snd (Primrec.const 2)))) Primrec.snd)
214 
215 /-- The natural numerator residue casts to the literal normalized Rat numerator. -/
216 theorem numeratorResidue_cast (c : RationalCode) {p : ℕ} (hp : p.Prime) :
217     (numeratorResidue c p : ZMod p) = ((value c).num : ZMod p) := by
218   have hb : normalNegative c % p ≤ normalPositive c % p + p :=
219     (Nat.mod_lt _ hp.pos).le.trans (Nat.le_add_left _ _)
220   rw [numeratorResidue, ZMod.natCast_mod, Nat.cast_sub hb, Nat.cast_add,
221     ZMod.natCast_mod, ZMod.natCast_mod, ZMod.natCast_self, add_zero, (normal_fields c).1]
222   push_cast
223   rfl
224 
225 /-- Exact reduction at every prime, including all exceptional denominator primes.
226 This is the total Rat cast used by the actual rational sections. -/
227 theorem primeResidue_cast (c : RationalCode) {p : ℕ} [Fact p.Prime] :
228     (primeResidue c p : ZMod p) = ((value c : ℚ) : ZMod p) := by
229   have hp : p.Prime := Fact.out
230   by_cases hden : normalDenominator c % p = 0
231   · have hd : (normalDenominator c : ZMod p) = 0 :=
232       (ZMod.natCast_eq_zero_iff _ _).mpr (Nat.dvd_of_mod_eq_zero hden)
233     simp only [primeResidue, hden, ↓reduceIte, Nat.cast_zero]
234     rw [Rat.cast_def, (normal_fields c).2, hd, div_zero]
235   · have hd : (normalDenominator c : ZMod p) ≠ 0 := by
236       intro hz
237       exact hden (Nat.mod_eq_zero_of_dvd ((ZMod.natCast_eq_zero_iff _ _).mp hz))
238     have hinv : (normalDenominator c : ZMod p) ^ (p - 2) =
239         (normalDenominator c : ZMod p)⁻¹ := by
240       symm
241       apply ZMod.inv_eq_of_mul_eq_one
242       have he : p - 2 + 1 = p - 1 := by have := hp.two_le; omega
243       rw [mul_comm, ← pow_succ, he]
244       exact ZMod.pow_card_sub_one_eq_one hd
245     rw [primeResidue, ite_eq_right hden, ZMod.natCast_mod, Nat.cast_mul,
246       numeratorResidue_cast c hp, Nat.cast_pow, hinv, Rat.cast_def, (normal_fields c).2,
247       div_eq_mul_inv]
248 
249 /-- The emitted representative is canonical in the finite field. -/
250 theorem primeResidue_lt (c : RationalCode) {p : ℕ} (hp : p.Prime) : primeResidue c p < p := by
251   unfold primeResidue
252   split
253   · exact hp.pos
254   · exact Nat.mod_lt _ hp.pos
255 
256 theorem primeResidue_eq_val (c : RationalCode) {p : ℕ} [Fact p.Prime] :
257     primeResidue c p = (((value c : ℚ) : ZMod p)).val := by
258   rw [← primeResidue_cast c, ZMod.val_natCast_of_lt (primeResidue_lt c (Fact.out : p.Prime))]
259 
260 end RationalCode
261 end RealizationBoundaries

Back to source controls