open Lacaml.D let () = let open Command_line in begin set_header_doc (Sys.argv.(0) ^ " - QCaml command"); set_description_doc "Runs a F12-Full CI calculation"; set_specs [ { short='b' ; long="basis" ; opt=Mandatory; arg=With_arg ""; doc="Name of the file containing the basis set"; } ; { short='a' ; long="aux_basis" ; opt=Mandatory; arg=With_arg ""; doc="Name of the file containing the auxiliary basis set"; } ; { short='x' ; long="xyz" ; opt=Mandatory; arg=With_arg ""; doc="Name of the file containing the nuclear coordinates in xyz format"; } ; { short='m' ; long="multiplicity" ; opt=Optional; arg=With_arg ""; doc="Spin multiplicity (2S+1). Default is singlet"; } ; { short='c' ; long="charge" ; opt=Optional; arg=With_arg ""; doc="Total charge of the molecule. Default is 0"; } ; (* { short='e' ; long="expo" ; opt=Optional; arg=With_arg ""; doc="Exponent of the Gaussian geminal"; } ; *) { short='s' ; long="state" ; opt=Optional; arg=With_arg ""; doc="Requested state. Default is 1."; } ; ] end; (* Handle options *) let basis_file = Util.of_some @@ Command_line.get "basis" in let aux_basis_filename = Util.of_some @@ Command_line.get "aux_basis" in let nuclei_file = Util.of_some @@ Command_line.get "xyz" in let charge = match Command_line.get "charge" with | Some x -> ( if x.[0] = 'm' then ~- (int_of_string (String.sub x 1 (String.length x - 1))) else int_of_string x ) | None -> 0 in (* let expo = match Command_line.get "expo" with | Some x -> float_of_string x | None -> 1.0 in *) let state = match Command_line.get "state" with | Some x -> int_of_string x | None -> 1 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 f12 = F12factor.gaussian_geminal expo in *) let f12 = F12factor.gaussian_geminal 1.0 in let simulation = Simulation.of_filenames ~f12 ~charge ~multiplicity ~nuclei:nuclei_file basis_file in let hf = HartreeFock.make simulation in Format.fprintf ppf "@[%a@]@." HartreeFock.pp hf; let mo_basis = MOBasis.of_hartree_fock hf in let fcif12 = F12CI.make ~simulation ~frozen_core:false ~mo_basis ~aux_basis_filename ~state () in let ci = F12CI.ci fcif12 in Format.fprintf ppf "FCI energy : "; Vec.iteri (fun i x -> if i <= state then Format.fprintf ppf "%20.16f@; " (x +. Simulation.nuclear_repulsion simulation) ) (CI.eigenvalues ci); Format.fprintf ppf "@."; let _, e_cif12 = F12CI.eigensystem fcif12 in Format.fprintf ppf "FCI-F12 energy : "; Vec.iteri (fun i x -> if i <= state then Format.fprintf ppf "%20.16f@; " (x +. Simulation.nuclear_repulsion simulation) ) e_cif12; Format.fprintf ppf "@."