(** * Factorisation de grands entiers à l'aide de courbes elliptiques * * @author Samuel Mimram *) (* * Copyright (C) 2002-2006 Samuel Mimram * * ECFacto is free software; you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by * the Free Software Foundation; either version 2 of the License, or * (at your option) any later version. * * ECFacto is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * GNU General Public License for more details. * * You should have received a copy of the GNU General Public License * along with Ocaml-mad; if not, write to the Free Software * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA *) open Big_int let y_init = 1 let bb_incr = unit_big_int and a_incr = big_int_of_int 17 and nb_essais = 1000 and x_init = unit_big_int and y_init = unit_big_int let a = ref unit_big_int and b = ref unit_big_int and p = ref zero_big_int exception Error_factor of big_int let rec bezout a b k = if lt_big_int a b then bezout b a k else ( let q, r = quomod_big_int a b in if eq_big_int r zero_big_int then k (zero_big_int, unit_big_int) else bezout b r (fun (u, v) -> k (v, sub_big_int u (mult_big_int q v))) ) let inv_big_int b n = let d, u, v = bezout b n ( fun (u, v) -> if gt_big_int b n then (add_big_int (mult_big_int b u) (mult_big_int n v), u, v) else (add_big_int (mult_big_int n u) (mult_big_int b v), v, u) ) in if eq_big_int d unit_big_int then u else raise (Error_factor d) type ec_point = O | ECP of big_int * big_int let ecp_print = function | O -> Printf.printf "Point : O\n" | ECP(x, y) -> Printf.printf "Point : %s, %s\n" (string_of_big_int x) (string_of_big_int y) let ecp_eq p1 p2 = match (p1, p2) with | (O, O) -> true | (O, _) | (_, O) -> false | (ECP(x1, y1), ECP(x2, y2)) -> eq_big_int x1 x2 && eq_big_int y1 y2 let ecp_add p1 p2 = match (p1, p2) with | (O, _) -> p2 | (_, O) -> p1 | (ECP(x1, y1), ECP(x2, y2)) -> if eq_big_int x1 x2 && eq_big_int y1 (mod_big_int (minus_big_int y2) !p) then O else let lambda = if ecp_eq p1 p2 then ( mod_big_int (mult_big_int (add_big_int (mult_int_big_int 3 (square_big_int x1)) !a) (inv_big_int (mult_int_big_int 2 y1) !p)) !p ) else ( mod_big_int (mult_big_int (sub_big_int y2 y1) (inv_big_int (sub_big_int x2 x1) !p)) !p ) in let x_ret = mod_big_int (sub_big_int (sub_big_int (square_big_int lambda) x1) x2) !p in ECP (x_ret, (mod_big_int (sub_big_int (mult_big_int lambda (sub_big_int x1 x_ret)) y1) !p)) let ecp_mult m p = let two = succ_big_int unit_big_int in let rec aux m pow res = if eq_big_int m zero_big_int then res else aux (div_big_int m two) (ecp_add pow pow) ( if eq_big_int (mod_big_int m two) unit_big_int then ecp_add res pow else res ) in aux m p O exception Result of big_int let calc_k b = div_big_int (mult_big_int b (pred_big_int b)) (gcd_big_int b (pred_big_int b)) let find_factor n bb' = let bb = ref bb' in let pt = ECP(x_init, y_init) and qt = ref O in let k = ref zero_big_int and tmp = ref zero_big_int in let a_essais = ref 0 and new_curve = ref false in if eq_big_int (mod_big_int n (succ_big_int unit_big_int)) zero_big_int then big_int_of_int 2 else if eq_big_int (mod_big_int n (big_int_of_int 3)) zero_big_int then big_int_of_int 3 else ( p := n; k := calc_k !bb; a := unit_big_int; b := unit_big_int; try while(0 = 0) do let ECP(x, y) = pt in b := sub_big_int (sub_big_int (square_big_int y) (power_big_int_positive_int x 3)) (mult_big_int !a x); tmp := gcd_big_int (add_big_int (mult_int_big_int 4 (power_big_int_positive_int !a 3)) (mult_int_big_int 27 (square_big_int !b))) n; while (not (eq_big_int !tmp unit_big_int)) do if eq_big_int !tmp n then ( a := add_big_int !a a_incr; tmp := gcd_big_int (add_big_int (mult_int_big_int 4 (power_big_int_positive_int !a 3)) (mult_int_big_int 27 (square_big_int !b))) n; ) else raise (Result !tmp) done; ( try qt := ecp_mult !k pt; with Error_factor f -> let r = gcd_big_int f n in if (lt_big_int r n) && (gt_big_int r unit_big_int) then raise (Result r) else a_essais := nb_essais; ); if (!a_essais < nb_essais) then ( a := add_big_int !a a_incr; incr a_essais; ) else ( a := unit_big_int; a_essais := 0; bb := add_big_int !bb bb_incr; while (eq_big_int !k (calc_k !bb)) do bb := add_big_int !bb bb_incr; done; k := calc_k !bb; ) done; pred_big_int zero_big_int; with | Result r -> r | _ -> pred_big_int zero_big_int ) let _ = let n = big_int_of_string "109849677793909" and bb = big_int_of_string "10487395" in Printf.printf "%s\n" (string_of_big_int (find_factor n bb))