Gabriel Pichot avatar Gabriel Pichot committed 8dc79ce

Exponentiation modulaire ajoutée, corrections dans matrixFactory

Comments (0)

Files changed (6)

 # }}}
 
 DIRS =
-SOURCES = log.ml arithmetik.ml matrixFactory.ml quantum.ml main.ml
+SOURCES = log.ml matrixFactory.ml testMatrix.ml quantum.ml arithmetik.ml main.ml
 EXEC = shor
 
 
  *)
 
 open Log;;
+open Quantum;;
 
 (* primeTest {{{1
  * Détermine si un nombre est premier ou non (exact).
   | _ -> let r = pow x (n / 2) in
     if n mod 2 = 0 then r * r else r * r * x
 (* }}} *)
+(* Exponentiation modulaire {{{1 *)
+(* Connaissant x a n on calcule x ^ a mod n *)
+let expModulaire x a n =
+  let rec expMod_aux xi ai r = match ai with
+    | 0 -> r
+    | _ -> let r = if ai land 1 = 1 then r * xi mod n else r
+      in expMod_aux (xi * xi mod n) (ai lsr 1) r
+  in expMod_aux x a 1
+(* }}} *)
 
 (* Recherche de l'ordre de p dans Z/nZ {{{1 *)
 (* On créé plusieurs petites fonctions auxiliaires *)
   let l = getQ n in
   let q = pow 2 l in
   log "q = %i = 2 ^ %i\n" q l;
-  l
+  let reg1 = new register l
+  and reg2 = new register l in
+  (* On met le premier registre dans un état de superposition uniforme *)
+  reg1#setUniformSuperposition q;
+  (* On calcule x ^ i mod n dans le second... *)
+  for i = 0 to reg2#nbStates () - 1 do
+    reg2#setStateProbability i (complex_of_int(expModulaire p i n))
+  done;
+  (* ... mais on veut des probabilités ! *)
+  reg2#normalize ();
+
+  (* Maintenant on applique la transformée de Fourier *)
+
+(*  reg1#dump ();
+  reg2#dump (); *)
+   l
+  
+
+(* }}} *)
 
 (** Vars **)
 (* Open the logfile *)
-let logfile = stdout (* open_out "shor.log" *)
+let logfile =  open_out "shor.log" (*stdout * open_out "shor.log" *)
 (* Enabled debug mod *)
 let mod_debug = ref false
 
  *  -? ne pas obliger l'utilisateur à rentrer l'entier p. 
  *)
 open Log;;        (* Log, debug *)
-open Quantum;;  
 open Sys;;        (* Arguments *)
 open Printf;;     (* Printf *)
-open Arithmetik;; (* primeTest *)
-
+open Arithmetik;; (* primeTest, order ...*)
 
 (* Lire un entier depuis l'entrée *)
 let read_int () = Scanf.scanf " %i" (fun x -> x);;
 (* }}} *)
 
 
-
   | _ -> iterator ~i:(i+1) n f (f i acc)
 let rec viterator ?(i=0) n f = match i with
   | _ when i > n || i = n -> ()
-  | _ -> viterator ~i:(i+1) n f
+  | _ -> f i; viterator ~i:(i+1) n f
 
 (* Ev Module {{{1 *)
 module type K = sig
   end
   (* }}} *)
   (* Vector Class {{{2 *)
-  class vector ?data n = object(self)
-    inherit matrix ~dims:(n,1) ?data () as super
+  class vector ?data ?(rows=0) () = object(self)
+    inherit matrix ~dims:(rows,1) ?data () as super
     method row n = super#get n 1
     method rowset n = super#set n 1
-    method norm () = sqrt (iterator n begin fun i sum -> 
-     (Ev.norm2 (self#row i)) +. sum
+    method norm () = sqrt (iterator rows begin fun i sum -> 
+     (Ev.norm2 (self#row (i + 1))) +. sum
     end 0.)
     method normalize () = let norm = self#norm () in
-      viterator n (fun i -> self#rowset i (Ev.mul (self#row i) (1. /. norm)))
+      viterator rows begin fun i -> 
+        self#rowset (i + 1) (Ev.mul (self#row (i + 1)) (1. /. norm))
+      end
   end
   (* }}} *)
   (* Operators Module {{{2 *)
   Complex.re = f;
   Complex.im = 0.;
 }
+and complex_of_int i = { 
+  Complex.re = float_of_int i;
+  Complex.im = 0.;
+}
 (* }}} *)
 
 
 
 
 (* Register Class {{{1 *)
-exception State_Does_Not_Exist of string;;
+exception Quant_Bad_Access of string;;
 class register n = object
   val size = n
-  val state = new vector (1 lsl n) (* Equivalent à 2^n *)
+  val state = new vector ~rows:(1 lsl n) () (* Equivalent à 2^n *)
   method size () = size
-  (* normalize permet d'être sûr de bien avoir des probabilités *)
+  method nbStates () = 1 lsl n
   method normalize () = state#normalize ()
-  method setStateProbability s v = state#rowset s v
-  (* Non destructif *)
+  method setStateProbability s v = state#rowset (s+1) v
   method getStateProbability s = 
-    if s > state#rows () then raise (State_Does_Not_Exist "getStateProbability")
+    if s > state#rows () then raise (Quant_Bad_Access "getStateProbability")
     else state#row (s+1)
   method dump () =
     print "Le registre est dans l'état :\n";
     for i = 1 to state#rows () do
-      printf "État %i de probabilité %f.\n" (i-1) (Complex.norm(state#row i));
+      printf "État %i : %f.\n" (i-1) (Complex.norm(state#row i));
     done
+  (* Met les n premiers états dans un état de superposition uniforme *)
+  method setUniformSuperposition n =
+    if n > state#rows () then raise (Quant_Bad_Access "setUniformSuperposition")
+    else begin
+      let prob = complex_of_float (sqrt(1. /. (float_of_int n))) in
+      for i = 1 to n do
+        state#rowset i prob
+      done
+    end
   (* Destructif *)
   method measureState () =
     (* On procède ainsi, imaginons l'état suivant :
Tip: Filter by directory path e.g. /media app.js to search for public/media/app.js.
Tip: Use camelCasing e.g. ProjME to search for ProjectModifiedEvent.java.
Tip: Filter by extension type e.g. /repo .js to search for all .js files in the /repo directory.
Tip: Separate your search with spaces e.g. /ssh pom.xml to search for src/ssh/pom.xml.
Tip: Use ↑ and ↓ arrow keys to navigate and return to view the file.
Tip: You can also navigate files with Ctrl+j (next) and Ctrl+k (previous) and view the file with Ctrl+o.
Tip: You can also navigate files with Alt+j (next) and Alt+k (previous) and view the file with Alt+o.