Initial version of Fermat factorization

diff --git a/numbertheory.txt b/numbertheory.txt
index 479d6d4..336fe7c 100644 (file)
--- a/numbertheory.txt
+++ b/numbertheory.txt
@@ -22,7 +22,7 @@ proc integerroot(a:nat): nat
 {
    var x:nat = 0;
    var y:ℕ[(K+1)^2] = 1;
-   var z:nat = 1;
+   var z:ℕ[(K+1)^2] = 1;
 
    while y ≤ a do
       invariant x^2 ≤ a ∧ y = x^2+z ∧ z = 2*x+1;
@@ -187,6 +187,35 @@ theorem lemma11(a:nat,p:nat) ⇔
                      (legendre(ap,pp) = legendre2(ap,pp));
 
 ////////////////////////////////////////////////////////////////////
+// Algorithm: Fermat Factorization
+////////////////////////////////////////////////////////////////////
+proc fermatfactor(t:ℕ[K]): nat
+   ensures divides(result,t);
+{
+   var factor:ℕ[2*K+2] = 1;
+   var a:ℕ[(2*K+1)^2] = 2*integerroot(t) + 1;
+   var b:ℕ[(2*K+1)^2] = 1;
+   var r:ℤ[-(K^3),(K^3)] = integerroot(t)^2-t;
+
+   if t ≠ 0 then
+   {
+      while r ≠ 0 do
+      {
+         r = r+a;
+         a = a+2;
+
+         while r>0 do
+         {
+            r = r-b;
+            b = b+2;
+         }
+      }
+      factor = (a-b)/2;
+   }
+   return factor;
+}
+
+////////////////////////////////////////////////////////////////////
 // Algorithm: Calculate x^n by Right-To-Left Exponentation
 ////////////////////////////////////////////////////////////////////
 proc r2lexponentation(x:Base,n:Exp): Result
@@ -227,3 +256,5 @@ proc xpowern(x:Base,n:Exp): Result
    return res;
 }
 
+
+