From 92c23f042fb3a4f149d13d5c83ce8e04dd83253f Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Fri, 1 Jan 2021 11:46:11 +0100 Subject: [PATCH] Added simulation --- common/charge.org | 5 +- common/lib/charge.ml | 12 ++-- common/lib/charge.mli | 11 +-- particles/lib/nuclei.ml | 22 +++++- particles/lib/nuclei.mli | 3 +- particles/nuclei.org | 27 ++++++- qcaml/README.org | 38 ++++++---- qcaml/lib/qcaml.ml | 1 + simulation/lib/simulation.ml | 54 ++++++++++---- simulation/lib/simulation.mli | 50 ++++++++----- simulation/simulation.org | 131 ++++++++++++++++++++++++++++++++++ top/lib/install_printers.ml | 2 + 12 files changed, 300 insertions(+), 56 deletions(-) create mode 100644 simulation/simulation.org diff --git a/common/charge.org b/common/charge.org index fcc36f9..c704d49 100644 --- a/common/charge.org +++ b/common/charge.org @@ -65,6 +65,7 @@ val ( + ) : t -> t -> t val ( - ) : t -> t -> t val ( * ) : t -> float -> t val ( / ) : t -> float -> t +val is_null : t -> bool #+end_src #+begin_src ocaml :tangle (eval ml) :exports none @@ -77,6 +78,8 @@ let ( + ) = gen_op ( +. ) let ( - ) = gen_op ( -. ) let ( * ) = gen_op ( *. ) let ( / ) = gen_op ( /. ) + +let is_null t = t == 0. #+end_src ** Printers @@ -87,6 +90,6 @@ val pp : Format.formatter -> t -> unit #+begin_src ocaml :tangle (eval ml) :exports none let pp ppf x = - Format.fprintf ppf "@[+%s@]" (to_string x) + Format.fprintf ppf "@[%s@]" (to_string x) #+end_src diff --git a/common/lib/charge.ml b/common/lib/charge.ml index 9312540..a3cb3b2 100644 --- a/common/lib/charge.ml +++ b/common/lib/charge.ml @@ -3,11 +3,11 @@ (* This type should be used for all charges in the program (electrons, nuclei,...). *) -(* [[file:~/QCaml/common/charge.org::*Type][Type:2]] *) +(* [[file:../charge.org::*Type][Type:2]] *) type t = float (* Type:2 ends here *) -(* [[file:~/QCaml/common/charge.org::*Conversions][Conversions:2]] *) +(* [[file:../charge.org::*Conversions][Conversions:2]] *) external of_float : float -> t = "%identity" external to_float : t -> float = "%identity" @@ -25,7 +25,7 @@ let to_string x = "0.0" (* Conversions:2 ends here *) -(* [[file:~/QCaml/common/charge.org::*Simple%20operations][Simple operations:2]] *) +(* [[file:../charge.org::*Simple operations][Simple operations:2]] *) let gen_op op = fun a b -> op (to_float a) (to_float b) @@ -35,9 +35,11 @@ let ( + ) = gen_op ( +. ) let ( - ) = gen_op ( -. ) let ( * ) = gen_op ( *. ) let ( / ) = gen_op ( /. ) + +let is_null t = t == 0. (* Simple operations:2 ends here *) -(* [[file:~/QCaml/common/charge.org::*Printers][Printers:2]] *) +(* [[file:../charge.org::*Printers][Printers:2]] *) let pp ppf x = - Format.fprintf ppf "@[+%s@]" (to_string x) + Format.fprintf ppf "@[%s@]" (to_string x) (* Printers:2 ends here *) diff --git a/common/lib/charge.mli b/common/lib/charge.mli index 6a8dba7..762b08a 100644 --- a/common/lib/charge.mli +++ b/common/lib/charge.mli @@ -1,14 +1,14 @@ (* Type *) -(* [[file:~/QCaml/common/charge.org::*Type][Type:1]] *) +(* [[file:../charge.org::*Type][Type:1]] *) type t (* Type:1 ends here *) (* Conversions *) -(* [[file:~/QCaml/common/charge.org::*Conversions][Conversions:1]] *) +(* [[file:../charge.org::*Conversions][Conversions:1]] *) val of_float : float -> t val to_float : t -> float @@ -22,16 +22,17 @@ val to_string: t -> string (* Simple operations *) -(* [[file:~/QCaml/common/charge.org::*Simple%20operations][Simple operations:1]] *) +(* [[file:../charge.org::*Simple operations][Simple operations:1]] *) val ( + ) : t -> t -> t val ( - ) : t -> t -> t val ( * ) : t -> float -> t -val ( / ) : t -> float -> t +val ( / ) : t -> float -> t +val is_null : t -> bool (* Simple operations:1 ends here *) (* Printers *) -(* [[file:~/QCaml/common/charge.org::*Printers][Printers:1]] *) +(* [[file:../charge.org::*Printers][Printers:1]] *) val pp : Format.formatter -> t -> unit (* Printers:1 ends here *) diff --git a/particles/lib/nuclei.ml b/particles/lib/nuclei.ml index 0408799..1b1839f 100644 --- a/particles/lib/nuclei.ml +++ b/particles/lib/nuclei.ml @@ -119,13 +119,33 @@ let to_xyz_string t = -(* | ~repulsion~ | Nuclear repulsion energy, in atomic units | +(* | ~formula~ | Returns the chemical formula | + * | ~repulsion~ | Nuclear repulsion energy, in atomic units | * | ~charge~ | Sum of the charges of the nuclei | * | ~small_core~ | Number of core electrons in the small core model | * | ~large_core~ | Number of core electrons in the large core model | *) (* [[file:../nuclei.org::*Query][Query:2]] *) +let formula t = + let dict = Hashtbl.create 67 in + Array.iter (fun (e,_) -> + let e = Element.to_string e in + let value = + try (Hashtbl.find dict e) + 1 + with Not_found -> 1 + in + Hashtbl.replace dict e value + ) t; + Hashtbl.to_seq_keys dict + |> List.of_seq + |> List.sort String.compare + |> List.fold_left (fun accu key -> + let x = Hashtbl.find dict key in + accu ^ key ^ "_{" ^ (string_of_int x) ^ "}") "" + + + let repulsion nuclei = let get_charge e = Element.to_charge e diff --git a/particles/lib/nuclei.mli b/particles/lib/nuclei.mli index 74bc3d9..98eb078 100644 --- a/particles/lib/nuclei.mli +++ b/particles/lib/nuclei.mli @@ -24,10 +24,11 @@ val to_string : t -> string val of_filename : string -> t (* Conversion:1 ends here *) -(* Query *) +(* TODO Query *) (* [[file:../nuclei.org::*Query][Query:1]] *) +val formula : t -> string val repulsion : t -> float val charge : t -> Charge.t val small_core : t -> int diff --git a/particles/nuclei.org b/particles/nuclei.org index 6ba7f1a..2b3e278 100644 --- a/particles/nuclei.org +++ b/particles/nuclei.org @@ -24,7 +24,9 @@ type t = (Element.t * Coordinate.t) array #+end_src #+begin_src ocaml :tangle (eval ml) :exports none -<> +open Common + +type t = (Element.t * Coordinate.t) array open Xyz_ast #+end_src @@ -312,21 +314,42 @@ let to_xyz_string t = |> String.concat "\n" #+end_src -** Query +** TODO Query #+begin_src ocaml :tangle (eval mli) +val formula : t -> string val repulsion : t -> float val charge : t -> Charge.t val small_core : t -> int val large_core : t -> int #+end_src + | ~formula~ | Returns the chemical formula | | ~repulsion~ | Nuclear repulsion energy, in atomic units | | ~charge~ | Sum of the charges of the nuclei | | ~small_core~ | Number of core electrons in the small core model | | ~large_core~ | Number of core electrons in the large core model | #+begin_src ocaml :tangle (eval ml) :exports none +let formula t = + let dict = Hashtbl.create 67 in + Array.iter (fun (e,_) -> + let e = Element.to_string e in + let value = + try (Hashtbl.find dict e) + 1 + with Not_found -> 1 + in + Hashtbl.replace dict e value + ) t; + Hashtbl.to_seq_keys dict + |> List.of_seq + |> List.sort String.compare + |> List.fold_left (fun accu key -> + let x = Hashtbl.find dict key in + accu ^ key ^ "_{" ^ (string_of_int x) ^ "}") "" + + + let repulsion nuclei = let get_charge e = Element.to_charge e diff --git a/qcaml/README.org b/qcaml/README.org index cd5f77d..edefbfd 100644 --- a/qcaml/README.org +++ b/qcaml/README.org @@ -10,14 +10,15 @@ : Main QCaml entry point -* Dune files :noexport: +* File generation :noexport: - - Use [C-c C-c] on the code below to create the output for the dune files + Use [C-c C-c] on the code below. This will create the dune files + and the =qcaml.ml= file. #+header: :noweb yes #+header: :var name=(file-name-directory buffer-file-name) #+header: :var dune="lib/dune" + #+header: :var mlfile="lib/qcaml.ml" #+header: :var dunetest="test/dune" #+header: :var dependencies=dependencies #+begin_src python :exports none :results output none @@ -26,13 +27,15 @@ synopsis = """ <> """ -dependencies = '\n '.join(map(lambda x: x[0], dependencies)) +excluded_modules = [ "top" ] + +dependencies = '\n '.join([x[0] for x in dependencies if x[0].replace("qcaml.","") not in excluded_modules]) with open(dune,'w') as f: f.write(f""" (library - (name {name}) - (public_name qcaml.{name}) + (name qcaml) + (public_name qcaml) (synopsis {synopsis} ) (libraries {dependencies} @@ -41,19 +44,29 @@ with open(dune,'w') as f: <> ) <> -""".replace("qcaml.qcaml","qcaml")) +""") with open(dunetest,'w') as f: f.write(f""" (library - (name test_{name}) - (synopsis "Test for {name} library") + (name test_qcaml) + (synopsis "Test for qcaml library") (libraries alcotest - qcaml.{name} + qcaml ) ) -""".replace("qcaml.qcaml","qcaml")) +""") + +with open(mlfile,'w') as f: + s = dependencies.split() + s = [x.replace("qcaml.","") for x in s ] + s = [x.capitalize() for x in s] + s = [ f"module {x} = {x}\n" for x in s ] + s = ''.join(s) + f.write("(* Auto-generated by qcaml/README.org *)\n") + f.write(s) + #+end_src @@ -63,7 +76,6 @@ grep "public_name" ../*/lib/dune | grep -v "qcaml)" | cut -d ' ' -f 3 | tr ')' ' #+end_src #+RESULTS: dependencies - | (libraries | | qcaml.ao | | qcaml.common | | qcaml.gaussian | @@ -74,4 +86,4 @@ grep "public_name" ../*/lib/dune | grep -v "qcaml)" | cut -d ' ' -f 3 | tr ')' ' | qcaml.particles | | qcaml.perturbation | | qcaml.simulation | - | ) | + | qcaml.top | diff --git a/qcaml/lib/qcaml.ml b/qcaml/lib/qcaml.ml index cbf0c0c..e48edc4 100644 --- a/qcaml/lib/qcaml.ml +++ b/qcaml/lib/qcaml.ml @@ -1,3 +1,4 @@ +(* Auto-generated by qcaml/README.org *) module Ao = Ao module Common = Common module Gaussian = Gaussian diff --git a/simulation/lib/simulation.ml b/simulation/lib/simulation.ml index dc8e314..c15e245 100644 --- a/simulation/lib/simulation.ml +++ b/simulation/lib/simulation.ml @@ -1,7 +1,10 @@ +(* [[file:../simulation.org::*Simulation][Simulation:2]] *) open Common open Particles open Operators - +(* Simulation:2 ends here *) + +(* [[file:../simulation.org::*Type][Type:2]] *) type t = { charge : Charge.t; electrons : Electrons.t; @@ -9,19 +12,42 @@ type t = { ao_basis : Ao.Basis.t; operators : Operator.t list; } +(* Type:2 ends here *) + + +(* | ~nuclei~ | Nuclear coordinates used in the smiulation | + * | ~charge~ | Total charge (electrons + nuclei) | + * | ~electrons~ | Electrons used in the simulation | + * | ~ao_basis~ | Atomic basis set | + * | ~nuclear_repulsion~ | Nuclear repulsion energy | + * | ~operators~ | List of extra operators (range-separation, f12, etc) | *) + + +(* [[file:../simulation.org::*Access][Access:2]] *) let nuclei t = t.nuclei let charge t = t.charge let electrons t = t.electrons let ao_basis t = t.ao_basis let nuclear_repulsion t = Nuclei.repulsion @@ nuclei t let operators t = t.operators +(* Access:2 ends here *) -let make ?(multiplicity=1) - ?(charge=0) - ?(operators=[]) - ~nuclei - ao_basis + + +(* Defaults: + * - multiplicity : ~1~ + * - charge : ~0~ + * - operators : ~[]~ *) + + +(* [[file:../simulation.org::*Creation][Creation:2]] *) +let make + ?(multiplicity=1) + ?(charge=0) + ?(operators=[]) + ~nuclei + ao_basis = (* Tune Garbage Collector *) @@ -36,10 +62,14 @@ let make ?(multiplicity=1) Charge.(Nuclei.charge nuclei + Electrons.charge electrons) in - { - charge ; nuclei ; electrons ; ao_basis ; - operators - } - - + { charge ; nuclei ; electrons ; ao_basis ; operators} +(* Creation:2 ends here *) +(* [[file:../simulation.org::*Printers][Printers:2]] *) +let pp ppf t = + let formula = Nuclei.formula t.nuclei in + let n_aos = Ao.Basis.size t.ao_basis in + let n_ops = List.length t.operators in + Format.fprintf ppf "@[@[%s@], @[%a@], @[%d AOs@], @[%d operators@]@]" + formula Electrons.pp t.electrons n_aos n_ops +(* Printers:2 ends here *) diff --git a/simulation/lib/simulation.mli b/simulation/lib/simulation.mli index 285489e..6c55e95 100644 --- a/simulation/lib/simulation.mli +++ b/simulation/lib/simulation.mli @@ -1,32 +1,50 @@ -(* Contains the state of a simulation *) +(* Simulation + * :PROPERTIES: + * :header-args: :noweb yes :comments both + * :END: + * + * Contains the state of the simulation. + * + * #+NAME: open *) +(* [[file:../simulation.org::open][open]] *) open Common open Particles open Operators - +(* open ends here *) + +(* Type + * + * #+NAME: types *) + +(* [[file:../simulation.org::types][types]] *) type t +(* types ends here *) -val nuclei : t -> Nuclei.t -(** Nuclear coordinates used in the smiulation *) +(* Access *) -val charge : t -> Charge.t -(** Total charge (electrons + nuclei) *) - -val electrons : t -> Electrons.t -(** Electrons used in the simulation *) - -val ao_basis : t -> Ao.Basis.t -(** Atomic basis set *) +(* [[file:../simulation.org::*Access][Access:1]] *) +val nuclei : t -> Nuclei.t +val charge : t -> Charge.t +val electrons : t -> Electrons.t +val ao_basis : t -> Ao.Basis.t val nuclear_repulsion : t -> float -(** Nuclear repulsion energy *) +val operators : t -> Operator.t list +(* Access:1 ends here *) -val operators : t -> Operator.t list -(** List of extra operators (range-separation, f12, etc) *) +(* Creation *) -(** {1 Creation} *) +(* [[file:../simulation.org::*Creation][Creation:1]] *) val make : ?multiplicity:int -> ?charge:int -> ?operators:Operator.t list-> nuclei:Nuclei.t -> Ao.Basis.t -> t +(* Creation:1 ends here *) +(* Printers *) + + +(* [[file:../simulation.org::*Printers][Printers:1]] *) +val pp : Format.formatter -> t -> unit +(* Printers:1 ends here *) diff --git a/simulation/simulation.org b/simulation/simulation.org new file mode 100644 index 0000000..2176052 --- /dev/null +++ b/simulation/simulation.org @@ -0,0 +1,131 @@ +#+begin_src elisp tangle: no :results none :exports none +(setq pwd (file-name-directory buffer-file-name)) +(setq name (file-name-nondirectory (substring buffer-file-name 0 -4))) +(setq lib (concat pwd "lib/")) +(setq testdir (concat pwd "test/")) +(setq mli (concat lib name ".mli")) +(setq ml (concat lib name ".ml")) +(setq test-ml (concat testdir name ".ml")) +(org-babel-tangle) +#+end_src + +* Simulation + :PROPERTIES: + :header-args: :noweb yes :comments both + :END: + + Contains the state of the simulation. + + #+NAME: open + #+begin_src ocaml :tangle (eval mli) +open Common +open Particles +open Operators + #+end_src + + #+begin_src ocaml :tangle (eval ml) :exports none +<> + #+end_src + +** Type + + #+NAME: types + #+begin_src ocaml :tangle (eval mli) +type t + #+end_src + + #+begin_src ocaml :tangle (eval ml) :exports none +type t = { + charge : Charge.t; + electrons : Electrons.t; + nuclei : Nuclei.t; + ao_basis : Ao.Basis.t; + operators : Operator.t list; +} + #+end_src + +** Access + + #+begin_src ocaml :tangle (eval mli) +val nuclei : t -> Nuclei.t +val charge : t -> Charge.t +val electrons : t -> Electrons.t +val ao_basis : t -> Ao.Basis.t +val nuclear_repulsion : t -> float +val operators : t -> Operator.t list + #+end_src + + | ~nuclei~ | Nuclear coordinates used in the smiulation | + | ~charge~ | Total charge (electrons + nuclei) | + | ~electrons~ | Electrons used in the simulation | + | ~ao_basis~ | Atomic basis set | + | ~nuclear_repulsion~ | Nuclear repulsion energy | + | ~operators~ | List of extra operators (range-separation, f12, etc) | + + #+begin_src ocaml :tangle (eval ml) :exports none +let nuclei t = t.nuclei +let charge t = t.charge +let electrons t = t.electrons +let ao_basis t = t.ao_basis +let nuclear_repulsion t = Nuclei.repulsion @@ nuclei t +let operators t = t.operators + #+end_src + +** Creation + + #+begin_src ocaml :tangle (eval mli) +val make : ?multiplicity:int -> ?charge:int -> + ?operators:Operator.t list-> nuclei:Nuclei.t -> + Ao.Basis.t -> t + #+end_src + + Defaults: + - multiplicity : ~1~ + - charge : ~0~ + - operators : ~[]~ + + #+begin_src ocaml :tangle (eval ml) :exports none +let make + ?(multiplicity=1) + ?(charge=0) + ?(operators=[]) + ~nuclei + ao_basis + = + + (* Tune Garbage Collector *) + let gc = Gc.get () in + Gc.set { gc with space_overhead = 1000 }; + + let electrons = + Electrons.of_atoms ~multiplicity ~charge nuclei + in + + let charge = + Charge.(Nuclei.charge nuclei + Electrons.charge electrons) + in + + { charge ; nuclei ; electrons ; ao_basis ; operators} + #+end_src + +** Printers + + #+begin_src ocaml :tangle (eval mli) +val pp : Format.formatter -> t -> unit + #+end_src + + #+begin_src ocaml :tangle (eval ml) :exports none +let pp ppf t = + let formula = Nuclei.formula t.nuclei in + let n_aos = Ao.Basis.size t.ao_basis in + let n_ops = List.length t.operators in + Format.fprintf ppf "@[@[%s@], @[%a@], @[%d AOs@], @[%d operators@]@]" + formula Electrons.pp t.electrons n_aos n_ops + #+end_src + + #+RESULTS: + : Line 2, characters 16-30: + : 2 | let formula = Nuclei.formula t.nuclei in + : ^^^^^^^^^^^^^^ + : Error: Unbound module Nuclei + diff --git a/top/lib/install_printers.ml b/top/lib/install_printers.ml index ae25f4f..894ee15 100644 --- a/top/lib/install_printers.ml +++ b/top/lib/install_printers.ml @@ -11,7 +11,9 @@ let printers = "Common.Zkey.pp" ; "Particles.Electrons.pp" ; "Particles.Element.pp" ; + "Particles.Nuclei.pp" ; "Particles.Zmatrix.pp" ; + "Simulation.Simulation.pp" ; ]