From beb5e5548141e3638fc9bb5a3a5fe274623140e3 Mon Sep 17 00:00:00 2001
From: Anthony Scemama <scemama@irsamc.ups-tlse.fr>
Date: Wed, 29 Oct 2014 22:13:03 +0100
Subject: [PATCH] Introduces Rst_string type

---
 ocaml/Atom.ml               |  2 +-
 ocaml/Input_ao_basis.ml     |  7 +++--
 ocaml/Input_bi_integrals.ml | 27 ++++++++++++----
 ocaml/Input_cis.ml          | 22 +++++++++++++
 ocaml/Input_cisd_sc2.ml     | 14 +++++++++
 ocaml/Input_determinants.ml |  7 +++--
 ocaml/Input_electrons.ml    |  6 ++--
 ocaml/Input_full_ci.ml      | 14 +++++++++
 ocaml/Input_hartree_fock.ml | 13 ++++++++
 ocaml/Input_mo_basis.ml     | 11 ++++---
 ocaml/Input_nuclei.ml       | 62 ++++++++++++++++++++++++++++++++++---
 ocaml/qp_edit.ml            | 35 +++++++++++----------
 ocaml/qptypes_generator.ml  |  3 ++
 ocaml/test_input.ml         | 35 ++++++++++++++++-----
 14 files changed, 209 insertions(+), 49 deletions(-)

diff --git a/ocaml/Atom.ml b/ocaml/Atom.ml
index 3b34951c..2aafe644 100644
--- a/ocaml/Atom.ml
+++ b/ocaml/Atom.ml
@@ -33,6 +33,6 @@ let to_string u a =
   [ Element.to_string a.element ;
     Charge.to_string  a.charge ;
     Point3d.to_string u a.coord ]
-  |> String.concat ?sep:(Some "   ")
+  |> String.concat ~sep:"   "
 ;;
 
diff --git a/ocaml/Input_ao_basis.ml b/ocaml/Input_ao_basis.ml
index fc4012bd..a094abd1 100644
--- a/ocaml/Input_ao_basis.ml
+++ b/ocaml/Input_ao_basis.ml
@@ -17,7 +17,7 @@ module Ao_basis : sig
   val read : unit -> t
   val to_string : t -> string
   val to_md5 : t -> MD5.t
-  val debug : t -> string
+  val to_rst : t -> Rst_string.t
 end = struct
   type t = 
     { ao_basis        : string ;
@@ -140,7 +140,7 @@ end = struct
     |> Long_basis.to_basis
   ;;
       
-  let to_string b =
+  let to_rst b =
     let short_basis = to_basis b in
     Printf.sprintf "
 Name of the AO basis ::
@@ -157,6 +157,7 @@ Basis set ::
      |> List.map ~f:(fun x-> "  "^x)
      |> String.concat ~sep:"\n"
     )
+  |> Rst_string.of_string
   ;;
 
   let to_md5 b =
@@ -164,7 +165,7 @@ Basis set ::
     Basis.to_md5 short_basis
   ;;
     
-  let debug b =
+  let to_string b =
     Printf.sprintf "
 ao_basis        = %s
 ao_num          = %s
diff --git a/ocaml/Input_bi_integrals.ml b/ocaml/Input_bi_integrals.ml
index d0a16288..b06d1d28 100644
--- a/ocaml/Input_bi_integrals.ml
+++ b/ocaml/Input_bi_integrals.ml
@@ -15,7 +15,7 @@ module Bielec_integrals : sig
   ;;
   val read : unit -> t
   val to_string : t -> string
