From f07dc33c73f5fd2b6b54b055d6daf9fc64f3b22f Mon Sep 17 00:00:00 2001
From: Anthony Scemama <scemama@irsamc.ups-tlse.fr>
Date: Mon, 18 Mar 2019 19:17:31 +0100
Subject: [PATCH] PT2 is not working yet

---
 run_cas.ml | 82 ++++++++++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 82 insertions(+)
 create mode 100644 run_cas.ml

diff --git a/run_cas.ml b/run_cas.ml
new file mode 100644
index 0000000..05b7f90
--- /dev/null
+++ b/run_cas.ml
@@ -0,0 +1,82 @@
+let () =
+  let open Command_line in
+  begin
+    set_header_doc (Sys.argv.(0) ^ " - QCaml command");
+    set_description_doc "Runs a CAS-CI calculation";
+    set_specs
+      [ { short='b' ; long="basis" ; opt=Mandatory;
+          arg=With_arg "<string>";
+          doc="Name of the file containing the basis set"; } ;
+
+        { short='x' ; long="xyz" ; opt=Mandatory;
+          arg=With_arg "<string>";
+          doc="Name of the file containing the nuclear coordinates in xyz format"; } ;
+
+        { short='m' ; long="multiplicity" ; opt=Optional;
+          arg=With_arg "<int>";
+          doc="Spin multiplicity (2S+1). Default is singlet"; } ;
+
+        { short='c' ; long="charge" ; opt=Optional;
+          arg=With_arg "<int>";
+          doc="Total charge of the molecule. Default is 0"; } ;
+
+        { short='n' ; long="n_elec" ; opt=Mandatory;
+          arg=With_arg "<int>";
+          doc="Number of electrons in the active space."; } ;
+
+        { short='o' ; long="n_orb" ; opt=Mandatory;
+          arg=With_arg "<int>";
+          doc="Number of MOs in the active space."; } ;
+      ]
+  end;
+
+  (* Handle options *)
+  let basis_file  = Util.of_some @@ Command_line.get "basis" in
+
+  let nuclei_file = Util.of_some @@ Command_line.get "xyz" in
+
+  let n = Util.of_some @@ Command_line.get "n_elec" |> int_of_string in
+  let m = Util.of_some @@ Command_line.get "n_orb" |> int_of_string in
+
+  let charge = 
+    match Command_line.get "charge" with
+    | Some x -> int_of_string x
+    | None   -> 0
+  in
+
+  let multiplicity =
+    match Command_line.get "multiplicity" with
+    | Some x -> int_of_string x
+    | None -> 1
+  in
+
+  let ppf = 
+    if Parallel.master then Format.std_formatter
+    else Printing.ppf_dev_null
+  in
+
+  let s =
+    Simulation.of_filenames ~charge ~multiplicity ~nuclei:nuclei_file basis_file
+  in
+
+  let hf = HartreeFock.make s in
+  Format.fprintf ppf "@[%a@]@." HartreeFock.pp_hf hf;
+
+  let mos =
+    MOBasis.of_hartree_fock hf
+  in
+
+  let space =
+    DeterminantSpace.cas_of_mo_basis mos n m
+  in
+  let ci = CI.make space in
+  Format.fprintf ppf "CAS-CI energy : %20.16f@." ((CI.eigenvalues ci).{1} +. Simulation.nuclear_repulsion s);
+
+  let pt2 = CI.pt2 ci in
+  Format.fprintf ppf "CAS-EN2 energy : %20.16f@." ((CI.eigenvalues ci).{1} +. Simulation.nuclear_repulsion s +. pt2)
+  (*
+  let s2 = Util.xt_o_x ~o:(CI.s2_matrix ci) ~x:(CI.eigenvectors ci) in
+  Util.list_range 1 (DeterminantSpace.size space)
+  |> List.iter (fun i -> Format.printf "@[%f@]@;" s2.{i,i});
+  *)
+