#open "sys";;
#open "graphics";;

let dir = "." ;;  (* répertoire où est installé pi.txt *)
chdir dir;;

(* ancienne version plus soft à comprendre *

let eratos n =
        let t = make_vect ((n-1)/2) true
        and i = ref 0 and r = int_of_float (sqrt (float_of_int n)) in
        while !i*2+3 <= r do
                while not t.(!i) do incr i done;
                for j = !i to (n/(!i*2+3) -1) /2 -1 do t.((!i*2+3)*(2*j+3) /2 -1) <- false done;
                incr i;
        done;
        let l = ref [] in
        for j = (n-1)/2 -1 downto 0 do if t.(j) then l:= (2*j+3) :: !l done;
        t,(2::!l);
;;

*)

let eratos n =
        let f k = if k mod 2 = 0 then 6*(k/2+1)-1 else 6*(k/2+1)+1
        and g n = n /12 *4 + n mod 12 /3 -1 in
        let fog n = let m = f (g n) in if m > n then f ((g n)-1) else m in
        let nn = fog n in
        let t = make_vect ((g nn)+1) true
        and i = ref 0 and r = int_of_float (sqrt (float_of_int nn)) in
        while (f !i) <= r do
                while not t.(!i) do incr i done;
                for j = !i to g (fog (n/(f !i))) do t.(g ((f !i)*(f j))) <- false done;
                incr i;
        done;
        let l = ref [] in
        for j = g nn downto 0 do if t.(j) then l:= (f j) :: !l done;
        t,(2::3::!l);
;;

let t = ref [||] and l = ref [] and nm = ref 0 and nl = ref 0 ;;

let init0 y n =
        if y || !nm < n then (
                let a,b = eratos n in
                t:= a ; l:= b ; nm:= n ; nl:= list_length b;
        );
;;

let init = init0 true ;;

init 32771;;


let factors0 b n =
        let rec pval p i = fun
                | n when n mod p = 0 -> pval p (i+1) (n/p)
                | n -> i,n
        in
        let rec aux = fun
                | 1 _ -> []
                | n (a::_) when a*a > n || a*a < 0 -> [n,1]
                | n (a::q) -> (let i,nn = pval a 0 n in
                               if i = 0 then aux n q else (if b then raise Exit else (a,i)::(aux nn q)))
                | n _ -> failwith ("décomposition inconnue de " ^ (string_of_int n))
        in
        aux n !l;
;;

let factors = factors0 false ;;

let isprime n =
        if n < 0 then failwith "borne max dépassée";
        try
                if !nm < n then factors0 true n = [n,1] else (
                if n = 1 then false else (if n = 2 || n = 3 then true else (
                if n mod 2 = 0 || n mod 3 = 0 then false
                else !t.(n /12 *4 + n mod 12 /3 -1))));
        with Exit -> false
;;

let pretty_string n =
        let aux1 p d =
                (string_of_int p) ^ match d with
                        | 1 -> ""
                        | 2 -> "²"   (* "\253" pour le bin mode *)
                        | 3 -> "³"   (* "\252" pour le bin mode *)
                        | d -> "^" ^ string_of_int d
        in
        let rec aux = fun
                | [] -> ""
                | [p,d] -> aux1 p d
                | ((p,d)::q) -> (aux1 p d) ^ "x" ^ (aux q)
        in
        if n > 1 then (
                let f = factors n in
                let nn = string_of_int n in
                if f = [n,1] then nn ^ " est premier."
                else nn ^ " = " ^ aux f;
        ) else "";
;;

let pretty_dec n = print_endline (pretty_string n) ;;


