diff --git a/Makefile b/Makefile index 8d4cd8e..5e1b13a 100644 --- a/Makefile +++ b/Makefile @@ -11,16 +11,11 @@ DOCS=$(patsubst %, docs/%.html, $(DIRS)) docs/index.html default: doc build -docs/index.html: docs/index.org - - ./bin/build_doc.sh docs - docs/%.html: %/*.org %/lib/*.ml %/lib/*.mli %/test/*.ml - ./bin/tangle.sh $* - - ./bin/build_doc.sh $* docs/top.html: */*.org */lib/*.ml */lib/*.mli - ./bin/tangle.sh top - - ./bin/build_doc.sh top doc: $(DOCS) diff --git a/ao/basis.org b/ao/basis.org deleted file mode 100644 index 01a3a57..0000000 --- a/ao/basis.org +++ /dev/null @@ -1,304 +0,0 @@ -#+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 - -* Basis - :PROPERTIES: - :header-args: :noweb yes :comments both - :END: - - Data structure for Atomic Orbitals. - -** Dimensions :noexport: - - #+begin_src ocaml :tangle lib/ao_dim.mli :exports none -type t - #+end_src - -** Polymorphic types - - <<<~Basis.t~>>> - #+NAME: types - #+begin_src ocaml :tangle lib/basis_poly.mli -type t = - | Unknown - | Gaussian of Basis_gaussian.t - #+end_src - - #+begin_src ocaml :tangle lib/basis_poly.ml :exports none -<> - #+end_src - -** Types - - <<<~Basis.t~>>> - #+begin_src ocaml :tangle (eval mli) -type t -type ao = Ao_dim.t - -open Common -open Particles -open Operators -open Linear_algebra - #+end_src - - #+begin_src ocaml :tangle (eval ml) :exports none -type t = - { ao_basis : Basis_poly.t ; - cartesian : bool - } - -type ao = Ao_dim.t - -open Linear_algebra -open Common - #+end_src - -** Conversions - - #+begin_src ocaml :tangle (eval mli) -val of_nuclei_and_basis_filename : - ?kind:[> `Gaussian ] -> - ?operators:Operator.t list -> - ?cartesian:bool -> - nuclei:Nuclei.t -> - string -> - t - #+end_src - - |--------------------------------+------------------------------------------------------------------------------------------------------------------------| - | ~of_nuclei_and_basis_filename~ | Creates the data structure for the atomic orbitals basis from a molecule ~Nuclei.t~ and the name of the basis-set file | - |--------------------------------+------------------------------------------------------------------------------------------------------------------------| - - Defaults: - - ~kind~ : ~`Gaussian~ - - ~operators~ : ~[]~ - - ~cartesian~ : ~false~ - - Example: - #+begin_example -let b = Ao.Basis.of_nuclei_and_basis_filename ~nuclei filename;; -val b : Ao.Basis.t = Gaussian Basis, spherical, 15 AOs - #+end_example - - #+begin_src ocaml :tangle (eval ml) :exports none -let not_implemented () = - Util.not_implemented "Only Gaussian is implemented" - -let of_nuclei_and_basis_filename ?(kind=`Gaussian) ?operators ?(cartesian=false) - ~nuclei filename = - match kind with - | `Gaussian -> - let basis = - Gaussian.Basis.of_nuclei_and_basis_filename ~nuclei filename - in - let ao_basis = - Basis_poly.Gaussian (Basis_gaussian.make ~basis ?operators ~cartesian nuclei ) - in - { ao_basis ; cartesian } - | _ -> not_implemented () - #+end_src - -** Access - - #+begin_src ocaml :tangle (eval mli) -val size : t -> int -val ao_basis : t -> Basis_poly.t -val overlap : t -> (ao,ao) Matrix.t -val multipole : t -> string -> (ao,ao) Matrix.t -val ortho : t -> (ao,'a) Matrix.t -val eN_ints : t -> (ao,ao) Matrix.t -val kin_ints : t -> (ao,ao) Matrix.t -val ee_ints : t -> ao Four_idx_storage.t -val ee_lr_ints : t -> ao Four_idx_storage.t -val f12_ints : t -> ao Four_idx_storage.t -val f12_over_r12_ints : t -> ao Four_idx_storage.t -val cartesian : t -> bool -val values : t -> Coordinate.t -> ao Vector.t - #+end_src - - |---------------------+--------------------------------------------------| - | ~size~ | Number of atomic orbitals in the AO basis set | - | ~ao_basis~ | One-electron basis set | - | ~overlap~ | Overlap matrix | - | ~multipole~ | Multipole matrices | - | ~ortho~ | Orthonormalization matrix of the overlap | - | ~eN_ints~ | Electron-nucleus potential integrals | - | ~kin_ints~ | Kinetic energy integrals | - | ~ee_ints~ | Electron-electron potential integrals | - | ~ee_lr_ints~ | Electron-electron long-range potential integrals | - | ~f12_ints~ | Electron-electron potential integrals | - | ~f12_over_r12_ints~ | Electron-electron potential integrals | - | ~cartesian~ | If true, use cartesian Gaussians (6d, 10f, ...) | - | ~values~ | Values of the AOs evaluated at a given point | - |---------------------+--------------------------------------------------| - - - #+begin_src ocaml :tangle (eval ml) :exports none -let ao_basis t = t.ao_basis - -let size t = - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.size b - | _ -> not_implemented () - -let overlap t = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.overlap b - | _ -> not_implemented () - end - |> Matrix.relabel - -let multipole t = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> - let m = Basis_gaussian.multipole b in - fun s -> - Gaussian_integrals.Multipole.matrix m s - |> Matrix.relabel - | _ -> not_implemented () - end - -let ortho t = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.ortho b - | _ -> not_implemented () - end - |> Matrix.relabel - -let eN_ints t = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.eN_ints b - | _ -> not_implemented () - end - |> Matrix.relabel - -let kin_ints t = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.kin_ints b - | _ -> not_implemented () - end - |> Matrix.relabel - -let ee_ints t = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.ee_ints b - | _ -> not_implemented () - end - |> Four_idx_storage.relabel - -let ee_lr_ints t = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.ee_lr_ints b - | _ -> not_implemented () - end - |> Four_idx_storage.relabel - -let f12_ints t = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.f12_ints b - | _ -> not_implemented () - end - |> Four_idx_storage.relabel - -let f12_over_r12_ints t = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.f12_over_r12_ints b - | _ -> not_implemented () - end - |> Four_idx_storage.relabel - -let cartesian t = t.cartesian - - -let values t point = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.values b point - | _ -> not_implemented () - end - |> Vector.relabel - - #+end_src - -** TREXIO - -*** Read - - #+begin_src ocaml :tangle (eval mli) -(* -val of_trexio : Trexio.trexio_file -> t -*) - #+end_src - - #+begin_src ocaml :tangle (eval ml) :exports none -(* -let of_trexio f = - - match Trexio.read_basis_type f with - | "G" - | "Gaussian" - | "gaussian" -> - let basis = - Gaussian.Basis.of_trexio f - in - let ao_basis = - Basis_poly.Gaussian (Basis_gaussian.make ~basis ?operators ~cartesian nuclei ) - in - { ao_basis ; cartesian } - | _ -> not_implemented () - -*) - #+end_src - -*** Write - - #+begin_src ocaml :tangle (eval mli) -(* -val to_trexio : Trexio.trexio_file -> t -> unit -*) - #+end_src - - #+begin_src ocaml :tangle (eval ml) :exports none -(* -let to_trexio f t = - - match t.ao_basis with - | Basis_poly.Gaussian b -> Gaussian.Basis.to_trexio f (Basis_gaussian.basis b)) - | _ -> not_implemented () - -*) - - #+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 = - begin - match t.ao_basis with - | Basis_poly.Gaussian b -> Basis_gaussian.pp ppf b - | _ -> not_implemented () - end - #+end_src - diff --git a/ao/basis_gaussian.org b/ao/basis_gaussian.org deleted file mode 100644 index 26bd7bc..0000000 --- a/ao/basis_gaussian.org +++ /dev/null @@ -1,232 +0,0 @@ -#+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 - -* Gaussian basis - :PROPERTIES: - :header-args: :noweb yes :comments both - :END: - - Data structure for Gaussian Atomic Orbitals: - - \[ - \chi_i(\mathbf{r}) = P_i(\mathbf{r}) \sum_k c_k \exp\left( -\alpha_k (\mathbf{r-R_A})^2 \right) - \] - - where the polynomial $P_i$ and the Gaussian part are both centered on - nucleus $A$. - -** Type - - <<<~Gaussian_basis.t~>>> - #+NAME: types - #+begin_src ocaml :tangle (eval mli) -open Common -open Particles -open Linear_algebra -open Gaussian_integrals -open Operators - -type t - #+end_src - - #+begin_src ocaml :tangle (eval ml) :exports none -open Linear_algebra -open Gaussian -open Gaussian_integrals -open Operators - -module Basis = Gaussian.Basis - -type t = -{ - basis : Basis.t ; - overlap : Overlap.t lazy_t; - multipole : Multipole.t lazy_t; - ortho : Orthonormalization.t lazy_t; - eN_ints : Electron_nucleus.t lazy_t; - kin_ints : Kinetic.t lazy_t; - ee_ints : Eri.t lazy_t; - ee_lr_ints : Eri_long_range.t lazy_t; - f12_ints : F12.t lazy_t; - f12_over_r12_ints : Screened_eri.t lazy_t; - cartesian : bool ; -} - #+end_src - -** Access - - #+begin_src ocaml :tangle (eval mli) -val basis : t -> Gaussian.Basis.t -val cartesian : t -> bool -val ee_ints : t -> Eri.t -val ee_lr_ints : t -> Eri_long_range.t -val eN_ints : t -> Electron_nucleus.t -val f12_ints : t -> F12.t -val f12_over_r12_ints : t -> Screened_eri.t -val kin_ints : t -> Kinetic.t -val multipole : t -> Multipole.t -val ortho : t -> Orthonormalization.t -val overlap : t -> Overlap.t -val size : t -> int - #+end_src - - |---------------------+--------------------------------------------------| - | ~basis~ | One-electron basis set | - | ~cartesian~ | If true, use cartesian Gaussians (6d, 10f, ...) | - | ~ee_ints~ | Electron-electron potential integrals | - | ~ee_lr_ints~ | Electron-electron long-range potential integrals | - | ~eN_ints~ | Electron-nucleus potential integrals | - | ~f12_ints~ | Electron-electron potential integrals | - | ~f12_over_r12_ints~ | Electron-electron potential integrals | - | ~kin_ints~ | Kinetic energy integrals | - | ~multipole~ | Multipole matrices | - | ~ortho~ | Orthonormalization matrix of the overlap | - | ~overlap~ | Overlap matrix | - | ~size~ | Number of atomic orbitals | - |---------------------+--------------------------------------------------| - - #+begin_src ocaml :tangle (eval ml) :exports none -let basis t = t.basis -let cartesian t = t.cartesian -let ee_ints t = Lazy.force t.ee_ints -let ee_lr_ints t = Lazy.force t.ee_lr_ints -let eN_ints t = Lazy.force t.eN_ints -let f12_ints t = Lazy.force t.f12_ints -let f12_over_r12_ints t = Lazy.force t.f12_over_r12_ints -let kin_ints t = Lazy.force t.kin_ints -let multipole t = Lazy.force t.multipole -let ortho t = Lazy.force t.ortho -let overlap t = Lazy.force t.overlap -let size t = Matrix.dim1 (Lazy.force t.overlap) - #+end_src - -** Computation - - #+begin_src ocaml :tangle (eval mli) -val values : t -> - Coordinate.t -> - Gaussian.Basis.t Vector.t - #+end_src - - |----------+--------------------------------------------------------------| - | ~values~ | Returns the values of all the AOs evaluated at a given point | - |----------+--------------------------------------------------------------| - - #+begin_src ocaml :tangle (eval ml) :exports none -module Cs = Contracted_shell - -let values t point = - let result = Vector.create (Basis.size t.basis) in - let resultx = Vector.to_bigarray_inplace result in - Array.iter (fun shell -> - Cs.values shell point - |> Array.iteri - (fun i_c value -> - let i = Cs.index shell + i_c + 1 in - resultx.{i} <- value) - ) (Basis.contracted_shells t.basis); - result - #+end_src - -** Creation - - #+begin_src ocaml :tangle (eval mli) -val make : basis:Gaussian.Basis.t -> - ?operators:Operator.t list -> - ?cartesian:bool -> - Nuclei.t -> - t - #+end_src - - Creates the data structure for atomic orbitals from a Gaussian basis and the - molecular geometry ~Nuclei.t~. - - Defaults: - - ~operators~ : ~[]~ - - ~cartesian~ : ~false~ - - Example: - #+begin_example -let b = Ao.Basis_gaussian.make ~basis nuclei ;; -val b : Ao.Basis_gaussian.t = Gaussian Basis, spherical, 15 AOs - #+end_example - - #+begin_src ocaml :tangle (eval ml) :exports none -let make ~basis ?(operators=[]) ?(cartesian=false) nuclei = - - let overlap = lazy ( - Overlap.of_basis basis - ) in - - let ortho = lazy ( - Lazy.force overlap - |> Orthonormalization.make ~cartesian ~basis - ) in - - let eN_ints = lazy ( - Electron_nucleus.of_basis_nuclei ~basis nuclei - ) in - - let kin_ints = lazy ( - Kinetic.of_basis basis - ) in - - let ee_ints = lazy ( - Eri.of_basis basis - ) in - - let rec get_f12 = function - | (Operator.F12 _ as f) :: _ -> f - | [] -> failwith "Missing F12 operator" - | _ :: rest -> get_f12 rest - in - - let rec get_rs = function - | (Operator.Range_sep _ as r) :: _ -> r - | [] -> failwith "Missing range-separation operator" - | _ :: rest -> get_rs rest - in - - let ee_lr_ints = lazy ( - Eri_long_range.of_basis basis~operator:(get_rs operators) - ) in - - let f12_ints = lazy ( - F12.of_basis basis ~operator:(get_f12 operators) - ) in - - let f12_over_r12_ints = lazy ( - Screened_eri.of_basis basis ~operator:(get_f12 operators) - ) in - - let multipole = lazy ( - Multipole.of_basis basis - ) in - - { basis ; overlap ; multipole ; ortho ; eN_ints ; kin_ints ; ee_ints ; - ee_lr_ints ; f12_ints ; f12_over_r12_ints ; cartesian } - #+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 cart = if t.cartesian then "cartesian" else "spherical" in - let nao = size t in - Format.fprintf ppf "@[@[Gaussian Basis@], @[%s@], @[%d AOs@]@]" - cart nao - #+end_src - diff --git a/ao/lib/basis.ml b/ao/lib/basis.ml index 8b4ad48..d1faf93 100644 --- a/ao/lib/basis.ml +++ b/ao/lib/basis.ml @@ -1,4 +1,3 @@ -(* [[file:~/QCaml/ao/basis.org::*Types][Types:2]] *) type t = { ao_basis : Basis_poly.t ; cartesian : bool @@ -8,30 +7,13 @@ type ao = Ao_dim.t open Linear_algebra open Common -(* Types:2 ends here *) -(* |--------------------------------+------------------------------------------------------------------------------------------------------------------------| - * | ~of_nuclei_and_basis_filename~ | Creates the data structure for the atomic orbitals basis from a molecule ~Nuclei.t~ and the name of the basis-set file | - * |--------------------------------+------------------------------------------------------------------------------------------------------------------------| - * - * Defaults: - * - ~kind~ : ~`Gaussian~ - * - ~operators~ : ~[]~ - * - ~cartesian~ : ~false~ - * - * Example: - * #+begin_example - * let b = Ao.Basis.of_nuclei_and_basis_filename ~nuclei filename;; - * val b : Ao.Basis.t = Gaussian Basis, spherical, 15 AOs - * #+end_example *) - - -(* [[file:~/QCaml/ao/basis.org::*Conversions][Conversions:2]] *) let not_implemented () = Util.not_implemented "Only Gaussian is implemented" + let of_nuclei_and_basis_filename ?(kind=`Gaussian) ?operators ?(cartesian=false) ~nuclei filename = match kind with @@ -44,29 +26,10 @@ let of_nuclei_and_basis_filename ?(kind=`Gaussian) ?operators ?(cartesian=false) in { ao_basis ; cartesian } | _ -> not_implemented () -(* Conversions:2 ends here *) -(* |---------------------+--------------------------------------------------| - * | ~size~ | Number of atomic orbitals in the AO basis set | - * | ~ao_basis~ | One-electron basis set | - * | ~overlap~ | Overlap matrix | - * | ~multipole~ | Multipole matrices | - * | ~ortho~ | Orthonormalization matrix of the overlap | - * | ~eN_ints~ | Electron-nucleus potential integrals | - * | ~kin_ints~ | Kinetic energy integrals | - * | ~ee_ints~ | Electron-electron potential integrals | - * | ~ee_lr_ints~ | Electron-electron long-range potential integrals | - * | ~f12_ints~ | Electron-electron potential integrals | - * | ~f12_over_r12_ints~ | Electron-electron potential integrals | - * | ~cartesian~ | If true, use cartesian Gaussians (6d, 10f, ...) | - * | ~values~ | Values of the AOs evaluated at a given point | - * |---------------------+--------------------------------------------------| *) - - -(* [[file:~/QCaml/ao/basis.org::*Access][Access:2]] *) let ao_basis t = t.ao_basis let size t = @@ -159,9 +122,8 @@ let values t point = | _ -> not_implemented () end |> Vector.relabel -(* Access:2 ends here *) -(* [[file:~/QCaml/ao/basis.org::*Read][Read:2]] *) + (* let of_trexio f = @@ -178,11 +140,8 @@ let of_trexio f = { ao_basis ; cartesian } | _ -> not_implemented () -*) -(* Read:2 ends here *) -(* [[file:~/QCaml/ao/basis.org::*Write][Write:2]] *) -(* + let to_trexio f t = match t.ao_basis with @@ -190,13 +149,12 @@ let to_trexio f t = | _ -> not_implemented () *) -(* Write:2 ends here *) -(* [[file:~/QCaml/ao/basis.org::*Printers][Printers:2]] *) + + let pp ppf t = begin match t.ao_basis with | Basis_poly.Gaussian b -> Basis_gaussian.pp ppf b | _ -> not_implemented () end -(* Printers:2 ends here *) diff --git a/ao/lib/basis.mli b/ao/lib/basis.mli index ae2c37c..c632d5e 100644 --- a/ao/lib/basis.mli +++ b/ao/lib/basis.mli @@ -1,8 +1,5 @@ -(* Types - * - * <<<~Basis.t~>>> *) +(** Data structure for Atomic Orbitals. *) -(* [[file:~/QCaml/ao/basis.org::*Types][Types:1]] *) type t type ao = Ao_dim.t @@ -10,61 +7,85 @@ open Common open Particles open Operators open Linear_algebra -(* Types:1 ends here *) -(* Conversions *) +(** Conversions *) - -(* [[file:~/QCaml/ao/basis.org::*Conversions][Conversions:1]] *) val of_nuclei_and_basis_filename : ?kind:[> `Gaussian ] -> ?operators:Operator.t list -> ?cartesian:bool -> nuclei:Nuclei.t -> - string -> - t -(* Conversions:1 ends here *) + string -> t +(** Creates the data structure for the atomic orbitals basis from a + molecule ~Nuclei.t~ and the name of the basis-set file. -(* Access *) + Defaults: + - ~kind~ : ~`Gaussian~ + - ~operators~ : ~[]~ + - ~cartesian~ : ~false~ + + Example: + #+begin_example +let b = Ao.Basis.of_nuclei_and_basis_filename ~nuclei filename;; +val b : Ao.Basis.t = Gaussian Basis, spherical, 15 AOs + #+end_example + +*) -(* [[file:~/QCaml/ao/basis.org::*Access][Access:1]] *) -val size : t -> int -val ao_basis : t -> Basis_poly.t -val overlap : t -> (ao,ao) Matrix.t -val multipole : t -> string -> (ao,ao) Matrix.t -val ortho : t -> (ao,'a) Matrix.t -val eN_ints : t -> (ao,ao) Matrix.t -val kin_ints : t -> (ao,ao) Matrix.t -val ee_ints : t -> ao Four_idx_storage.t -val ee_lr_ints : t -> ao Four_idx_storage.t -val f12_ints : t -> ao Four_idx_storage.t +(** Access *) + +val size : t -> int +(** Number of atomic orbitals in the AO basis set *) + +val ao_basis : t -> Basis_poly.t +(** One-electron basis set *) + +val overlap : t -> (ao,ao) Matrix.t +(** Overlap matrix *) + +val multipole : t -> string -> (ao,ao) Matrix.t +(** Multipole matrices *) + +val ortho : t -> (ao,'a) Matrix.t +(** Orthonormalization matrix of the overlap *) + +val eN_ints : t -> (ao,ao) Matrix.t +(** Electron-nucleus potential integrals *) + +val kin_ints : t -> (ao,ao) Matrix.t +(** Kinetic energy integrals *) + +val ee_ints : t -> ao Four_idx_storage.t +(** Electron-electron potential integrals *) + +val ee_lr_ints : t -> ao Four_idx_storage.t +(** Electron-electron long-range potential integrals *) + +val f12_ints : t -> ao Four_idx_storage.t +(** Electron-electron potential integrals *) + val f12_over_r12_ints : t -> ao Four_idx_storage.t -val cartesian : t -> bool -val values : t -> Coordinate.t -> ao Vector.t -(* Access:1 ends here *) +(** ectron-electron potential integrals *) -(* Read *) +val cartesian : t -> bool +(** If true, use cartesian Gaussians (6d, 10f, ...) *) + +val values : t -> Coordinate.t -> ao Vector.t +(** Values of the AOs evaluated at a given point *) -(* [[file:~/QCaml/ao/basis.org::*Read][Read:1]] *) +(** TREXIO *) + (* val of_trexio : Trexio.trexio_file -> t -*) -(* Read:1 ends here *) +(** Reads the basis set from a TREXIO file *) -(* Write *) - - -(* [[file:~/QCaml/ao/basis.org::*Write][Write:1]] *) -(* val to_trexio : Trexio.trexio_file -> t -> unit +(** Writes the basis set to a TREXIO file *) *) -(* Write:1 ends here *) + (* Printers *) - -(* [[file:~/QCaml/ao/basis.org::*Printers][Printers:1]] *) val pp : Format.formatter -> t -> unit -(* Printers:1 ends here *) diff --git a/ao/lib/basis_gaussian.ml b/ao/lib/basis_gaussian.ml index 8095509..47ee3a9 100644 --- a/ao/lib/basis_gaussian.ml +++ b/ao/lib/basis_gaussian.ml @@ -1,4 +1,3 @@ -(* [[file:~/QCaml/ao/basis_gaussian.org::*Type][Type:2]] *) open Linear_algebra open Gaussian open Gaussian_integrals @@ -20,27 +19,7 @@ type t = f12_over_r12_ints : Screened_eri.t lazy_t; cartesian : bool ; } -(* Type:2 ends here *) - - -(* |---------------------+--------------------------------------------------| - * | ~basis~ | One-electron basis set | - * | ~cartesian~ | If true, use cartesian Gaussians (6d, 10f, ...) | - * | ~ee_ints~ | Electron-electron potential integrals | - * | ~ee_lr_ints~ | Electron-electron long-range potential integrals | - * | ~eN_ints~ | Electron-nucleus potential integrals | - * | ~f12_ints~ | Electron-electron potential integrals | - * | ~f12_over_r12_ints~ | Electron-electron potential integrals | - * | ~kin_ints~ | Kinetic energy integrals | - * | ~multipole~ | Multipole matrices | - * | ~ortho~ | Orthonormalization matrix of the overlap | - * | ~overlap~ | Overlap matrix | - * | ~size~ | Number of atomic orbitals | - * |---------------------+--------------------------------------------------| *) - - -(* [[file:~/QCaml/ao/basis_gaussian.org::*Access][Access:2]] *) let basis t = t.basis let cartesian t = t.cartesian let ee_ints t = Lazy.force t.ee_ints @@ -53,7 +32,6 @@ let multipole t = Lazy.force t.multipole let ortho t = Lazy.force t.ortho let overlap t = Lazy.force t.overlap let size t = Matrix.dim1 (Lazy.force t.overlap) -(* Access:2 ends here *) @@ -62,7 +40,6 @@ let size t = Matrix.dim1 (Lazy.force t.overlap) * |----------+--------------------------------------------------------------| *) -(* [[file:~/QCaml/ao/basis_gaussian.org::*Computation][Computation:2]] *) module Cs = Contracted_shell let values t point = @@ -76,25 +53,10 @@ let values t point = resultx.{i} <- value) ) (Basis.contracted_shells t.basis); result -(* Computation:2 ends here *) -(* Creates the data structure for atomic orbitals from a Gaussian basis and the - * molecular geometry ~Nuclei.t~. - * - * Defaults: - * - ~operators~ : ~[]~ - * - ~cartesian~ : ~false~ - * - * Example: - * #+begin_example - * let b = Ao.Basis_gaussian.make ~basis nuclei ;; - * val b : Ao.Basis_gaussian.t = Gaussian Basis, spherical, 15 AOs - * #+end_example *) - -(* [[file:~/QCaml/ao/basis_gaussian.org::*Creation][Creation:2]] *) let make ~basis ?(operators=[]) ?(cartesian=false) nuclei = let overlap = lazy ( @@ -148,12 +110,11 @@ let make ~basis ?(operators=[]) ?(cartesian=false) nuclei = { basis ; overlap ; multipole ; ortho ; eN_ints ; kin_ints ; ee_ints ; ee_lr_ints ; f12_ints ; f12_over_r12_ints ; cartesian } -(* Creation:2 ends here *) -(* [[file:~/QCaml/ao/basis_gaussian.org::*Printers][Printers:2]] *) + + let pp ppf t = let cart = if t.cartesian then "cartesian" else "spherical" in let nao = size t in Format.fprintf ppf "@[@[Gaussian Basis@], @[%s@], @[%d AOs@]@]" cart nao -(* Printers:2 ends here *) diff --git a/ao/lib/basis_gaussian.mli b/ao/lib/basis_gaussian.mli index 7c177bd..99aa10a 100644 --- a/ao/lib/basis_gaussian.mli +++ b/ao/lib/basis_gaussian.mli @@ -1,9 +1,16 @@ -(* Type - * - * <<<~Gaussian_basis.t~>>> - * #+NAME: types *) +(** Gaussian basis + + Data structure for Gaussian Atomic Orbitals: + + \[ + \chi_i(\mathbf{r}) = P_i(\mathbf{r}) \sum_k c_k \exp\left( -\alpha_k (\mathbf{r-R_A})^2 \r$ + \] + + where the polynomial $P_i$ and the Gaussian part are both centered on + nucleus $A$. + +*) -(* [[file:~/QCaml/ao/basis_gaussian.org::types][types]] *) open Common open Particles open Linear_algebra @@ -11,49 +18,74 @@ open Gaussian_integrals open Operators type t -(* types ends here *) -(* Access *) +(** Access *) +val basis : t -> Gaussian.Basis.t +(** One-electron basis set *) + +val cartesian : t -> bool +(** If true, use cartesian Gaussians (6d, 10f, ...) *) + +val ee_ints : t -> Eri.t +(* Electron-electron potential integrals *) + +val ee_lr_ints : t -> Eri_long_range.t +(** Electron-electron long-range potential integrals *) + +val eN_ints : t -> Electron_nucleus.t +(** Electron-nucleus potential integrals *) + +val f12_ints : t -> F12.t +(** Electron-electron potential integrals *) -(* [[file:~/QCaml/ao/basis_gaussian.org::*Access][Access:1]] *) -val basis : t -> Gaussian.Basis.t -val cartesian : t -> bool -val ee_ints : t -> Eri.t -val ee_lr_ints : t -> Eri_long_range.t -val eN_ints : t -> Electron_nucleus.t -val f12_ints : t -> F12.t val f12_over_r12_ints : t -> Screened_eri.t -val kin_ints : t -> Kinetic.t -val multipole : t -> Multipole.t -val ortho : t -> Orthonormalization.t +(**lectron-electron potential integrals *) + +val kin_ints : t -> Kinetic.t +(** Kinetic energy integrals *) + +val multipole : t -> Multipole.t +(** Multipole matrices *) + +val ortho : t -> Orthonormalization.t +(** Orthonormalization matrix of the overlap *) + val overlap : t -> Overlap.t -val size : t -> int -(* Access:1 ends here *) +(** Overlap matrix *) -(* Computation *) +val size : t -> int +(** Number of atomic orbitals *) -(* [[file:~/QCaml/ao/basis_gaussian.org::*Computation][Computation:1]] *) -val values : t -> - Coordinate.t -> - Gaussian.Basis.t Vector.t -(* Computation:1 ends here *) +(** Computation *) -(* Creation *) +val values : t -> Coordinate.t -> Gaussian.Basis.t Vector.t +(** Returns the values of all the AOs evaluated at a given point *) -(* [[file:~/QCaml/ao/basis_gaussian.org::*Creation][Creation:1]] *) +(** Creation *) + val make : basis:Gaussian.Basis.t -> ?operators:Operator.t list -> ?cartesian:bool -> - Nuclei.t -> - t -(* Creation:1 ends here *) + Nuclei.t -> t -(* Printers *) +(** Creates the data structure for atomic orbitals from a Gaussian basis and the + molecular geometry ~Nuclei.t~. + + Defaults: + - ~operators~ : ~[]~ + - ~cartesian~ : ~false~ + + Example: + #+begin_example +let b = Ao.Basis_gaussian.make ~basis nuclei ;; +val b : Ao.Basis_gaussian.t = Gaussian Basis, spherical, 15 AOs + #+end_example +*) -(* [[file:~/QCaml/ao/basis_gaussian.org::*Printers][Printers:1]] *) +(** Printers *) + val pp : Format.formatter -> t -> unit -(* Printers:1 ends here *) diff --git a/bin/build_doc.sh b/bin/build_doc.sh deleted file mode 100755 index 3c2a531..0000000 --- a/bin/build_doc.sh +++ /dev/null @@ -1,34 +0,0 @@ -#!/bin/bash -# Usage: $0 [DIR] - -if [[ -z $1 ]] ; then - echo "Usage: $0 [DIR]" - exit -1 -fi - - -if [[ $(basename $PWD) != "QCaml" ]] ; then - echo "This script needs to be run in the QCaml directory" - exit -1 -fi - -DIR=${1%/} - -CONFIG="--load docs/htmlize.el --load docs/config.el" -if [[ ${DIR} == "docs" ]] ; then - - emacs --debug-init --batch $CONFIG docs/index.org -f org-html-export-to-html - -else - rm -f docs/${DIR}.org - for i in ${DIR}/README.org ${DIR}/[a-z]*.org - do - cat $i >> docs/${DIR}.org - done - - emacs --debug-init --batch $CONFIG docs/${DIR}.org -f org-html-export-to-html \ - && rm docs/${DIR}.org -fi - - - diff --git a/common/angular_momentum.org b/common/angular_momentum.org deleted file mode 100644 index 2b7313e..0000000 --- a/common/angular_momentum.org +++ /dev/null @@ -1,288 +0,0 @@ -#+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 - - -* Angular Momentum - :PROPERTIES: - :header-args: :noweb yes :comments both - :END: - - Azimuthal quantum number, repsesented as \( s,p,d,\dots \) . - -** Type - - <<<~Angular_momentum.t~>>> - #+NAME: types - #+begin_src ocaml :tangle (eval mli) -type t = - | S | P | D | F | G | H | I | J | K | L | M | N | O - | Int of int - -exception Angular_momentum_error of string - -type kind = - Singlet of t - | Doublet of (t * t) - | Triplet of (t * t * t) - | Quartet of (t * t * t * t) - - #+end_src - - An exception is raised when the ~Angular_momentum.t~ element can't - be created. - - The ~kind~ is used to build shells, shell doublets, triplets or - quartets, use in the two-electron operators. - - #+begin_src ocaml :tangle (eval ml) :exports none -<> -open Powers - #+end_src - -** Conversions - - #+begin_src ocaml :tangle (eval mli) -val of_char : char -> t -val to_char : t -> char - -val to_int : t -> int -val of_int : int -> t - -val to_string : t -> string - #+end_src - - | ~of_char~ | Returns an ~Angular_momentum.t~ when a shell is given as a character (case insensitive) | - | ~to_char~ | Converts the angular momentum into a char | - | ~of_int~ | Returns a shell given an $l$ value. | - | ~to_int~ | Returns the $l_{max}$ value of the shell | - | ~to_string~ | Converts the angular momentum into a string | - - Example: - #+begin_example -Angular_momentum.of_char 'p';; -- : Angular_momentum.t = P - -Angular_momentum.(to_char P);; -- : char = 'P' - -Angular_momentum.of_int 2;; -- : Angular_momentum.t = D - -Angular_momentum.(to_int D);; -- : int = 2 - -Angular_momentum.(to_string D);; -- : string = "D" - #+end_example - - #+begin_src ocaml :tangle (eval ml) :exports none -let of_char = function - | 's' | 'S' -> S | 'p' | 'P' -> P - | 'd' | 'D' -> D | 'f' | 'F' -> F - | 'g' | 'G' -> G | 'h' | 'H' -> H - | 'i' | 'I' -> I | 'j' | 'J' -> J - | 'k' | 'K' -> K | 'l' | 'L' -> L - | 'm' | 'M' -> M | 'n' | 'N' -> N - | 'o' | 'O' -> O - | c -> raise (Angular_momentum_error (String.make 1 c)) - -let to_string = function - | S -> "S" | P -> "P" - | D -> "D" | F -> "F" - | G -> "G" | H -> "H" - | I -> "I" | J -> "J" - | K -> "K" | L -> "L" - | M -> "M" | N -> "N" - | O -> "O" | Int i -> string_of_int i - - -let to_char = function - | S -> 'S' | P -> 'P' - | D -> 'D' | F -> 'F' - | G -> 'G' | H -> 'H' - | I -> 'I' | J -> 'J' - | K -> 'K' | L -> 'L' - | M -> 'M' | N -> 'N' - | O -> 'O' | Int _ -> '_' - - -let to_int = function - | S -> 0 | P -> 1 - | D -> 2 | F -> 3 - | G -> 4 | H -> 5 - | I -> 6 | J -> 7 - | K -> 8 | L -> 9 - | M -> 10 | N -> 11 - | O -> 12 | Int i -> i - - -let of_int = function - | 0 -> S | 1 -> P - | 2 -> D | 3 -> F - | 4 -> G | 5 -> H - | 6 -> I | 7 -> J - | 8 -> K | 9 -> L - | 10 -> M | 11 -> N - | 12 -> O | i -> Int i - #+end_src - -** Shell functions - - #+begin_src ocaml :tangle (eval mli) -val n_functions : t -> int -val zkey_array : kind -> Zkey.t array - #+end_src - - | ~n_functions~ | Returns the number of cartesian functions in a shell. | - | ~zkey_array~ | Array of ~Zkey.t~, where each element is a a key associated with the the powers of $x,y,z$. | - - Example: - #+begin_example -Angular_momentum.(n_functions D) ;; -- : int = 6 - -Angular_momentum.( zkey_array (Doublet (P,S)) );; -- : Zkey.t array = -[|< 01125899906842624 | 1, 0, 0, 0, 0, 0 >; - < 01099511627776 | 0, 1, 0, 0, 0, 0 >; - < 01073741824 | 0, 0, 1, 0, 0, 0 >|] - #+end_example - - #+begin_src ocaml :tangle (eval ml) :exports none -let n_functions a = - let a = - to_int a - in - (a*a + 3*a + 2)/2 - - -let zkey_array_memo : (kind, Zkey.t array) Hashtbl.t = - Hashtbl.create 13 - -let zkey_array a = - - let keys_1d l = - let create_z { x ; y ; _ } = - Powers.of_int_tuple (x,y,l-(x+y)) - in - let rec create_y accu xyz = - let { x ; y ; z ;_ } = xyz in - match y with - | 0 -> (create_z xyz)::accu - | _ -> let ynew = y-1 in - (create_y [@tailcall]) ( (create_z xyz)::accu) (Powers.of_int_tuple (x,ynew,z)) - in - let rec create_x accu xyz = - let { x ; z ;_ } = xyz in - match x with - | 0 -> (create_y [] xyz)@accu - | _ -> let xnew = x-1 in - let ynew = l-xnew in - (create_x [@tailcall]) ((create_y [] xyz)@accu) (Powers.of_int_tuple (xnew, ynew, z)) - in - create_x [] (Powers.of_int_tuple (l,0,0)) - |> List.rev - in - - try - Hashtbl.find zkey_array_memo a - - with Not_found -> - - let result = - begin - match a with - | Singlet l1 -> - List.rev_map (fun x -> Zkey.of_powers_three x) (keys_1d @@ to_int l1) - - | Doublet (l1, l2) -> - List.rev_map (fun a -> - List.rev_map (fun b -> Zkey.of_powers_six a b) (keys_1d @@ to_int l2) - ) (keys_1d @@ to_int l1) - |> List.concat - - | Triplet (l1, l2, l3) -> - - List.rev_map (fun a -> - List.rev_map (fun b -> - List.rev_map (fun c -> - Zkey.of_powers_nine a b c) (keys_1d @@ to_int l3) - ) (keys_1d @@ to_int l2) - |> List.concat - ) (keys_1d @@ to_int l1) - |> List.concat - - | Quartet (l1, l2, l3, l4) -> - - List.rev_map (fun a -> - List.rev_map (fun b -> - List.rev_map (fun c -> - List.rev_map (fun d -> - Zkey.of_powers_twelve a b c d) (keys_1d @@ to_int l4) - ) (keys_1d @@ to_int l3) - |> List.concat - ) (keys_1d @@ to_int l2) - |> List.concat - ) (keys_1d @@ to_int l1) - |> List.concat - end - |> List.rev - |> Array.of_list - in - Hashtbl.add zkey_array_memo a result; - result - #+end_src - -** Arithmetic - - #+begin_src ocaml :tangle (eval mli) -val ( + ) : t -> t -> t -val ( - ) : t -> t -> t - #+end_src - - Example: - #+begin_example -Angular_momentum.(D + P);; -- : Angular_momentum.t = F - -Angular_momentum.(F - P);; -- : Angular_momentum.t = D - #+end_example - - #+begin_src ocaml :tangle (eval ml) :exports none -let ( + ) a b = - of_int ( (to_int a) + (to_int b) ) - -let ( - ) a b = - of_int ( (to_int a) - (to_int b) ) - #+end_src - -** Printers - - Printers can print as a string (default) or as an integer. - - #+begin_src ocaml :tangle (eval mli) -val pp : Format.formatter -> t -> unit -val pp_string : Format.formatter -> t -> unit -val pp_int : Format.formatter -> t -> unit - #+end_src - - #+begin_src ocaml :tangle (eval ml) :exports none -let pp_string ppf x = - Format.fprintf ppf "@[%s@]" (to_string x) - -let pp_int ppf x = - Format.fprintf ppf "@[%d@]" (to_int x) - -let pp = pp_string - #+end_src - -** TODO Tests diff --git a/common/bitstring.org b/common/bitstring.org deleted file mode 100644 index 52f5b02..0000000 --- a/common/bitstring.org +++ /dev/null @@ -1,460 +0,0 @@ -#+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 - -* Bit string - :PROPERTIES: - :header-args: :noweb yes :comments both - :END: - - We define here a data type to handle bit strings efficiently. When - the bit string contains less than 64 bits, it is stored internally - in a 63-bit integer and uses bitwise instructions. When more than - 63 bits are required, the =zarith= library is used to consider the - bit string as a multi-precision integer. - - -** Single-integer implementation :noexport: - - #+begin_src ocaml :tangle (eval ml) :exports none -module One = struct - - let of_int x = - assert (x > 0); x - - let numbits _ = 63 - let zero = 0 - let is_zero x = x = 0 - let shift_left x i = x lsl i - let shift_right x i = x lsr i - let shift_left_one i = 1 lsl i - let testbit x i = ( (x lsr i) land 1 ) = 1 - let logor a b = a lor b - let neg a = - a - let logxor a b = a lxor b - let logand a b = a land b - let lognot a = lnot a - let minus_one a = a - 1 - let plus_one a = a + 1 - - let popcount = function - | 0 -> 0 - | r -> Util.popcnt (Int64.of_int r) - - let trailing_zeros r = - Util.trailz (Int64.of_int r) - - let hamdist a b = - a lxor b - |> popcount - - let pp ppf s = - Format.fprintf ppf "@[@[%a@]@]" (Util.pp_bitstring 64) - (Z.of_int s) - -end - #+end_src - -** Zarith implementation :noexport: - - #+begin_src ocaml :tangle (eval ml) :exports none -module Many = struct - - let of_z x = x - let zero = Z.zero - let is_zero x = x = Z.zero - let shift_left = Z.shift_left - let shift_right = Z.shift_right - let shift_left_one i = Z.shift_left Z.one i - let testbit = Z.testbit - let logor = Z.logor - let logxor = Z.logxor - let logand = Z.logand - let lognot = Z.lognot - let neg = Z.neg - let minus_one = Z.pred - let plus_one = Z.succ - let trailing_zeros = Z.trailing_zeros - let hamdist = Z.hamdist - let numbits i = max (Z.numbits i) 64 - - let popcount z = - if z = Z.zero then 0 else Z.popcount z - - let pp ppf s = - Format.fprintf ppf "@[@[%a@]@]" (Util.pp_bitstring (Z.numbits s)) s - -end - #+end_src - -** Type - - <<<~Bitstring.t~>>> - #+begin_src ocaml :tangle (eval mli) -type t - #+end_src - - #+begin_src ocaml :tangle (eval ml) :exports none -type t = - | One of int - | Many of Z.t - #+end_src - -** Tests header :noexport: - - #+begin_src ocaml :tangle (eval test-ml) -open Common.Bitstring -let check_bool = Alcotest.(check bool) -let check msg x = check_bool msg true x -let test_all () = - let x = 8745687 in - let one_x = of_int x in - let z = Z.shift_left (Z.of_int x) 64 in - let many_x = of_z z in - #+end_src - -** General implementation - - #+begin_src ocaml :tangle (eval mli) -val of_int : int -> t -val of_z : Z.t -> t -val zero : int -> t - -val is_zero : t -> bool -val numbits : t -> int -val testbit : t -> int -> bool - -val neg : t -> t -val shift_left : t -> int -> t -val shift_right : t -> int -> t -val shift_left_one : int -> int -> t - -val logor : t -> t -> t -val logxor : t -> t -> t -val logand : t -> t -> t -val lognot : t -> t - -val plus_one : t -> t -val minus_one : t -> t - -val hamdist : t -> t -> int -val trailing_zeros : t -> int -val popcount : t -> int - -val to_list : ?accu:(int list) -> t -> int list -val permutations : int -> int -> t list - #+end_src - - | ~of_int~ | Creates a bit string from an ~int~ | - | ~of_z~ | Creates a bit string from an ~Z.t~ multi-precision integer | - | ~zero~ | ~zero n~ creates a zero bit string with ~n~ bits | - | ~is_zero~ | True if all the bits of the bit string are zero. | - | ~numbits~ | Returns the number of bits used to represent the bit string | - | ~testbit~ | ~testbit t n~ is true if the ~n~-th bit of the bit string ~t~ is set to ~1~ | - | ~neg~ | Returns the negative of the integer interpretation of the bit string | - | ~shift_left~ | ~shift_left t n~ returns a new bit strings with all the bits shifted ~n~ positions to the left | - | ~shift_right~ | ~shift_right t n~ returns a new bit strings with all the bits shifted ~n~ positions to the right | - | ~shift_left_one~ | ~shift_left_one size n~ returns a new bit strings with the ~n~-th bit set to one. It is equivalent as shifting ~1~ by ~n~ bits to the left, ~size~ is the total number of bits of the bit string | - | ~logor~ | Bitwise logical or | - | ~logxor~ | Bitwise logical exclusive or | - | ~logand~ | Bitwise logical and | - | ~lognot~ | Bitwise logical negation | - | ~plus_one~ | Takes the integer representation of the bit string and adds one | - | ~minus_one~ | Takes the integer representation of the bit string and removes one | - | ~hamdist~ | Returns the Hamming distance, i.e. the number of bits differing between two bit strings | - | ~trailing_zeros~ | Returns the number of trailing zeros in the bit string | - | ~permutations~ | ~permutations m n~ generates the list of all possible ~n~-bit strings with ~m~ bits set to ~1~. Algorithm adapted from [[https://graphics.stanford.edu/~seander/bithacks.html#NextBitPermutation][Bit twiddling hacks]] | - | ~popcount~ | Returns the number of bits set to one in the bit string | - | ~to_list~ | Converts a bit string into a list of integers indicating the positions where the bits are set to ~1~. The first value for the position is not ~0~ but ~1~ | - - #+begin_src ocaml :tangle (eval ml) :exports none -let of_int x = - One (One.of_int x) - -let of_z x = - if Z.numbits x < 64 then One (Z.to_int x) else Many (Many.of_z x) - #+end_src - - #+begin_src ocaml :tangle (eval test-ml) :exports none -check_bool "of_x" true (one_x = (of_int x)); -check_bool "of_z" true (one_x = (of_z (Z.of_int x))); - #+end_src - - Example: - #+begin_example -Bitstring.of_int 15;; -- : Bitstring.t = -++++------------------------------------------------------------ - #+end_example - - #+begin_src ocaml :tangle (eval ml) :exports none -let zero = function - | n when n < 64 -> One (One.zero) - | _ -> Many (Many.zero) - - -let numbits = function - | One x -> One.numbits x - | Many x -> Many.numbits x - - -let is_zero = function - | One x -> One.is_zero x - | Many x -> Many.is_zero x - - -let neg = function - | One x -> One (One.neg x) - | Many x -> Many (Many.neg x) - - -let shift_left x i = match x with - | One x -> One (One.shift_left x i) - | Many x -> Many (Many.shift_left x i) - - -let shift_right x i = match x with - | One x -> One (One.shift_right x i) - | Many x -> Many (Many.shift_right x i) - -let shift_left_one = function - | n when n < 64 -> fun i -> One (One.shift_left_one i) - | _ -> fun i -> Many (Many.shift_left_one i) - -let testbit = function - | One x -> One.testbit x - | Many x -> Many.testbit x - #+end_src - - Example: - #+begin_example -Bitstring.(shift_left (of_int 15) 2);; -- : Bitstring.t = ---++++---------------------------------------------------------- - -Bitstring.(shift_right (of_int 15) 2);; -- : Bitstring.t = -++-------------------------------------------------------------- - -Bitstring.shift_left_one 32 4;; -- : Bitstring.t = -----+----------------------------------------------------------- - -Bitstring.(testbit (of_int 15) 3);; -- : bool = true - -Bitstring.(testbit (of_int 15) 4);; -- : bool = false - #+end_example - - #+begin_src ocaml :tangle (eval test-ml) :exports none -check_bool "shift_left1" true (of_int (x lsl 3) = shift_left one_x 3); -check_bool "shift_left2" true (of_z (Z.shift_left z 3) = shift_left many_x 3); -check_bool "shift_left3" true (of_z (Z.shift_left z 100) = shift_left many_x 100); -check_bool "shift_right1" true (of_int (x lsr 3) = shift_right one_x 3); -check_bool "shift_right2" true (of_z (Z.shift_right z 3) = shift_right many_x 3); -check_bool "shift_left_one1" true (of_int (1 lsl 3) = shift_left_one 4 3); -check_bool "shift_left_one2" true (of_z (Z.shift_left Z.one 200) = shift_left_one 300 200); -check_bool "testbit1" true (testbit (of_int 8) 3); -check_bool "testbit2" false (testbit (of_int 8) 2); -check_bool "testbit3" false (testbit (of_int 8) 4); -check_bool "testbit4" true (testbit (of_z (Z.of_int 8)) 3); -check_bool "testbit5" false (testbit (of_z (Z.of_int 8)) 2); -check_bool "testbit6" false (testbit (of_z (Z.of_int 8)) 4); - #+end_src - - - #+begin_src ocaml :tangle (eval ml) :exports none -let logor a b = - match a,b with - | One a, One b -> One (One.logor a b) - | Many a, Many b -> Many (Many.logor a b) - | _ -> invalid_arg "Bitstring.logor" - - -let logxor a b = - match a,b with - | One a, One b -> One (One.logxor a b) - | Many a, Many b -> Many (Many.logxor a b) - | _ -> invalid_arg "Bitstring.logxor" - - -let logand a b = - match a,b with - | One a, One b -> One (One.logand a b) - | Many a, Many b -> Many (Many.logand a b) - | _ -> invalid_arg "Bitstring.logand" - - -let lognot = function - | One x -> One (One.lognot x) - | Many x -> Many (Many.lognot x) - #+end_src - - - Example: - #+begin_example -Bitstring.(logor (of_int 15) (of_int 73));; -- : Bitstring.t = -++++--+--------------------------------------------------------- - -Bitstring.(logand (of_int 15) (of_int 10));; -- : Bitstring.t = --+-+------------------------------------------------------------ - -Bitstring.(logxor (of_int 15) (of_int 73));; -- : Bitstring.t = --++---+--------------------------------------------------------- - #+end_example - - - #+begin_src ocaml :tangle (eval test-ml) :exports none -check_bool "logor1" true (of_int (1 lor 2) = logor (of_int 1) (of_int 2)); -check_bool "logor2" true (of_z (Z.of_int (1 lor 2)) = logor (of_z Z.one) (of_z (Z.of_int 2))); -check_bool "logxor1" true (of_int (1 lxor 2) = logxor (of_int 1) (of_int 2)); -check_bool "logxor2" true (of_z (Z.of_int (1 lxor 2)) = logxor (of_z Z.one) (of_z (Z.of_int 2))); -check_bool "logand1" true (of_int (1 land 3) = logand (of_int 1) (of_int 3)); -check_bool "logand2" true (of_z (Z.of_int (1 land 3)) = logand (of_z Z.one) (of_z (Z.of_int 3))); - #+end_src - - - #+begin_src ocaml :tangle (eval ml) :exports none -let minus_one = function - | One x -> One (One.minus_one x) - | Many x -> Many (Many.minus_one x) - - -let plus_one = function - | One x -> One (One.plus_one x) - | Many x -> Many (Many.plus_one x) - #+end_src - - - Example: - #+begin_example -Bitstring.(plus_one (of_int 15));; -- : Bitstring.t = -----+----------------------------------------------------------- - -Bitstring.(minus_one (of_int 15));; -- : Bitstring.t = --+++------------------------------------------------------------ - #+end_example - - - #+begin_src ocaml :tangle (eval ml) :exports none -let trailing_zeros = function - | One x -> One.trailing_zeros x - | Many x -> Many.trailing_zeros x - - -let hamdist a b = match a, b with - | One a, One b -> One.hamdist a b - | Many a, Many b -> Many.hamdist a b - | _ -> invalid_arg "Bitstring.hamdist" - - -let popcount = function - | One x -> One.popcount x - | Many x -> Many.popcount x - #+end_src - - - Example: - #+begin_example -Bitstring.(trailing_zeros (of_int 12));; -- : int = 2 - -Bitstring.(hamdist (of_int 15) (of_int 73));; -- : int = 3 - -Bitstring.(popcount (of_int 15));; -- : int = 4 - #+end_example - - #+begin_src ocaml :tangle (eval ml) :exports none -let rec to_list ?(accu=[]) = function - | t when (is_zero t) -> List.rev accu - | t -> let newlist = - (trailing_zeros t + 1)::accu - in - let b = logand t @@ minus_one t in - (to_list [@tailcall]) ~accu:newlist b - #+end_src - - #+begin_src ocaml :tangle (eval test-ml) :exports none -check_bool "to_list" true ([ 1 ; 3 ; 4 ; 6 ] = (to_list (of_int 45))); - #+end_src - - Example: - #+begin_example -Bitstring.(to_list (of_int 45));; -- : int list = [1; 3; 4; 6] - #+end_example - - #+begin_src ocaml :tangle (eval ml) :exports none -let permutations m n = - - let rec aux k u rest = - if k=1 then - List.rev (u :: rest) - else - let t = logor u @@ minus_one u in - let t' = plus_one t in - let not_t = lognot t in - let neg_not_t = neg not_t in - let t'' = shift_right (minus_one @@ logand not_t neg_not_t) (trailing_zeros u + 1) in - (* - let t'' = shift_right (minus_one (logand (lognot t) t')) (trailing_zeros u + 1) in - ,*) - (aux [@tailcall]) (k-1) (logor t' t'') (u :: rest) - in - aux (Util.binom n m) (minus_one (shift_left_one n m)) [] - #+end_src - - Example: - #+begin_example - Bitstring.permutations 2 4;; -- : Bitstring.t list = -[++--------------------------------------------------------------; - +-+-------------------------------------------------------------; - -++-------------------------------------------------------------; - +--+------------------------------------------------------------; - -+-+------------------------------------------------------------; - --++------------------------------------------------------------] - #+end_example - - #+begin_src ocaml :tangle (eval test-ml) :exports none -check "permutations" - (permutations 2 4 = List.map of_int - [ 3 ; 5 ; 6 ; 9 ; 10 ; 12 ]); - #+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 = function - | One x -> One.pp ppf x - | Many x -> Many.pp ppf x - #+end_src - -** Tests :noexport: - - #+begin_src ocaml :tangle (eval test-ml) :exports none -() - -let tests = [ - "all", `Quick, test_all; -] - #+end_src diff --git a/common/lib/angular_momentum.ml b/common/lib/angular_momentum.ml index d76486d..4f622b9 100644 --- a/common/lib/angular_momentum.ml +++ b/common/lib/angular_momentum.ml @@ -1,13 +1,3 @@ - - -(* An exception is raised when the ~Angular_momentum.t~ element can't - * be created. - * - * The ~kind~ is used to build shells, shell doublets, triplets or - * quartets, use in the two-electron operators. *) - - -(* [[file:~/QCaml/common/angular_momentum.org::*Type][Type:2]] *) type t = | S | P | D | F | G | H | I | J | K | L | M | N | O | Int of int @@ -15,42 +5,14 @@ type t = exception Angular_momentum_error of string type kind = - Singlet of t + | Singlet of t | Doublet of (t * t) | Triplet of (t * t * t) | Quartet of (t * t * t * t) open Powers -(* Type:2 ends here *) - -(* | ~of_char~ | Returns an ~Angular_momentum.t~ when a shell is given as a character (case insensitive) | - * | ~to_char~ | Converts the angular momentum into a char | - * | ~of_int~ | Returns a shell given an $l$ value. | - * | ~to_int~ | Returns the $l_{max}$ value of the shell | - * | ~to_string~ | Converts the angular momentum into a string | - * - * Example: - * #+begin_example - * Angular_momentum.of_char 'p';; - * - : Angular_momentum.t = P - * - * Angular_momentum.(to_char P);; - * - : char = 'P' - * - * Angular_momentum.of_int 2;; - * - : Angular_momentum.t = D - * - * Angular_momentum.(to_int D);; - * - : int = 2 - * - * Angular_momentum.(to_string D);; - * - : string = "D" - * #+end_example *) - - -(* [[file:~/QCaml/common/angular_momentum.org::*Conversions][Conversions:2]] *) let of_char = function | 's' | 'S' -> S | 'p' | 'P' -> P | 'd' | 'D' -> D | 'f' | 'F' -> F @@ -99,39 +61,22 @@ let of_int = function | 8 -> K | 9 -> L | 10 -> M | 11 -> N | 12 -> O | i -> Int i -(* Conversions:2 ends here *) -(* | ~n_functions~ | Returns the number of cartesian functions in a shell. | - * | ~zkey_array~ | Array of ~Zkey.t~, where each element is a a key associated with the the powers of $x,y,z$. | - * - * Example: - * #+begin_example - * Angular_momentum.(n_functions D) ;; - * - : int = 6 - * - * Angular_momentum.( zkey_array (Doublet (P,S)) );; - * - : Zkey.t array = - * [|< 01125899906842624 | 1, 0, 0, 0, 0, 0 >; - * < 01099511627776 | 0, 1, 0, 0, 0, 0 >; - * < 01073741824 | 0, 0, 1, 0, 0, 0 >|] - * #+end_example *) -(* [[file:~/QCaml/common/angular_momentum.org::*Shell functions][Shell functions:2]] *) let n_functions a = - let a = - to_int a - in + let a = to_int a in (a*a + 3*a + 2)/2 -let zkey_array_memo : (kind, Zkey.t array) Hashtbl.t = - Hashtbl.create 13 +(** Keeps data for function ~zkey_array~ *) +let zkey_array_memo : (kind, Zkey.t array) Hashtbl.t = Hashtbl.create 13 + let zkey_array a = - + let keys_1d l = let create_z { x ; y ; _ } = Powers.of_int_tuple (x,y,l-(x+y)) @@ -141,15 +86,17 @@ let zkey_array a = match y with | 0 -> (create_z xyz)::accu | _ -> let ynew = y-1 in - (create_y [@tailcall]) ( (create_z xyz)::accu) (Powers.of_int_tuple (x,ynew,z)) + (create_y [@tailcall]) ((create_z xyz)::accu) + (Powers.of_int_tuple (x,ynew,z)) in let rec create_x accu xyz = let { x ; z ;_ } = xyz in match x with | 0 -> (create_y [] xyz)@accu | _ -> let xnew = x-1 in - let ynew = l-xnew in - (create_x [@tailcall]) ((create_y [] xyz)@accu) (Powers.of_int_tuple (xnew, ynew, z)) + let ynew = l-xnew in + (create_x [@tailcall]) ((create_y [] xyz)@accu) + (Powers.of_int_tuple (xnew, ynew, z)) in create_x [] (Powers.of_int_tuple (l,0,0)) |> List.rev @@ -164,72 +111,59 @@ let zkey_array a = begin match a with | Singlet l1 -> - List.rev_map (fun x -> Zkey.of_powers_three x) (keys_1d @@ to_int l1) + List.rev_map (fun x -> Zkey.of_powers_three x) (keys_1d @@ to_int l1) | Doublet (l1, l2) -> - List.rev_map (fun a -> - List.rev_map (fun b -> Zkey.of_powers_six a b) (keys_1d @@ to_int l2) - ) (keys_1d @@ to_int l1) - |> List.concat - + List.rev_map (fun a -> + List.rev_map (fun b -> Zkey.of_powers_six a b) (keys_1d @@ to_int l2) + ) (keys_1d @@ to_int l1) + |> List.concat + | Triplet (l1, l2, l3) -> - - List.rev_map (fun a -> - List.rev_map (fun b -> - List.rev_map (fun c -> - Zkey.of_powers_nine a b c) (keys_1d @@ to_int l3) - ) (keys_1d @@ to_int l2) - |> List.concat - ) (keys_1d @@ to_int l1) - |> List.concat - + + List.rev_map (fun a -> + List.rev_map (fun b -> + List.rev_map (fun c -> + Zkey.of_powers_nine a b c) (keys_1d @@ to_int l3) + ) (keys_1d @@ to_int l2) + |> List.concat + ) (keys_1d @@ to_int l1) + |> List.concat + | Quartet (l1, l2, l3, l4) -> - - List.rev_map (fun a -> - List.rev_map (fun b -> - List.rev_map (fun c -> - List.rev_map (fun d -> - Zkey.of_powers_twelve a b c d) (keys_1d @@ to_int l4) - ) (keys_1d @@ to_int l3) - |> List.concat - ) (keys_1d @@ to_int l2) - |> List.concat - ) (keys_1d @@ to_int l1) - |> List.concat + + List.rev_map (fun a -> + List.rev_map (fun b -> + List.rev_map (fun c -> + List.rev_map (fun d -> + Zkey.of_powers_twelve a b c d) (keys_1d @@ to_int l4) + ) (keys_1d @@ to_int l3) + |> List.concat + ) (keys_1d @@ to_int l2) + |> List.concat + ) (keys_1d @@ to_int l1) + |> List.concat end |> List.rev |> Array.of_list in Hashtbl.add zkey_array_memo a result; result -(* Shell functions:2 ends here *) -(* Example: - * #+begin_example - * Angular_momentum.(D + P);; - * - : Angular_momentum.t = F - * - * Angular_momentum.(F - P);; - * - : Angular_momentum.t = D - * #+end_example *) - - -(* [[file:~/QCaml/common/angular_momentum.org::*Arithmetic][Arithmetic:2]] *) let ( + ) a b = of_int ( (to_int a) + (to_int b) ) let ( - ) a b = of_int ( (to_int a) - (to_int b) ) -(* Arithmetic:2 ends here *) -(* [[file:~/QCaml/common/angular_momentum.org::*Printers][Printers:2]] *) + + let pp_string ppf x = Format.fprintf ppf "@[%s@]" (to_string x) - + let pp_int ppf x = Format.fprintf ppf "@[%d@]" (to_int x) - + let pp = pp_string -(* Printers:2 ends here *) diff --git a/common/lib/angular_momentum.mli b/common/lib/angular_momentum.mli index 527ab18..2ff3b5c 100644 --- a/common/lib/angular_momentum.mli +++ b/common/lib/angular_momentum.mli @@ -1,58 +1,89 @@ -(* Type - * - * <<<~Angular_momentum.t~>>> - * #+NAME: types *) +(** Azimuthal quantum number, repsesented as \( s,p,d,\dots \) *) -(* [[file:~/QCaml/common/angular_momentum.org::types][types]] *) type t = | S | P | D | F | G | H | I | J | K | L | M | N | O | Int of int - + exception Angular_momentum_error of string +(** An exception is raised when the ~Angular_momentum.t~ element can't + be created. *) type kind = - Singlet of t + | Singlet of t | Doublet of (t * t) | Triplet of (t * t * t) | Quartet of (t * t * t * t) -(* types ends here *) -(* Conversions *) +(** Type ~kind~ is used to build shells, shell doublets, triplets or + quartets, use in the two-electron operators. *) -(* [[file:~/QCaml/common/angular_momentum.org::*Conversions][Conversions:1]] *) +(** Conversions *) + val of_char : char -> t +(** Returns an ~Angular_momentum.t~ when a shell is given as a character + (case insensitive) + +Angular_momentum.of_char 'p';; +- : Angular_momentum.t = P +*) + val to_char : t -> char +(** Converts the angular momentum into a char + +Angular_momentum.(to_char P);; +- : char = 'P' +*) val to_int : t -> int +(** Returns the $l_{max}$ value of the shell + +Angular_momentum.(to_int D);; +- : int = 2 +*) + val of_int : int -> t +(** Returns a shell given an $l$ value + +Angular_momentum.of_int 2;; +- : Angular_momentum.t = D +*) val to_string : t -> string -(* Conversions:1 ends here *) +(** Converts the angular momentum into a string -(* Shell functions *) +Angular_momentum.(to_string D);; +- : string = "D" +*) -(* [[file:~/QCaml/common/angular_momentum.org::*Shell functions][Shell functions:1]] *) +(** Shell functions *) + val n_functions : t -> int +(** Returns the number of cartesian functions in a shell. *) + val zkey_array : kind -> Zkey.t array -(* Shell functions:1 ends here *) - -(* Arithmetic *) +(** Array of ~Zkey.t~, where each element is a a key associated with the + the powers of $x,y,z$. *) -(* [[file:~/QCaml/common/angular_momentum.org::*Arithmetic][Arithmetic:1]] *) +(** Arithmetic *) + val ( + ) : t -> t -> t +(** +Angular_momentum.(D + P);; +- : Angular_momentum.t = F +*) + val ( - ) : t -> t -> t -(* Arithmetic:1 ends here *) - -(* Printers - * - * Printers can print as a string (default) or as an integer. *) +(** +Angular_momentum.(F - P);; +- : Angular_momentum.t = D +*) -(* [[file:~/QCaml/common/angular_momentum.org::*Printers][Printers:1]] *) +(** Printers *) + val pp : Format.formatter -> t -> unit val pp_string : Format.formatter -> t -> unit val pp_int : Format.formatter -> t -> unit -(* Printers:1 ends here *) diff --git a/common/lib/bitstring.ml b/common/lib/bitstring.ml index 179eb16..51dfa51 100644 --- a/common/lib/bitstring.ml +++ b/common/lib/bitstring.ml @@ -1,12 +1,8 @@ -(* Single-integer implementation :noexport: *) +(** Single-integer implementation *) - -(* [[file:~/QCaml/common/bitstring.org::*Single-integer implementation][Single-integer implementation:1]] *) module One = struct - let of_int x = - assert (x > 0); x - + let of_int x = assert (x > 0); x let numbits _ = 63 let zero = 0 let is_zero x = x = 0 @@ -38,12 +34,9 @@ module One = struct (Z.of_int s) end -(* Single-integer implementation:1 ends here *) - -(* Zarith implementation :noexport: *) -(* [[file:~/QCaml/common/bitstring.org::*Zarith implementation][Zarith implementation:1]] *) +(** Zarith implementation *) module Many = struct let of_z x = x @@ -71,58 +64,21 @@ module Many = struct Format.fprintf ppf "@[@[%a@]@]" (Util.pp_bitstring (Z.numbits s)) s end -(* Zarith implementation:1 ends here *) -(* [[file:~/QCaml/common/bitstring.org::*Type][Type:2]] *) type t = | One of int | Many of Z.t -(* Type:2 ends here *) -(* | ~of_int~ | Creates a bit string from an ~int~ | - * | ~of_z~ | Creates a bit string from an ~Z.t~ multi-precision integer | - * | ~zero~ | ~zero n~ creates a zero bit string with ~n~ bits | - * | ~is_zero~ | True if all the bits of the bit string are zero. | - * | ~numbits~ | Returns the number of bits used to represent the bit string | - * | ~testbit~ | ~testbit t n~ is true if the ~n~-th bit of the bit string ~t~ is set to ~1~ | - * | ~neg~ | Returns the negative of the integer interpretation of the bit string | - * | ~shift_left~ | ~shift_left t n~ returns a new bit strings with all the bits shifted ~n~ positions to the left | - * | ~shift_right~ | ~shift_right t n~ returns a new bit strings with all the bits shifted ~n~ positions to the right | - * | ~shift_left_one~ | ~shift_left_one size n~ returns a new bit strings with the ~n~-th bit set to one. It is equivalent as shifting ~1~ by ~n~ bits to the left, ~size~ is the total number of bits of the bit string | - * | ~logor~ | Bitwise logical or | - * | ~logxor~ | Bitwise logical exclusive or | - * | ~logand~ | Bitwise logical and | - * | ~lognot~ | Bitwise logical negation | - * | ~plus_one~ | Takes the integer representation of the bit string and adds one | - * | ~minus_one~ | Takes the integer representation of the bit string and removes one | - * | ~hamdist~ | Returns the Hamming distance, i.e. the number of bits differing between two bit strings | - * | ~trailing_zeros~ | Returns the number of trailing zeros in the bit string | - * | ~permutations~ | ~permutations m n~ generates the list of all possible ~n~-bit strings with ~m~ bits set to ~1~. Algorithm adapted from [[https://graphics.stanford.edu/~seander/bithacks.html#NextBitPermutation][Bit twiddling hacks]] | - * | ~popcount~ | Returns the number of bits set to one in the bit string | - * | ~to_list~ | Converts a bit string into a list of integers indicating the positions where the bits are set to ~1~. The first value for the position is not ~0~ but ~1~ | *) - - -(* [[file:~/QCaml/common/bitstring.org::*General implementation][General implementation:2]] *) let of_int x = One (One.of_int x) let of_z x = if Z.numbits x < 64 then One (Z.to_int x) else Many (Many.of_z x) -(* General implementation:2 ends here *) - -(* Example: - * #+begin_example - * Bitstring.of_int 15;; - * - : Bitstring.t = - * ++++------------------------------------------------------------ - * #+end_example *) - -(* [[file:~/QCaml/common/bitstring.org::*General implementation][General implementation:4]] *) let zero = function | n when n < 64 -> One (One.zero) | _ -> Many (Many.zero) @@ -159,9 +115,9 @@ let shift_left_one = function let testbit = function | One x -> One.testbit x | Many x -> Many.testbit x -(* General implementation:4 ends here *) -(* [[file:~/QCaml/common/bitstring.org::*General implementation][General implementation:6]] *) + +(** General implementation *) let logor a b = match a,b with | One a, One b -> One (One.logor a b) @@ -186,9 +142,8 @@ let logand a b = let lognot = function | One x -> One (One.lognot x) | Many x -> Many (Many.lognot x) -(* General implementation:6 ends here *) -(* [[file:~/QCaml/common/bitstring.org::*General implementation][General implementation:8]] *) + let minus_one = function | One x -> One (One.minus_one x) | Many x -> Many (Many.minus_one x) @@ -197,25 +152,8 @@ let minus_one = function let plus_one = function | One x -> One (One.plus_one x) | Many x -> Many (Many.plus_one x) -(* General implementation:8 ends here *) - - -(* Example: - * #+begin_example - * Bitstring.(plus_one (of_int 15));; - * - : Bitstring.t = - * ----+----------------------------------------------------------- - * - * Bitstring.(minus_one (of_int 15));; - * - : Bitstring.t = - * -+++------------------------------------------------------------ - * #+end_example *) - - - -(* [[file:~/QCaml/common/bitstring.org::*General implementation][General implementation:9]] *) let trailing_zeros = function | One x -> One.trailing_zeros x | Many x -> Many.trailing_zeros x @@ -230,24 +168,10 @@ let hamdist a b = match a, b with let popcount = function | One x -> One.popcount x | Many x -> Many.popcount x -(* General implementation:9 ends here *) -(* Example: - * #+begin_example - * Bitstring.(trailing_zeros (of_int 12));; - * - : int = 2 - * - * Bitstring.(hamdist (of_int 15) (of_int 73));; - * - : int = 3 - * - * Bitstring.(popcount (of_int 15));; - * - : int = 4 - * #+end_example *) - - (* [[file:~/QCaml/common/bitstring.org::*General implementation][General implementation:10]] *) let rec to_list ?(accu=[]) = function | t when (is_zero t) -> List.rev accu @@ -267,7 +191,6 @@ let rec to_list ?(accu=[]) = function * #+end_example *) -(* [[file:~/QCaml/common/bitstring.org::*General implementation][General implementation:12]] *) let permutations m n = let rec aux k u rest = @@ -285,10 +208,9 @@ let permutations m n = (aux [@tailcall]) (k-1) (logor t' t'') (u :: rest) in aux (Util.binom n m) (minus_one (shift_left_one n m)) [] -(* General implementation:12 ends here *) -(* [[file:~/QCaml/common/bitstring.org::*Printers][Printers:2]] *) + + let pp ppf = function | One x -> One.pp ppf x | Many x -> Many.pp ppf x -(* Printers:2 ends here *) diff --git a/common/lib/bitstring.mli b/common/lib/bitstring.mli index 10c19ed..fa85dd9 100644 --- a/common/lib/bitstring.mli +++ b/common/lib/bitstring.mli @@ -1,47 +1,166 @@ -(* Type - * - * <<<~Bitstring.t~>>> *) +(** Bitstring + We define here a data type to handle bit strings efficiently. When + the bit string contains less than 64 bits, it is stored internally + in a 63-bit integer and uses bitwise instructions. When more than + 63 bits are required, the =zarith= library is used to consider the + bit string as a multi-precision integer. + +*) -(* [[file:~/QCaml/common/bitstring.org::*Type][Type:1]] *) type t -(* Type:1 ends here *) -(* General implementation *) - - -(* [[file:~/QCaml/common/bitstring.org::*General implementation][General implementation:1]] *) val of_int : int -> t +(** Creates a bit string from an ~int~ + +Bitstring.of_int 15;; +- : Bitstring.t = +++++------------------------------------------------------------ + *) + val of_z : Z.t -> t +(** Creates a bit string from an ~Z.t~ multi-precision integer *) + val zero : int -> t +(** ~zero n~ creates a zero bit string with ~n~ bits *) val is_zero : t -> bool +(** True if all the bits of the bit string are zero. *) + val numbits : t -> int +(** Returns the number of bits used to represent the bit string *) + val testbit : t -> int -> bool +(** ~testbit t n~ is true if the ~n~-th bit of the bit string ~t~ is set to ~1~ + +Bitstring.(testbit (of_int 15) 3);; +- : bool = true + +Bitstring.(testbit (of_int 15) 4);; +- : bool = false +*) val neg : t -> t +(** Returns the negative of the integer interpretation of the bit string *) + val shift_left : t -> int -> t +(** ~shift_left t n~ returns a new bit strings with all the bits shifted ~n~ + positions to the left + +Bitstring.(shift_left (of_int 15) 2);; +- : Bitstring.t = +--++++---------------------------------------------------------- +*) + val shift_right : t -> int -> t +(** ~shift_right t n~ returns a new bit strings with all the bits shifted ~n~ + positions to the right + +Bitstring.(shift_right (of_int 15) 2);; +- : Bitstring.t = +++-------------------------------------------------------------- +*) + val shift_left_one : int -> int -> t +(** shift_left_one size n~ returns a new bit strings with the ~n~-th bit set to + one. It is equivalent as shifting ~1~ by ~n~ bits to the left, ~size~ is the + total number of bits of the bit string + +Bitstring.shift_left_one 32 4;; +- : Bitstring.t = +----+----------------------------------------------------------- +*) val logor : t -> t -> t +(** Bitwise logical or + +Bitstring.(logor (of_int 15) (of_int 73));; +- : Bitstring.t = +++++--+--------------------------------------------------------- +*) + val logxor : t -> t -> t +(** Bitwise logical exclusive or + +Bitstring.(logxor (of_int 15) (of_int 73));; +- : Bitstring.t = +-++---+--------------------------------------------------------- +*) + val logand : t -> t -> t +(** Bitwise logical and + +Bitstring.(logand (of_int 15) (of_int 10));; +- : Bitstring.t = +-+-+------------------------------------------------------------ +*) + val lognot : t -> t +(** Bitwise logical negation *) val plus_one : t -> t +(** Takes the integer representation of the bit string and adds one + +Bitstring.(plus_one (of_int 15));; +- : Bitstring.t = +----+----------------------------------------------------------- +*) + val minus_one : t -> t +(** Takes the integer representation of the bit string and removes one + +Bitstring.(minus_one (of_int 15));; +- : Bitstring.t = +-+++------------------------------------------------------------ +*) + val hamdist : t -> t -> int +(** Returns the Hamming distance, i.e. the number of bits differing between two + bit strings + +Bitstring.(hamdist (of_int 15) (of_int 73));; +- : int = 3 +*) + val trailing_zeros : t -> int +(** Returns the number of trailing zeros in the bit string + +Bitstring.(trailing_zeros (of_int 12));; +- : int = 2 +*) + val popcount : t -> int +(** Returns the number of bits set to one in the bit string + +Bitstring.(popcount (of_int 15));; +- : int = 4 +*) + val to_list : ?accu:(int list) -> t -> int list +(** Converts a bit string into a list of integers indicating the positions where + the bits are set to ~1~. The first value for the position is not ~0~ but ~1~ + +Bitstring.(to_list (of_int 45));; +- : int list = [1; 3; 4; 6] +*) + val permutations : int -> int -> t list -(* General implementation:1 ends here *) +(** ~permutations m n~ generates the list of all possible ~n~-bit strings with + ~m~ bits set to ~1~. Algorithm adapted from + [[https://graphics.stanford.edu/~seander/bithacks.html#NextBitPermutation][Bit twiddling hacks]] -(* Printers *) +Bitstring.permutations 2 4;; +- : Bitstring.t list = +[++--------------------------------------------------------------; + +-+-------------------------------------------------------------; + -++-------------------------------------------------------------; + +--+------------------------------------------------------------; + -+-+------------------------------------------------------------; + --++------------------------------------------------------------] +*) -(* [[file:~/QCaml/common/bitstring.org::*Printers][Printers:1]] *) +(** Printers *) + val pp : Format.formatter -> t -> unit -(* Printers:1 ends here *) diff --git a/common/lib/non_negative_float.ml b/common/lib/non_negative_float.ml index 0a83876..0b025f1 100644 --- a/common/lib/non_negative_float.ml +++ b/common/lib/non_negative_float.ml @@ -1,14 +1,5 @@ -(* [[file:~/QCaml/common/non_negative_float.org::*Type][Type:2]] *) type t = float -(* Type:2 ends here *) - - -(* The ~of_float~ function checks that the float is non-negative. - * The unsafe variant doesn't do this check. *) - - -(* [[file:~/QCaml/common/non_negative_float.org::*Conversions][Conversions:2]] *) let of_float x = if x < 0. then invalid_arg (__FILE__^": of_float"); x @@ -21,4 +12,3 @@ let to_string x = let of_string x = let f = float_of_string x in of_float f -(* Conversions:2 ends here *) diff --git a/common/lib/non_negative_float.mli b/common/lib/non_negative_float.mli index efdfb63..6872b69 100644 --- a/common/lib/non_negative_float.mli +++ b/common/lib/non_negative_float.mli @@ -1,19 +1,14 @@ -(* Type - * - * <<<~Non_negative_float.t~>> *) - -(* [[file:~/QCaml/common/non_negative_float.org::*Type][Type:1]] *) +(** $x \ge 0$ *) type t = private float -(* Type:1 ends here *) -(* Conversions *) +(** Conversions *) - -(* [[file:~/QCaml/common/non_negative_float.org::*Conversions][Conversions:1]] *) val of_float : float -> t -val unsafe_of_float : float -> t -val to_float : t -> float +(** Checks that the float is non-negative. *) +val unsafe_of_float : float -> t +(** Fast conversion without checking that the float is non-negative. *) + +val to_float : t -> float val of_string : string -> t val to_string : t -> string -(* Conversions:1 ends here *) diff --git a/common/lib/qcaml.ml b/common/lib/qcaml.ml index 945977c..bb5cf06 100644 --- a/common/lib/qcaml.ml +++ b/common/lib/qcaml.ml @@ -1,11 +1,3 @@ - - -(* | ~root~ | Path to the QCaml source directory | - * | ~name~ | ~"QCaml"~ | - * | ~num_domains~ | Number of threads in parallel computations | *) - - -(* [[file:~/QCaml/common/qcaml.org::*QCaml][QCaml:2]] *) let name = "QCaml" let root = @@ -29,4 +21,3 @@ let num_domains = with Not_found -> 1 in result - 1 -(* QCaml:2 ends here *) diff --git a/common/lib/qcaml.mli b/common/lib/qcaml.mli index 0d21c3c..cc4a542 100644 --- a/common/lib/qcaml.mli +++ b/common/lib/qcaml.mli @@ -1,13 +1,11 @@ -(* QCaml - * :PROPERTIES: - * :header-args: :noweb yes :comments both - * :END: - * - * QCaml-specific parameters *) - +(** QCaml-specific parameters *) -(* [[file:~/QCaml/common/qcaml.org::*QCaml][QCaml:1]] *) val root : string +(** Path to the QCaml source directory *) + val name : string +(** "QCaml" *) + val num_domains : int -(* QCaml:1 ends here *) +(** Number of threads in parallel computations. Obtained from the environment variable + ~OMP_NUM_THREADS~. *) diff --git a/common/non_negative_float.org b/common/non_negative_float.org deleted file mode 100644 index ca79f3b..0000000 --- a/common/non_negative_float.org +++ /dev/null @@ -1,55 +0,0 @@ -#+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 - -* Non-negative float - :PROPERTIES: - :header-args: :noweb yes :comments both - :END: - -** Type - - <<<~Non_negative_float.t~>> - #+begin_src ocaml :tangle (eval mli) -type t = private float - #+end_src - - #+begin_src ocaml :tangle (eval ml) :exports none -type t = float - #+end_src - -** Conversions - - #+begin_src ocaml :tangle (eval mli) -val of_float : float -> t -val unsafe_of_float : float -> t -val to_float : t -> float - -val of_string : string -> t -val to_string : t -> string - #+end_src - - The ~of_float~ function checks that the float is non-negative. - The unsafe variant doesn't do this check. - - #+begin_src ocaml :tangle (eval ml) :exports none -let of_float x = - if x < 0. then invalid_arg (__FILE__^": of_float"); - x - -external to_float : t -> float = "%identity" -external unsafe_of_float : float -> t = "%identity" - -let to_string x = - let f = to_float x in string_of_float f - -let of_string x = - let f = float_of_string x in of_float f - #+end_src diff --git a/common/qcaml.org b/common/qcaml.org deleted file mode 100644 index a6a29a3..0000000 --- a/common/qcaml.org +++ /dev/null @@ -1,55 +0,0 @@ -#+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 - -* QCaml - :PROPERTIES: - :header-args: :noweb yes :comments both - :END: - - QCaml-specific parameters - - #+begin_src ocaml :tangle (eval mli) -val root : string -val name : string -val num_domains : int - #+end_src - - | ~root~ | Path to the QCaml source directory | - | ~name~ | ~"QCaml"~ | - | ~num_domains~ | Number of threads in parallel computations | - - #+begin_src ocaml :tangle (eval ml) :exports none -let name = "QCaml" - -let root = - let rec chop = function - | [] -> [] - | x :: _ as l when x = name -> l - | _ :: rest -> chop rest - in - String.split_on_char Filename.dir_sep.[0] (Sys.getcwd ()) - |> List.rev - |> chop - |> List.rev - |> String.concat Filename.dir_sep - - -let num_domains = - let result = - try - Unix.getenv "OMP_NUM_THREADS" - |> int_of_string - with Not_found -> 1 - in - result - 1 - - #+end_src - diff --git a/docs/ao.html b/docs/ao.html index aebe4ac..f345dff 100644 --- a/docs/ao.html +++ b/docs/ao.html @@ -3,7 +3,7 @@ "http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd"> - + Atomic Orbitals @@ -238,28 +238,6 @@ org_html_manager.setup(); // activate after the parameters are set /*]]>*///--> // @license-end - -
@@ -272,37 +250,13 @@ org_html_manager.setup(); // activate after the parameters are set