-  val of_string : string -> t
+  val to_rst : t -> Rst_string.t
 end = struct
   type t = 
     { read_ao_integrals  : bool;
@@ -116,6 +116,25 @@ end = struct
 
   let to_string b =
     Printf.sprintf "
+read_ao_integrals = %s
+read_mo_integrals = %s
+write_ao_integrals = %s
+write_mo_integrals = %s
+threshold_ao = %s
+threshold_mo = %s
+direct = %s
+"
+        (Bool.to_string b.read_ao_integrals)
+        (Bool.to_string b.read_mo_integrals)
+        (Bool.to_string b.write_ao_integrals)
+        (Bool.to_string b.write_mo_integrals)
+        (Threshold.to_string b.threshold_ao)
+        (Threshold.to_string b.threshold_mo)
+        (Bool.to_string b.direct)
+  ;;
+
+  let to_rst b =
+    Printf.sprintf "
 Read AO/MO integrals from disk ::
 
   read_ao_integrals = %s
@@ -143,11 +162,7 @@ Direct calculation of integrals ::
         (Threshold.to_string b.threshold_ao)
         (Threshold.to_string b.threshold_mo)
         (Bool.to_string b.direct)
-  ;;
-
-  let of_string s =
-    input_to_sexp s
-    |> t_of_sexp
+  |> Rst_string.of_string
   ;;
 
 end
diff --git a/ocaml/Input_cis.ml b/ocaml/Input_cis.ml
index 98af9c43..79cf03e8 100644
--- a/ocaml/Input_cis.ml
+++ b/ocaml/Input_cis.ml
@@ -14,6 +14,7 @@ module Cis_dressed : sig
   ;;
   val read : unit -> t
   val to_string : t -> string
+  val to_rst : t -> Rst_string.t
 end = struct
   type t = 
     { n_state_cis        : States_number.t;
@@ -95,6 +96,24 @@ end = struct
 
   let to_string b =
     Printf.sprintf "
+n_state_cis = %s
+n_core_cis = %s
+n_act_cis  = %s
+mp2_dressing       = %s
+standard_doubles   = %s
+en_2_2             = %s
+"
+        (States_number.to_string b.n_state_cis)
+        (Positive_int.to_string b.n_core_cis)
+        (Positive_int.to_string b.n_act_cis)
+        (Bool.to_string b.mp2_dressing)
+        (Bool.to_string b.standard_doubles)
+        (Bool.to_string b.en_2_2)
+  ;;
+
+
+  let to_rst b =
+    Printf.sprintf "
 Number of states ::
 
   n_state_cis = %s
@@ -123,6 +142,9 @@ Epstein-Nesbet 2x2 diagonalization ::
         (Bool.to_string b.mp2_dressing)
         (Bool.to_string b.standard_doubles)
         (Bool.to_string b.en_2_2)
+
+   |> Rst_string.of_string
+  ;;
 end
 
 
diff --git a/ocaml/Input_cisd_sc2.ml b/ocaml/Input_cisd_sc2.ml
index f02a46a5..b87dcc91 100644
--- a/ocaml/Input_cisd_sc2.ml
+++ b/ocaml/Input_cisd_sc2.ml
@@ -11,6 +11,7 @@ module Cisd_sc2 : sig
   ;;
   val read : unit -> t
   val to_string : t -> string
+  val to_rst : t -> Rst_string.t
 end = struct
   type t = 
     { n_det_max_cisd_sc2 : Det_number_max.t;
@@ -62,6 +63,17 @@ end = struct
 
   let to_string b =
     Printf.sprintf "
+n_det_max_cisd_sc2 = %s
+pt2_max = %s
+do_pt2_end = %s
+"
+        (Det_number_max.to_string b.n_det_max_cisd_sc2)
+        (PT2_energy.to_string b.pt2_max)
+        (Bool.to_string b.do_pt2_end)
+  ;;
+
+  let to_rst b =
+    Printf.sprintf "
 Stop when the `n_det` > `n_det_max_cisd_sc2` ::
 
   n_det_max_cisd_sc2 = %s
@@ -78,6 +90,8 @@ Compute E(PT2) at the end ::
         (Det_number_max.to_string b.n_det_max_cisd_sc2)
         (PT2_energy.to_string b.pt2_max)
         (Bool.to_string b.do_pt2_end)
+   |> Rst_string.of_string
+
 end
 
 
diff --git a/ocaml/Input_determinants.ml b/ocaml/Input_determinants.ml
index 4b958c94..ec4b380c 100644
--- a/ocaml/Input_determinants.ml
+++ b/ocaml/Input_determinants.ml
@@ -22,7 +22,7 @@ module Determinants : sig
   ;;
   val read : unit -> t
   val to_string : t -> string
-  val debug : t -> string
+  val to_rst : t -> Rst_string.t
 end = struct
   type t = 
     { n_int                  : N_int_number.t;
@@ -245,7 +245,7 @@ end = struct
     }
   ;;
 
-  let to_string b =
+  let to_rst b =
     let mo_tot_num = Ezfio.get_mo_basis_mo_tot_num () 
     |> MO_number.of_int in
     let det_text = 
@@ -309,9 +309,10 @@ Determinants ::
      (b.n_det_max_jacobi |> Strictly_positive_int.to_string)
      (b.n_det         |> Det_number.to_string)
      det_text
+     |> Rst_string.of_string
   ;;
 
-  let debug b =
+  let to_string b =
     let mo_tot_num = Ezfio.get_mo_basis_mo_tot_num () 
     |> MO_number.of_int in
     Printf.sprintf "
diff --git a/ocaml/Input_electrons.ml b/ocaml/Input_electrons.ml
index 601a066c..c95830e1 100644
--- a/ocaml/Input_electrons.ml
+++ b/ocaml/Input_electrons.ml
@@ -11,6 +11,7 @@ module Electrons : sig
   ;;
   val read : unit -> t
   val to_string : t -> string
+  val to_rst : t -> Rst_string.t
 end = struct
   type t = 
     { elec_alpha_num     : Elec_alpha_number.t;
@@ -46,7 +47,7 @@ end = struct
     }
   ;;
 
-  let to_string b =
+  let to_rst b =
     Printf.sprintf "
 Spin multiplicity is %s.
 
@@ -60,9 +61,10 @@ Number of alpha and beta electrons ::
          |> Multiplicity.to_string)
         (Elec_alpha_number.to_string b.elec_alpha_num)
         (Elec_beta_number.to_string b.elec_beta_num)
+    |> Rst_string.of_string
   ;;
 
-  let debug b =
+  let to_string b =
     Printf.sprintf "elec_alpha_num     = %s
 elec_beta_num      = %s
 elec_num           = %s
diff --git a/ocaml/Input_full_ci.ml b/ocaml/Input_full_ci.ml
index 1a8b76ff..abe4e959 100644
--- a/ocaml/Input_full_ci.ml
+++ b/ocaml/Input_full_ci.ml
@@ -11,6 +11,7 @@ module Full_ci : sig
   ;;
   val read : unit -> t
   val to_string : t -> string
+  val to_rst : t -> Rst_string.t
 end = struct
   type t = 
     { n_det_max_fci      : Det_number_max.t;
@@ -60,6 +61,17 @@ end = struct
 
   let to_string b =
     Printf.sprintf "
+n_det_max_fci = %s
+pt2_max = %s
+do_pt2_end = %s
+"
+        (Det_number_max.to_string b.n_det_max_fci)
+        (PT2_energy.to_string b.pt2_max)
+        (Bool.to_string b.do_pt2_end)
+    ;;
+
+  let to_rst b =
+    Printf.sprintf "
 Stop when the `n_det` > `n_det_max_fci` ::
 
   n_det_max_fci = %s
@@ -76,6 +88,8 @@ Compute E(PT2) at the end ::
         (Det_number_max.to_string b.n_det_max_fci)
         (PT2_energy.to_string b.pt2_max)
         (Bool.to_string b.do_pt2_end)
+    |> Rst_string.of_string
+    ;;
 end
 
 
diff --git a/ocaml/Input_hartree_fock.ml b/ocaml/Input_hartree_fock.ml
index cc56e7b9..2444f0a0 100644
--- a/ocaml/Input_hartree_fock.ml
+++ b/ocaml/Input_hartree_fock.ml
@@ -10,6 +10,7 @@ module Hartree_fock : sig
   ;;
   val read : unit -> t
   val to_string : t -> string
+  val to_rst : t -> Rst_string.t
 end = struct
   type t = 
     { n_it_scf_max       : Strictly_positive_int.t;
@@ -47,6 +48,15 @@ end = struct
 
   let to_string b =
     Printf.sprintf "
+n_it_scf_max = %s
+thresh_scf = %s
+"
+      (Strictly_positive_int.to_string b.n_it_scf_max)
+      (Threshold.to_string b.thresh_scf)
+  ;;
+
+  let to_rst b =
+    Printf.sprintf "
 Max number of SCF iterations ::
 
   n_it_scf_max = %s
@@ -58,6 +68,9 @@ SCF convergence criterion (on energy) ::
 "
         (Strictly_positive_int.to_string b.n_it_scf_max)
         (Threshold.to_string b.thresh_scf)
+    |> Rst_string.of_string
+  ;;
+
 end
 
 
diff --git a/ocaml/Input_mo_basis.ml b/ocaml/Input_mo_basis.ml
index 94024aba..1d1a9aee 100644
--- a/ocaml/Input_mo_basis.ml
+++ b/ocaml/Input_mo_basis.ml
@@ -12,7 +12,7 @@ module Mo_basis : sig
   ;;
   val read : unit -> t
   val to_string : t -> string
-  val debug : t -> string
+  val to_rst : t -> Rst_string.t
 end = struct
   type t = 
     { mo_tot_num      : MO_number.t ;
@@ -132,15 +132,15 @@ end = struct
       | _ -> assert false
     in
     let rec create_list accu i = 
-      if (i+5 < mo_tot_num) then
-        create_list ( (print_five i (i+4) |> String.concat ~sep:"\n")::accu ) (i+5)
+      if (i+4 < mo_tot_num) then
+        create_list ( (print_five i (i+3) |> String.concat ~sep:"\n")::accu ) (i+4)
       else
         (print_five i (mo_tot_num-1) |> String.concat ~sep:"\n")::accu |> List.rev
     in
     create_list [] 0 |> String.concat ~sep:"\n\n"
   ;;
 
-  let to_string b =
+  let to_rst b =
     Printf.sprintf "
 Label of the molecular orbitals ::
 
@@ -157,10 +157,11 @@ MO coefficients ::
     (MO_label.to_string b.mo_label)
     (MO_number.to_string b.mo_tot_num)
     (mo_coef_to_string b.mo_coef)
+    |> Rst_string.of_string
 
   ;;
 
-  let debug b =
+  let to_string b =
     Printf.sprintf "
 mo_label        = %s
 mo_tot_num      = \"%s\"
diff --git a/ocaml/Input_nuclei.ml b/ocaml/Input_nuclei.ml
index 9a285669..969a78b4 100644
--- a/ocaml/Input_nuclei.ml
+++ b/ocaml/Input_nuclei.ml
@@ -11,8 +11,9 @@ module Nuclei : sig
     } with sexp
   ;;
   val read : unit -> t
-  val debug : t -> string
   val to_string : t -> string
+  val to_rst : t -> Rst_string.t
+  val of_rst : Rst_string.t -> t
 end = struct
   type t = 
     { nucl_num        : Nucl_number.t ;
@@ -58,6 +59,58 @@ end = struct
     result
   ;;
 
+  let of_rst s =
+    let l = Rst_string.to_string s 
+    |> String.split ~on:'\n'
+    in
+    (* Find lines containing the xyz data *)
+    let rec extract_begin = function
+    | [] -> raise Not_found
+    | line::tail ->
+      let line = String.strip line in
+      if (String.length line > 3) &&
+        (String.sub line ~pos:((String.length line)-2)
+            ~len:2 = "::") then
+           tail
+      else
+         extract_begin tail
+    in
+    (* Read the xyz data *)
+    let rec read_line = function
+    | [] -> []
+    | line::tail -> 
+      if (String.strip line = "") then []
+      else 
+        (read_line tail)
+    in
+    (* Create a list of Atom.t *)
+    let atom_list = 
+      match (extract_begin l) with 
+      | _ :: nucl_num :: title :: lines ->
+        begin
+          let nucl_num = nucl_num
+          |> String.strip
+          |> Int.of_string
+          |> Nucl_number.of_int
+          and lines = Array.of_list lines
+          in
+          List.init (Nucl_number.to_int nucl_num) ~f:(fun i ->
+            Atom.of_string Units.Angstrom lines.(i))
+        end
+      | _ -> failwith "Error in xyz format"
+    in
+    (* Create the Nuclei.t data structure *)
+    { nucl_num = List.length atom_list
+        |> Nucl_number.of_int;
+      nucl_label = List.map atom_list ~f:(fun x ->
+        x.Atom.element) |> Array.of_list ;
+      nucl_charge = List.map atom_list ~f:(fun x ->
+        x.Atom.charge ) |> Array.of_list ;
+      nucl_coord = List.map atom_list ~f:(fun x ->
+        x.Atom.coord ) |> Array.of_list ;
+    }
+  ;;
+
   let read () =
     { nucl_num        = read_nucl_num ();
       nucl_label      = read_nucl_label () ;
@@ -66,7 +119,7 @@ end = struct
     }
   ;;
 
-  let debug b =
+  let to_string b =
     Printf.sprintf "
 nucl_num         = %s
 nucl_label       = %s
@@ -82,7 +135,7 @@ nucl_coord       = %s
       ~f:(Point3d.to_string Units.Bohr) |> String.concat ~sep:"\n" )
   ;;
 
-   let to_string b = 
+   let to_rst b = 
      let nucl_num = Nucl_number.to_int b.nucl_num in
      let text = 
        ( Printf.sprintf "  %d\n  "
@@ -101,8 +154,9 @@ Nuclear coordinates in xyz format (Angstroms) ::
 %s
 
 " text
-     
+     |> Rst_string.of_string
 ;;
+     
 end
 
 
diff --git a/ocaml/qp_edit.ml b/ocaml/qp_edit.ml
index 410f5fa1..622be6d0 100644
--- a/ocaml/qp_edit.ml
+++ b/ocaml/qp_edit.ml
@@ -43,27 +43,28 @@ let make_header kw =
 ;;
 
 let get s = 
-  let header = (make_header s)
-  in header^(match s with
-  | Bielec_integrals -> 
-    Input.Bielec_integrals.(to_string (read ()))
-  | Ao_basis ->
-    Input.Ao_basis.(to_string (read ()))
-  | Cisd_sc2 ->
-    Input.Cisd_sc2.(to_string (read ()))
-  | Determinants ->
-    Input.Determinants.(to_string (read ()))
-  | Electrons ->
-    Input.Electrons.(to_string (read ()))
+  let header = (make_header s) 
+  and rst = match s with
   | Full_ci ->
-    Input.Full_ci.(to_string (read ()))
+    Input.Full_ci.(to_rst (read ()))
   | Hartree_fock ->
-    Input.Hartree_fock.(to_string (read ()))
+    Input.Hartree_fock.(to_rst (read ()))
   | Mo_basis ->
-    Input.Mo_basis.(to_string (read ()))
+    Input.Mo_basis.(to_rst (read ()))
+  | Electrons ->
+    Input.Electrons.(to_rst (read ()))
+  | Determinants ->
+    Input.Determinants.(to_rst (read ()))
+  | Cisd_sc2 ->
+    Input.Cisd_sc2.(to_rst (read ()))
   | Nuclei ->
-    Input.Nuclei.(to_string (read ()))
-  )
+    Input.Nuclei.(to_rst (read ()))
+  | Ao_basis ->
+    Input.Ao_basis.(to_rst (read ()))
+  | Bielec_integrals -> 
+    Input.Bielec_integrals.(to_rst (read ()))
+ 
+  in header^(Rst_string.to_string rst)
 ;;
 
 (*
diff --git a/ocaml/qptypes_generator.ml b/ocaml/qptypes_generator.ml
index ced370cb..2633bdcc 100644
--- a/ocaml/qptypes_generator.ml
+++ b/ocaml/qptypes_generator.ml
@@ -131,6 +131,9 @@ let input_data = "
 
 * MD5 : string
   assert ((String.length x) = 32);
+
+* Rst_string : string
+
 "
 ;;
 
diff --git a/ocaml/test_input.ml b/ocaml/test_input.ml
index ee5bba60..7f200a92 100644
--- a/ocaml/test_input.ml
+++ b/ocaml/test_input.ml
@@ -1,9 +1,11 @@
+open Qptypes;;
+
 let test_ao () =
   Ezfio.set_file "F2.ezfio" ;
   let b = Input.Ao_basis.read ()
   in
-  print_endline (Input.Ao_basis.debug b);
   print_endline (Input.Ao_basis.to_string b);
+  print_endline (Input.Ao_basis.to_rst b |> Rst_string.to_string);
 ;;
 
 let test_bielec_intergals () =
@@ -13,12 +15,6 @@ let test_bielec_intergals () =
   let output = Input.Bielec_integrals.to_string b
   in
   print_endline output;
-  let b2 = Input.Bielec_integrals.of_string output in
-  if (b = b2) then
-   print_endline "OK"
-  else
-   print_endline "b <> b2"
-  ;
 ;;
 
 let test_bitmasks () =
@@ -85,6 +81,29 @@ let test_nucl () =
   print_endline (Input.Nuclei.to_string b);
 ;;
 
+let test_nucl_read () =
+  let rst_input = Rst_string.of_string "
+Molecule
+========
+
+
+Nuclear coordinates in xyz format (Angstroms) ::
+
+  2
+  
+  F    9 0.00000000  0.00000000 -0.70000000
+  F    9 0.00000000  0.00000000 0.70000000
+
+
+
+Electrons
+=========
+" in
+  let b = Input.Nuclei.of_rst rst_input
+  in
+  print_endline (Input.Nuclei.to_string b);
+;;
+
 (*
 test_ao ();;
 test_bielec_intergals ();;
@@ -96,5 +115,5 @@ test_hf ();;
 test_mo ();;
 test_nucl ();
 *)
-test_nucl();;
+test_nucl_read();;