(* lecture d'une base de données n,pi(n) *)

let lpi =
        let lpi0 = [12582893,823749] in
        try
                let c = open_in "pi.txt" and
                lpi = ref [] in
                try while true do
                        let s = input_line c in
                        let n = ref 13 in
                        while s.[!n] = ` ` do decr n done;
                        lpi:= ((int_of_string (sub_string s 0 (!n+1))),
                               (int_of_string (sub_string s 14 (string_length s - 14)))) :: !lpi;
                done;
                [];
                with End_of_file -> close_in c ; sort__merge (fun x y -> (fst x) > (fst y)) lpi0 !lpi
        with Sys_error _ -> lpi0
;;

let rec extraction ineq = fun
        | [] -> failwith "mauvaise liste d'initialisation"
        | (a::_) when ineq a -> a
        | (_::q) -> extraction ineq q
;;


(* k_ième nombre premier *)

let ithprime k =
        if !nl < k && k < 823750 then (
                if interactive then init0 false 12582912
                else (
                        let kk = float_of_int k in
                        init (min (int_of_float (kk*.(log kk)*.1.14260611623)) 12582912);
                );
        );
        let rec aux = fun
                | [] _ -> failwith "borne max dépassée"
                | (a::_) 1 -> a
                | (_::q) n -> aux q (n-1)
        in
        if k < 823750 then aux !l k else (
                let e = extraction (fun x -> snd x < k) lpi in
                let i = ref (snd e) and n = ref (fst e) in
                n:= !n - 1 + !n mod 2;
                while !i < k do
                        n:= !n + 2;
                        if isprime !n then incr i;
                done;
                !n;
        );
;;

(* la fonction inverse de ithprime, nb de nb premiers <= n *)

let pi n =
        let b = !nm < n && n < 12582893 in
        if b && not interactive then (
                init n;
                !nl;
        ) else (
                if b then init0 false 12582912;
                let rec aux n i = fun
                        | [] -> i
                        | (a::_) when a > n -> i
                        | (_::q) -> aux n (i+1) q
                in
                if n < 12582893 then aux n 0 !l else (
                        let e = extraction (fun x -> fst x <= n) lpi in
                        let i = ref (snd e) and a = ref (fst e) in
                        a:= !a + 1 + !a mod 2;
                        while !a <= n && !a > 0 do
                                if isprime !a then incr i;
                                a:= !a + 2;
                        done;
                        !i;
                );
        );
;;


let nextprime n =
        let i = ref (n+1 + n mod 2) in
        while not isprime !i do i:= !i+2 done;
        !i;
;;

let prevprime n =
        let i = ref (n-1 - n mod 2) in
        while not isprime !i do i:= !i-2 done;
        !i;
;;


let etalonne_pi s n pin step max =
        let i = ref pin and a = ref (n-1 + n mod 2) in
        if !a < 2 then i:= 1;
        while !a < max do
                a:= !a + 2;
                let b = isprime !a and aa = !a mod step in
                if aa <= 1 then (
                        let aaa = !a - aa in
                        if aaa > n && aaa <= max then (
                                let aaa = string_of_int aaa in
                                if b && aa = 0 then incr i;
                                print_endline (aaa ^ (make_string (14-(string_length aaa)) ` `) ^ (string_of_int !i));
                                if b && aa = 0 then decr i;
                                flush std_out;
                                if s then sound 200 400;
                        );
                );
                if b then incr i;
        done;
        if s then (sound 200 200 ; sound 400 200 ; sound 800 200;);
;;


let main args =
        let err = "prime [is|ith|pi|next|prev] <n> | prime etal [-s] <n> <pi(n)> <step> <max>" in
        try
                (try
                        if args.(1) = "etal" && args.(2) = "-s" then (
                                etalonne_pi true (int_of_string args.(3)) (int_of_string args.(4))
                                                 (int_of_string args.(5)) (int_of_string args.(6));
                                                 true;
                        ) else (
                                let nn = args.(2) in let n = int_of_string nn in
                                match args.(1) with
                                        | "is" -> if isprime n then print_endline (nn ^ " est premier.")
                                                  else print_endline (nn ^ " est composé.");
                                                  true
                                        | "ith" -> print_endline ("le " ^ nn ^ "e nombre premier est " ^ 
                                                        (string_of_int (ithprime n)));
                                                   true
                                        | "pi" -> print_endline ("pi(" ^ nn ^ ") = " ^ (string_of_int (pi n)));
                                                  true
                                        | "next" -> print_endline (" -> " ^ (string_of_int (nextprime n)));
                                                    true
                                        | "prev" -> print_endline (" -> " ^ (string_of_int (prevprime n)));
                                                    true
                                        | "etal" -> etalonne_pi false n (int_of_string args.(3)) (int_of_string args.(4))
                                                                        (int_of_string args.(5));
                                                    true
                                        | _ -> print_endline err ; false
                        );
                with Invalid_argument "vect_item" -> pretty_dec (int_of_string args.(1)) ; true)
        with | Invalid_argument "vect_item" -> print_endline err ; false
             | Failure "int_of_string" -> print_endline err ; false
             | Failure s -> print_endline s ; false
;;

if not interactive then (
        let b = main command_line in
        flush std_out;
        if b then exit 0 else exit 1;
);;