Table of Contents

-
-

1 Summary

+
+

1 Summary

This modules contains the data structures used to characterize the @@ -310,444 +264,10 @@ atomic basis set.

- -
-

2 Gaussian basis

-
-

-Data structure for Gaussian Atomic Orbitals: -

- -

-\[ - \chi_i(\mathbf{r}) = P_i(\mathbf{r}) \sum_k c_k \exp\left( -\alpha_k (\mathbf{r-R_A})^2 \right) - \] -

- -

-where the polynomial \(P_i\) and the Gaussian part are both centered on -nucleus \(A\). -

-
- -
-

2.1 Type

-
-

-Gaussian_basis.t -

-
-
open Common
-open Particles
-open Linear_algebra
-open Gaussian_integrals
-open Operators
-
-type t
-
-
-
-
- -
-

2.2 Access

-
-
-
val basis             : t -> Gaussian.Basis.t
-val cartesian         : t -> bool
-val ee_ints           : t -> Eri.t
-val ee_lr_ints        : t -> Eri_long_range.t
-val eN_ints           : t -> Electron_nucleus.t
-val f12_ints          : t -> F12.t
-val f12_over_r12_ints : t -> Screened_eri.t
-val kin_ints          : t -> Kinetic.t
-val multipole         : t -> Multipole.t
-val ortho             : t -> Orthonormalization.t
-val overlap           : t -> Overlap.t
-val size              : t -> int
-
-
- - - - --- -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
basisOne-electron basis set
cartesianIf true, use cartesian Gaussians (6d, 10f, …)
ee_intsElectron-electron potential integrals
ee_lr_intsElectron-electron long-range potential integrals
eN_intsElectron-nucleus potential integrals
f12_intsElectron-electron potential integrals
f12_over_r12_intsElectron-electron potential integrals
kin_intsKinetic energy integrals
multipoleMultipole matrices
orthoOrthonormalization matrix of the overlap
overlapOverlap matrix
sizeNumber of atomic orbitals
-
-
- -
-

