From d6b970e8a73de793e0a7cd5423a7931e88ea873c Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Tue, 22 Mar 2016 13:28:03 +0100 Subject: [PATCH] Auto-generation of libint required files --- install/scripts/install_ocaml.sh | 2 +- ocaml/Atom.ml | 12 ++- ocaml/Atom.mli | 1 + ocaml/Basis.ml | 22 ++++- ocaml/Basis.mli | 2 +- ocaml/Gto.ml | 35 +++++++- ocaml/Gto.mli | 6 +- ocaml/Input_ao_basis.ml | 1 + ocaml/Input_nuclei.ml | 17 ++++ ocaml/Libint.ml | 85 +++++++++++++++++++ ocaml/Molecule.ml | 11 ++- ocaml/Molecule.mli | 1 + ocaml/TaskServer.mli | 84 ++++++++++++++++++ ocaml/_tags | 2 +- ...debug_libinit.irp.f => debug_libint.irp.f} | 64 ++++++-------- scripts/compilation/qp_create_ninja.py | 2 +- scripts/ezfio_interface/qp_edit_template | 11 ++- src/Integrals_Bielec/integral_bielec.cc | 4 +- src/Integrals_Monoelec/shell.irp.f | 30 +++++++ 19 files changed, 332 insertions(+), 60 deletions(-) create mode 100644 ocaml/Libint.ml create mode 100644 ocaml/TaskServer.mli rename plugins/Hartree_Fock/{debug_libinit.irp.f => debug_libint.irp.f} (68%) create mode 100644 src/Integrals_Monoelec/shell.irp.f diff --git a/install/scripts/install_ocaml.sh b/install/scripts/install_ocaml.sh index 86e4e8b7..234cd09f 100755 --- a/install/scripts/install_ocaml.sh +++ b/install/scripts/install_ocaml.sh @@ -5,7 +5,7 @@ QP_ROOT=$PWD cd - # Normal installation -PACKAGES="core cryptokit ocamlfind sexplib ZMQ" +PACKAGES="core cryptokit ocamlfind sexplib pa_sexp_conv ZMQ" # Needed for ZeroMQ export C_INCLUDE_PATH="${QP_ROOT}"/lib:"${C_INCLUDE_PATH}" diff --git a/ocaml/Atom.ml b/ocaml/Atom.ml index 832cfa5b..72932b1f 100644 --- a/ocaml/Atom.ml +++ b/ocaml/Atom.ml @@ -1,4 +1,4 @@ -open Core.Std;; +open Core.Std exception AtomError of string @@ -27,12 +27,18 @@ let of_string ~units s = coord = Point3d.of_string ~units (String.concat [x; y; z] ~sep:" ") } | _ -> raise (AtomError s) -;; + let to_string ~units a = [ Element.to_string a.element ; Charge.to_string a.charge ; Point3d.to_string ~units a.coord ] |> String.concat ~sep:" " -;; + + +let to_xyz a = + Printf.sprintf "%-3s %s" + (Element.to_string a.element) + (Point3d.to_string ~units:Units.Angstrom a.coord) + diff --git a/ocaml/Atom.mli b/ocaml/Atom.mli index 28915993..4b1963d5 100644 --- a/ocaml/Atom.mli +++ b/ocaml/Atom.mli @@ -7,3 +7,4 @@ val sexp_of_t : t -> Sexplib.Sexp.t val of_string : units:Units.units -> string -> t val to_string : units:Units.units -> t -> string +val to_xyz : t -> string diff --git a/ocaml/Basis.ml b/ocaml/Basis.ml index 237e5547..c2675c57 100644 --- a/ocaml/Basis.ml +++ b/ocaml/Basis.ml @@ -35,11 +35,11 @@ let read_element in_channel at_number element = read in_channel at_number -let to_string b = + +let to_string_general ~fmt ~atom_sep b = let new_nucleus n = Printf.sprintf "Atom %d" n in - let rec do_work accu current_nucleus = function | [] -> List.rev accu | (g,n)::tail -> @@ -47,15 +47,29 @@ let to_string b = in let accu = if (n <> current_nucleus) then - (new_nucleus n)::""::accu + (new_nucleus n)::atom_sep::accu else accu in - do_work ((Gto.to_string g)::accu) n tail + do_work ((Gto.to_string ~fmt g)::accu) n tail in do_work [new_nucleus 1] 1 b |> String.concat ~sep:"\n" +let to_string_gamess = + to_string_general ~fmt:Gto.Gamess ~atom_sep:"" + +let to_string_gaussian b = + String.concat ~sep:"\n" [ "****" ; + to_string_general ~fmt:Gto.Gaussian ~atom_sep:"****" b ; + "****" + ] + +let to_string ?(fmt=Gto.Gamess) = + match fmt with + | Gto.Gamess -> to_string_gamess + | Gto.Gaussian -> to_string_gaussian + include To_md5 let to_md5 = to_md5 sexp_of_t diff --git a/ocaml/Basis.mli b/ocaml/Basis.mli index 4da99266..249c14f9 100644 --- a/ocaml/Basis.mli +++ b/ocaml/Basis.mli @@ -14,7 +14,7 @@ val read_element : in_channel -> Nucl_number.t -> Element.t -> (Gto.t * Nucl_number.t) list (** Convert the basis to a string *) -val to_string : (Gto.t * Nucl_number.t) list -> string +val to_string : ?fmt:Gto.fmt -> (Gto.t * Nucl_number.t) list -> string (** Convert the basis to an MD5 hash *) val to_md5 : (Gto.t * Nucl_number.t) list -> MD5.t diff --git a/ocaml/Gto.ml b/ocaml/Gto.ml index 69aeba37..fb576ee7 100644 --- a/ocaml/Gto.ml +++ b/ocaml/Gto.ml @@ -4,6 +4,10 @@ open Qptypes exception GTO_Read_Failure of string exception End_Of_Basis +type fmt = +| Gamess +| Gaussian + type t = { sym : Symmetry.t ; lc : ((Primitive.t * AO_coef.t) list) @@ -68,8 +72,8 @@ let read_one in_channel = -(** Transform the gto to a string *) -let to_string { sym = sym ; lc = lc } = +(** Write the GTO in Gamess format *) +let to_string_gamess { sym = sym ; lc = lc } = let result = Printf.sprintf "%s %3d" (Symmetry.to_string sym) (List.length lc) in @@ -88,3 +92,30 @@ let to_string { sym = sym ; lc = lc } = |> String.concat ~sep:"\n" +(** Write the GTO in Gaussian format *) +let to_string_gaussian { sym = sym ; lc = lc } = + let result = + Printf.sprintf "%s %3d 1.00" (Symmetry.to_string sym) (List.length lc) + in + let rec do_work accu i = function + | [] -> List.rev accu + | (p,c)::tail -> + let p = AO_expo.to_float p.Primitive.expo + and c = AO_coef.to_float c + in + let result = + Printf.sprintf "%15.7f %15.7f" p c + in + do_work (result::accu) (i+1) tail + in + (do_work [result] 1 lc) + |> String.concat ~sep:"\n" + + +(** Transform the gto to a string *) +let to_string ?(fmt=Gamess) = + match fmt with + | Gamess -> to_string_gamess + | Gaussian -> to_string_gaussian + + diff --git a/ocaml/Gto.mli b/ocaml/Gto.mli index fad133a3..753cd81a 100644 --- a/ocaml/Gto.mli +++ b/ocaml/Gto.mli @@ -1,5 +1,9 @@ exception GTO_Read_Failure of string exception End_Of_Basis +type fmt = +| Gamess +| Gaussian + type t = { sym : Symmetry.t ; lc : (Primitive.t * Qptypes.AO_coef.t) list; @@ -13,4 +17,4 @@ val of_prim_coef_list : val read_one : in_channel -> t (** Convert to string for printing *) -val to_string : t -> string +val to_string : ?fmt:fmt -> t -> string diff --git a/ocaml/Input_ao_basis.ml b/ocaml/Input_ao_basis.ml index 82bc4964..88e277ee 100644 --- a/ocaml/Input_ao_basis.ml +++ b/ocaml/Input_ao_basis.ml @@ -17,6 +17,7 @@ module Ao_basis : sig ;; val read : unit -> t option val to_string : t -> string + val to_basis : t -> Basis.t val write : t -> unit val to_md5 : t -> MD5.t val to_rst : t -> Rst_string.t diff --git a/ocaml/Input_nuclei.ml b/ocaml/Input_nuclei.ml index d050ded9..ca81629e 100644 --- a/ocaml/Input_nuclei.ml +++ b/ocaml/Input_nuclei.ml @@ -13,6 +13,7 @@ module Nuclei : sig val read : unit -> t option val write : t -> unit val to_string : t -> string + val to_atom_list : t -> Atom.t list val to_rst : t -> Rst_string.t val of_rst : Rst_string.t -> t option end = struct @@ -134,6 +135,22 @@ end = struct ;; + let to_atom_list b = + let rec loop accu (coord, charge, label) = function + | -1 -> accu + | i -> + let atom = + { Atom.element = label.(i) ; + Atom.charge = charge.(i) ; + Atom.coord = coord.(i) ; + } + in + loop (atom::accu) (coord, charge, label) (i-1) + in + loop [] (b.nucl_coord, b.nucl_charge, b.nucl_label) + ( (Nucl_number.to_int b.nucl_num) - 1) + ;; + let to_string b = Printf.sprintf " nucl_num = %s diff --git a/ocaml/Libint.ml b/ocaml/Libint.ml new file mode 100644 index 00000000..0b88ab41 --- /dev/null +++ b/ocaml/Libint.ml @@ -0,0 +1,85 @@ +open Core.Std +open Qptypes + +let molecule = lazy ( + let atom_list = + match Input.Nuclei.read () with + | Some data -> Input.Nuclei.to_atom_list data + | None -> failwith "No coordinate found" + and data = + match Input.Electrons.read () with + | Some data -> data + | None -> failwith "No electrons found" + in + { Molecule.nuclei = atom_list; + Molecule.elec_alpha = data.Input.Electrons.elec_alpha_num; + Molecule.elec_beta = data.Input.Electrons.elec_beta_num + } +) + +let write_xyz_file libint_dir = + assert (Sys.is_directory_exn libint_dir); + let filename = + Filename.concat libint_dir "xyz" + and molecule = + Lazy.force molecule + in + Out_channel.with_file filename ~f:(fun oc -> + Molecule.to_xyz molecule + |> Out_channel.output_string oc + ) + +let write_basis_file libint_dir = + assert (Sys.is_directory_exn libint_dir); + let filename = + Filename.concat libint_dir "basis.g94" + and molecule = + Lazy.force molecule + in + let text = + let rec substitute accu i = function + | ele :: tail -> + let pattern = + Printf.sprintf "Atom %d" i + in + let new_string = + String.substr_replace_first accu ~pattern ~with_:(Printf.sprintf "%s 0" ele) + in + substitute new_string (i+1) tail + | [] -> accu + in + let accu = + let b = + match Input.Ao_basis.read () with + | Some data -> Input.Ao_basis.to_basis data + | None -> failwith "No AO basis" + in + Basis.to_string ~fmt:Gto.Gaussian b + and atom_names = + List.map molecule.Molecule.nuclei ~f:(fun x -> Element.to_string x.Atom.element) + in + substitute accu 1 atom_names + in + Out_channel.with_file filename ~f:(fun oc -> + Out_channel.output_string oc text + ) + + +let write_files ezfio_filename = + assert (Sys.is_directory_exn ezfio_filename); + + let libint_dir = + Filename.concat ezfio_filename "libint" + in + + let () = + match Sys.is_directory libint_dir with + | `Yes -> () + | `No -> Unix.mkdir libint_dir + | `Unknown -> failwith ("Unable to tell if "^libint_dir^" exists.") + in + + write_xyz_file libint_dir; + write_basis_file libint_dir + + diff --git a/ocaml/Molecule.ml b/ocaml/Molecule.ml index f0800f7f..a9d73432 100644 --- a/ocaml/Molecule.ml +++ b/ocaml/Molecule.ml @@ -85,7 +85,7 @@ let name m = String.concat (result) -let to_string m = +let to_string_general ~f m = let { nuclei ; elec_alpha ; elec_beta } = m in let n = @@ -94,10 +94,15 @@ let to_string m = let title = name m in - [ Int.to_string n ; title ] @ - (List.map ~f:(fun x -> Atom.to_string Units.Angstrom x) nuclei) + [ Int.to_string n ; title ] @ (List.map ~f nuclei) |> String.concat ~sep:"\n" +let to_string = + to_string_general ~f:(fun x -> Atom.to_string Units.Angstrom x) + +let to_xyz = + to_string_general ~f:Atom.to_xyz + let of_xyz_string ?(charge=(Charge.of_int 0)) ?(multiplicity=(Multiplicity.of_int 1)) diff --git a/ocaml/Molecule.mli b/ocaml/Molecule.mli index 1a3d9715..f81f28a3 100644 --- a/ocaml/Molecule.mli +++ b/ocaml/Molecule.mli @@ -20,6 +20,7 @@ val name : t -> string (** Conversion for printing *) val to_string : t -> string +val to_xyz : t -> string (** Creates a molecule from an xyz file *) diff --git a/ocaml/TaskServer.mli b/ocaml/TaskServer.mli new file mode 100644 index 00000000..f16ddaab --- /dev/null +++ b/ocaml/TaskServer.mli @@ -0,0 +1,84 @@ +type t = +{ + queue : Queuing_system.t ; + state : Message.State.t option ; + address_tcp : Address.Tcp.t option ; + address_inproc : Address.Inproc.t option ; + psi : Message.Psi.t option; + progress_bar : Progress_bar.t option ; + running : bool; +} + + +(** {1} Debugging *) + +(** Fetch the QP_TASK_DEBUG environment variable *) +val debug_env : bool + +(** Print a debug message *) +val debug : string -> unit + +(** {1} ZMQ *) + +(** ZeroMQ context *) +val zmq_context : ZMQ.Context.t + +(** Bind a ZMQ socket *) +val bind_socket : + socket_type:string -> socket:'a ZMQ.Socket.t -> address:string -> unit + +(** Name of the host on which the server runs *) +val hostname : string lazy_t + +(** IP address of the current host *) +val ip_address : string lazy_t + +(** Standard messages *) +val reply_ok : [> `Req ] ZMQ.Socket.t -> unit +val reply_wrong_state : [> `Req ] ZMQ.Socket.t -> unit + +(** Stop server *) +val stop : port:int -> unit + +(** {1} Server functions *) + +(** Create a new job *) +val new_job : Message.Newjob_msg.t -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Finish a running job *) +val end_job : Message.Endjob_msg.t -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Connect a client *) +val connect: Message.Connect_msg.t -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Disconnect a client *) +val disconnect: Message.Disconnect_msg.t -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Add a task to the pool *) +val add_task: Message.AddTask_msg.t -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Mark the task as done by the client *) +val task_done: Message.TaskDone_msg.t -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Delete a task when it has been pulled by the collector *) +val del_task: Message.DelTask_msg.t -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** The client get a new task to execute *) +val get_task: Message.GetTask_msg.t -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Terminate server *) +val terminate : t -> [> `Req ] ZMQ.Socket.t -> t + +(** Put a wave function in the task server *) +val put_psi : + Message.PutPsi_msg.t -> string list -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Get the wave function stored in the task server *) +val get_psi : Message.GetPsi_msg.t -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Reply an Error message *) +val error : string -> t -> [> `Req ] ZMQ.Socket.t -> t + +(** Run server *) +val run : port:int -> unit + diff --git a/ocaml/_tags b/ocaml/_tags index fd4c4804..6570ddca 100644 --- a/ocaml/_tags +++ b/ocaml/_tags @@ -1,3 +1,3 @@ -true: package(core,sexplib.syntax,cryptokit,ZMQ) +true: package(core,sexplib,pa_sexp_conv,cryptokit,ZMQ) true: thread false: profile diff --git a/plugins/Hartree_Fock/debug_libinit.irp.f b/plugins/Hartree_Fock/debug_libint.irp.f similarity index 68% rename from plugins/Hartree_Fock/debug_libinit.irp.f rename to plugins/Hartree_Fock/debug_libint.irp.f index 586f393e..bfbb4027 100644 --- a/plugins/Hartree_Fock/debug_libinit.irp.f +++ b/plugins/Hartree_Fock/debug_libint.irp.f @@ -1,20 +1,10 @@ program debug_libint - use libint_module implicit none double precision :: ao_bielec_integral - integer, allocatable :: s2bf(:) double precision, allocatable :: buffer_int(:) - call init_libint(trim(ezfio_filename)//char(0)) - - integer :: nb_shell_f - nb_shell_f = get_nb_shell() - - allocate(s2bf(2*nb_shell_f)) - call map_shell_to_basis_function_interval(2*nb_shell_f,s2bf) - integer :: s1, s2,s3,s4 integer :: bf1,bf2,bf3,bf4 integer :: bf1_begin,bf2_begin,bf3_begin,bf4_begin @@ -26,30 +16,30 @@ ! Loop over the shell ! ! =================== ! - do s1 = 1, nb_shell_f + do s1 = 1, shell_num - print*, s1, "/", nb_shell_f + print*, s1, "/", shell_num - bf1_begin = s2bf(2*s1-1) - bf1_end = s2bf(2*s1) + bf1_begin = shell_idx(1,s1) + bf1_end = shell_idx(2,s1) n1 = 1 + bf1_end - bf1_begin - do s2 = 1, nb_shell_f + do s2 = 1, shell_num - bf2_begin = s2bf(2*s2-1) - bf2_end = s2bf(2*s2) + bf2_begin = shell_idx(1,s2) + bf2_end = shell_idx(2,s2) n2 = 1 + bf2_end - bf2_begin - do s3 = 1, nb_shell_f + do s3 = 1, shell_num - bf3_begin = s2bf(2*s3-1) - bf3_end = s2bf(2*s3) + bf3_begin = shell_idx(1,s3) + bf3_end = shell_idx(2,s3) n3 = 1 + bf3_end - bf3_begin - do s4 = 1, nb_shell_f + do s4 = 1, shell_num - bf4_begin = s2bf(2*s4-1) - bf4_end = s2bf(2*s4) + bf4_begin = shell_idx(1,s4) + bf4_end = shell_idx(2,s4) n4 = 1 + bf4_end - bf4_begin ! ========================== ! @@ -93,19 +83,19 @@ libint = buffer_int(f1234) * norm !Verify with the manu's one -! ref = ao_bielec_integral(bf1,bf2,bf3,bf4) -! -! if ( (ABS(ABS(ref) - ABS(libint)).GE.1e-6) ) THEN -! -! print*, bf1,bf2,bf3,bf4 -! print*,"r", ref -! print*,"l", libint -! print*,"r/l", ref/libint -! print*,"l/r", libint/ref -! print*,"n", norm -! -! call exit(1) -! end if + ref = ao_bielec_integral(bf1,bf2,bf3,bf4) + + if ( (ABS(ABS(ref) - ABS(libint)).GE.1e-6) ) THEN + + print*, bf1,bf2,bf3,bf4 + print*,"r", ref + print*,"l", libint + print*,"r/l", ref/libint + print*,"l/r", libint/ref + print*,"n", norm + + call exit(1) + end if enddo enddo @@ -121,4 +111,4 @@ enddo call finalize_libint() - end debug_libint + end diff --git a/scripts/compilation/qp_create_ninja.py b/scripts/compilation/qp_create_ninja.py index df2a1dff..946b184c 100755 --- a/scripts/compilation/qp_create_ninja.py +++ b/scripts/compilation/qp_create_ninja.py @@ -38,7 +38,7 @@ from qp_path import QP_ROOT, QP_SRC, QP_EZFIO LIB = "" # join(QP_ROOT, "lib", "rdtsc.o") EZFIO_LIB = join(QP_ROOT, "lib", "libezfio_irp.a") ZMQ_LIB = join(QP_ROOT, "lib", "libf77zmq.a") + " " + join(QP_ROOT, "lib", "libzmq.a") + " -lstdc++ -lrt" -INT_LIB = join(QP_ROOT, "libint","lib", ".libs", "libint2.a") +INT_LIB = join(QP_ROOT, "install", "libint","lib", ".libs", "libint2.a") ROOT_BUILD_NINJA = join(QP_ROOT, "config", "build.ninja") header = r"""# diff --git a/scripts/ezfio_interface/qp_edit_template b/scripts/ezfio_interface/qp_edit_template index 408ca3f7..b580355b 100644 --- a/scripts/ezfio_interface/qp_edit_template +++ b/scripts/ezfio_interface/qp_edit_template @@ -6,7 +6,7 @@ open Core.Std;; WARNING This file is autogenerad by -`${{QP_ROOT}}/script/ezfio_interface/ei_handler.py` +`${{QP_ROOT}}/scripts/ezfio_interface/ei_handler.py` *) @@ -120,7 +120,7 @@ let set str s = | Nuclei -> write Nuclei.(of_rst, write) s | Ao_basis -> () (* TODO *) | Mo_basis -> () (* TODO *) - end + end ;; @@ -169,7 +169,9 @@ let run check_only ezfio_filename = in (* Create the temp file *) - let temp_filename = create_temp_file ezfio_filename tasks in + let temp_filename = + create_temp_file ezfio_filename tasks + in (* Open the temp file with external editor *) let editor = @@ -191,9 +193,10 @@ let run check_only ezfio_filename = In_channel.input_all in_channel) in List.iter ~f:(fun x -> set temp_string x) tasks; + Libint.write_files (!Ezfio.ezfio_filename); (* Remove temp_file *) - Sys.remove temp_filename; + Sys.remove temp_filename ;; diff --git a/src/Integrals_Bielec/integral_bielec.cc b/src/Integrals_Bielec/integral_bielec.cc index 1ea9f160..e63345f9 100644 --- a/src/Integrals_Bielec/integral_bielec.cc +++ b/src/Integrals_Bielec/integral_bielec.cc @@ -18,7 +18,7 @@ */ #if __cplusplus <= 199711L -# error "Hartree-Fock test requires C++11 support" +# error "C++11 support is required" #endif // standard C++ headers @@ -244,4 +244,4 @@ void compute_ao_bielec_integrals_shell(int s1, int s2, int s3, int s4, int sze, for(auto i=0; i!=sze; i++) values[i] = buf_1234[i]; - }; +}; diff --git a/src/Integrals_Monoelec/shell.irp.f b/src/Integrals_Monoelec/shell.irp.f new file mode 100644 index 00000000..b91476e1 --- /dev/null +++ b/src/Integrals_Monoelec/shell.irp.f @@ -0,0 +1,30 @@ +use libint_module +BEGIN_PROVIDER [ logical, libint_initialized ] + implicit none + BEGIN_DOC + ! if true, libint is initialized + END_DOC + call init_libint(trim(ezfio_filename)//char((0))) +END_PROVIDER + +BEGIN_PROVIDER [ integer, shell_num ] + implicit none + BEGIN_DOC + ! Number of shells + END_DOC + PROVIDE libint_initialized + shell_num = get_nb_shell() + +END_PROVIDER + + +BEGIN_PROVIDER [ integer, shell_idx, (2,shell_num) ] + implicit none + BEGIN_DOC + ! Contains the 1st and last AO index in each shell + END_DOC + PROVIDE libint_initialized + call map_shell_to_basis_function_interval(2*shell_num,shell_idx) +END_PROVIDER + +