2.3 Computation

-
-
-
val values : t ->
-             Coordinate.t ->
-             Gaussian.Basis.t Vector.t
-
-
- - - - --- -- - - - - - - -
valuesReturns the values of all the AOs evaluated at a given point
-
-
- -
-

2.4 Creation

-
-
-
val make : basis:Gaussian.Basis.t ->
-           ?operators:Operator.t list ->
-           ?cartesian:bool ->
-           Nuclei.t ->
-           t
-
-
- -

-Creates the data structure for atomic orbitals from a Gaussian basis and the -molecular geometry Nuclei.t. -

- -

-Defaults: -

-
    -
  • operators : []
    • -
  • cartesian : false
    • -
- -

-Example: -

-
-let b = Ao.Basis_gaussian.make ~basis nuclei ;;
-val b : Ao.Basis_gaussian.t = Gaussian Basis, spherical, 15 AOs
-
-
-
- -
-

2.5 Printers

-
-
-
val pp : Format.formatter -> t -> unit
-
-
-
-
-
- -
-

3 Basis

-
-

-Data structure for Atomic Orbitals. -

-
- -
-

3.1 Polymorphic types

-
-

-Basis.t -

-
-
type t =
-  | Unknown
-  | Gaussian of Basis_gaussian.t
-
-
-
-
- -
-

3.2 Types

-
-

-Basis.t -

-
-
type t
-type ao = Ao_dim.t
-
-open Common
-open Particles
-open Operators
-open Linear_algebra
-
-
-
-
- -
-

3.3 Conversions

-
-
-
val of_nuclei_and_basis_filename :
-  ?kind:[> `Gaussian ] ->
-  ?operators:Operator.t list ->
-  ?cartesian:bool ->
-  nuclei:Nuclei.t ->
-  string ->
-  t
-
-
- - - - --- -- - - - - - - -
of_nuclei_and_basis_filenameCreates the data structure for the atomic orbitals basis from a molecule Nuclei.t and the name of the basis-set file
- -

-Defaults: -

-
    -
  • kind : `Gaussian
    • -
  • operators : []
    • -
  • cartesian : false
    • -
- -

-Example: -

-
-let b = Ao.Basis.of_nuclei_and_basis_filename ~nuclei filename;;
-val b : Ao.Basis.t = Gaussian Basis, spherical, 15 AOs
-
-
-
- -
-

3.4 Access

-
-
-
val size              : t -> int
-val ao_basis          : t -> Basis_poly.t
-val overlap           : t -> (ao,ao) Matrix.t
-val multipole         : t -> string -> (ao,ao) Matrix.t
-val ortho             : t -> (ao,'a) Matrix.t
-val eN_ints           : t -> (ao,ao) Matrix.t
-val kin_ints          : t -> (ao,ao) Matrix.t
-val ee_ints           : t -> ao Four_idx_storage.t
-val ee_lr_ints        : t -> ao Four_idx_storage.t
-val f12_ints          : t -> ao Four_idx_storage.t
-val f12_over_r12_ints : t -> ao Four_idx_storage.t
-val cartesian         : t -> bool
-val values            : t -> Coordinate.t -> ao Vector.t
-
-
- - - - --- -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
sizeNumber of atomic orbitals in the AO basis set
ao_basisOne-electron basis set
overlapOverlap matrix
multipoleMultipole matrices
orthoOrthonormalization matrix of the overlap
eN_intsElectron-nucleus potential integrals
kin_intsKinetic energy integrals
ee_intsElectron-electron potential integrals
ee_lr_intsElectron-electron long-range potential integrals
f12_intsElectron-electron potential integrals
f12_over_r12_intsElectron-electron potential integrals
cartesianIf true, use cartesian Gaussians (6d, 10f, …)
valuesValues of the AOs evaluated at a given point
-
-
- - -
-

3.5 TREXIO

-
-
-
-

3.5.1 Read

-
-
-
(*
-val of_trexio : Trexio.trexio_file -> t
-*)
-
-
-
-
- -
-

3.5.2 Write

-
-
-
(*
-val to_trexio : Trexio.trexio_file -> t -> unit
-*)
-
-
-
-
-
- -
-

3.6 Printers

-
-
-
val pp : Format.formatter -> t -> unit
-
-
-
-
-

Author: Anthony Scemama

-

Created: 2023-04-24 Mon 18:06

+

Created: 2023-06-27 Tue 09:39

Validate

diff --git a/docs/top.html b/docs/top.html index 3077cea..ebe657c 100644 --- a/docs/top.html +++ b/docs/top.html @@ -3,7 +3,7 @@ "http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd"> - + Top-level @@ -250,18 +250,18 @@ org_html_manager.setup(); // activate after the parameters are set

Table of Contents

-
-

1 Summary

+
+

1 Summary

Author: Anthony Scemama

-

Created: 2023-06-26 Mon 15:36

+

Created: 2023-06-27 Tue 09:39

Validate

diff --git a/top/lib/install_printers.ml b/top/lib/install_printers.ml index 09ffdd8..bf38de5 100644 --- a/top/lib/install_printers.ml +++ b/top/lib/install_printers.ml @@ -1,10 +1,6 @@ (* [[file:~/QCaml/top/install_printers.org::*Intall printers][Intall printers:3]] *) let printers = [ - "Ao.Basis.pp" ; - "Ao.Basis_gaussian.pp" ; - "Common.Angular_momentum.pp" ; - "Common.Bitstring.pp" ; "Common.Charge.pp" ; "Common.Coordinate.pp" ; "Common.Powers.pp" ;