LCPQ/QCaml
10
1
Fork 0
mirror of https://gitlab.com/scemama/QCaml.git synced 2024-09-01 09:13:40 +02:00
Code Releases Activity

Compare commits

base: LCPQ:1e778594ec7f3cdc29bb85a58822ab6e66b0a33c
Branches Tags
LCPQ:master
LCPQ:dev
LCPQ:smallarray
..
compare: LCPQ:701ac05e2361196735cd419c8c8b129d1f3cd4ae
Branches Tags
LCPQ:dev
LCPQ:smallarray
LCPQ:master

No commits in common. "1e778594ec7f3cdc29bb85a58822ab6e66b0a33c" and "701ac05e2361196735cd419c8c8b129d1f3cd4ae" have entirely different histories.
1e778594ec ... 701ac05e23

13 changed files with 105 additions and 110 deletions
Show Stats Expand all files Collapse all files

16
.merlin
View File

@ -1,8 +1,12 @@
EXCLUDE_QUERY_DIR
B _build/default/.run_hartree_fock.eobjs/byte
B _build/default/Basis/.Basis.objs/byte
B _build/default/Nuclei/.Nuclei.objs/byte
PKG str unix bigarray lacaml zarith mpi alcotest
S .
S Basis
S Test
S Nuclei
FLG -w @a-4-29-40-41-42-44-45-48-58-59-60-40 -strict-sequence -strict-formats -short-paths -keep-locs
S Parallel
S CI
S Utils
S Basis
S SCF
S MOBasis
S CI
B _build/**

6
Basis/TwoElectronRRVectorized.ml
View File

@ -756,14 +756,14 @@ let contracted_class_shell_pairs ~zero_m ?schwartz_p ?schwartz_q shell_p shell_q
in
let empty = Array.make (maxm+1) 0. in
let center_qc_tmp = Array.init nq (fun cd ->
Coordinate.make { Coordinate.
Bohr.make { Point.
x = (center_qc Co.X).(cd) ;
y = (center_qc Co.Y).(cd) ;
z = (center_qc Co.Z).(cd) ;
})
in
Array.iteri (fun ab shell_ab ->
let center_pa = Coordinate.make { Coordinate.
let center_pa = Bohr.make { Point.
x = (center_pa Co.X).(ab) ;
y = (center_pa Co.Y).(ab) ;
z = (center_pa Co.Z).(ab) ;
@ -789,7 +789,7 @@ let contracted_class_shell_pairs ~zero_m ?schwartz_p ?schwartz_q shell_p shell_q
in
zero_m Zp.{
maxm ; expo_p_inv ; expo_q_inv ; norm_pq_sq ;
center_pq = Coordinate.make Coordinate.{x ; y ; z} ;
center_pq = Bohr.make Point.{x ; y ; z} ;
center_pa ; center_qc = center_qc_tmp.(cd) ;
normalization ;
}

17
Nuclei/Nuclei.ml
View File

@ -28,9 +28,7 @@ let of_xyz_file filename =
let of_zmt_string buffer =
Zmatrix.of_string buffer
|> Zmatrix.to_xyz
|> Array.map (fun (e,x,y,z) ->
(e, Coordinate.(angstrom_to_bohr @@ make_angstrom { x ; y ; z} ))
)
|> Array.map (fun (e,x,y,z) -> (e, Coordinate.angstrom_to_bohr (Angstrom.make {Point.x ; y ; z} )))
let of_zmt_file filename =
@ -58,13 +56,12 @@ let to_string atoms =
-----------------------------------------------------------------------
" ^
(Array.mapi (fun i (e, coord) ->
let open Coordinate in
let coord =
bohr_to_angstrom coord
Coordinate.bohr_to_angstrom coord
in
Printf.sprintf " %5d %5d %5s %12.6f %12.6f %12.6f"
(i+1) (Element.to_int e) (Element.to_string e)
coord.x coord.y coord.z
coord.Angstrom.x coord.Angstrom.y coord.Angstrom.z
) atoms
|> Array.to_list
|> String.concat "\n" ) ^
@ -103,12 +100,11 @@ let charge nuclei =
let to_xyz_string t =
[ string_of_int (Array.length t) ; "" ] @
( Array.mapi (fun i (e, coord) ->
let open Coordinate in
let coord =
bohr_to_angstrom coord
Coordinate.bohr_to_angstrom coord
in
Printf.sprintf " %5s %12.6f %12.6f %12.6f"
(Element.to_string e) coord.x coord.y coord.z
(Element.to_string e) coord.Angstrom.x coord.Angstrom.y coord.Angstrom.z
) t
|> Array.to_list )
|> String.concat "\n"
@ -121,9 +117,8 @@ let to_t2_string t =
(2 * Array.length t);
"# Znuc x y z" ]
@ (Array.mapi (fun i (e, coord) ->
let open Coordinate in
Printf.sprintf " %5f %12.6f %12.6f %12.6f"
(Element.to_int e |> float_of_int) coord.x coord.y coord.z
(Element.to_int e |> float_of_int) coord.Bohr.x coord.Bohr.y coord.Bohr.z
) t
|> Array.to_list )
|> String.concat "\n"

2
Nuclei/Nuclei.mli
View File

@ -2,7 +2,7 @@
of tuples ({!Element.t}, {!Coordinate.t}).
*)
type t = (Element.t * Coordinate.t) array
type t = (Element.t * Bohr.t) array
val of_xyz_file : string -> t

2
Nuclei/Xyz_ast.mli
View File

@ -4,7 +4,7 @@
type nucleus =
{
element: Element.t ;
coord : Coordinate.angstrom Coordinate.point;
coord : Angstrom.t;
}
type xyz_file =

5
Nuclei/Xyz_parser.mly
View File

@ -3,9 +3,10 @@
%{
exception InputError of string
let make_angstrom x y z =
Coordinate.(make_angstrom {
Angstrom.make {
Point.
x ; y ; z
})
}
let output_of f x y z =
let a = make_angstrom x y z in

1
Utils/Angstrom.ml Normal file
View File

@ -0,0 +1 @@
include Bohr

4
Utils/Angstrom.mli Normal file
View File

@ -0,0 +1,4 @@
(** 3D point in cartesian coordinates, where units are Angstroms. *)
include module type of Bohr

11
Utils/Bohr.ml Normal file
View File

@ -0,0 +1,11 @@
type t = {
x : float ;
y : float ;
z : float ;
}
external make : Point.t -> t = "%identity"

12
Utils/Bohr.mli Normal file
View File

@ -0,0 +1,12 @@
(** 3D point in cartesian coordinates, where units are atomic units (Bohr). *)
type t = private {
x : float ;
y : float ;
z : float ;
}
val make : Point.t -> t

16
Utils/Cholesky.ml
View File

@ -3,16 +3,13 @@ open Lacaml.D
let full_ldl m_A =
let n = Mat.dim1 m_A in
assert (Mat.dim2 m_A = n);
let v_D = Vec.make0 n in
let m_Lt = Mat.identity n in
let v_D_inv = Vec.make0 n in
let compute_d j =
let l_jk =
Mat.col m_Lt j
@ -23,7 +20,6 @@ let full_ldl m_A =
m_A.{j,j} -. dot ~n:(j-1) l_jk l_jk__d_k
in
let compute_l i =
let l_ik__d_k =
Mat.col m_Lt i
@ -37,7 +33,6 @@ let full_ldl m_A =
v_D_inv.{j} *. ( m_A.{j,i} -. dot ~n:(j-1) l_ik__d_k l_jk )
in
for i=1 to n do
for j=1 to (i-1) do
m_Lt.{j,i} <- compute_l i j;
@ -90,7 +85,6 @@ let pivoted_ldl threshold m_A =
v_D_inv.{j} *. ( m_A.{pi.(j),pi.(i)} -. dot ~n:(j-1) l_ik__d_k l_jk )
in
let maxloc v i =
let rec aux pos value i =
if i > n then
@ -109,9 +103,7 @@ let pivoted_ldl threshold m_A =
let p_i, p_j = pi.(i), pi.(j) in
pi.(i) <- p_j;
pi.(j) <- p_i;
in
in
let () =
@ -132,11 +124,7 @@ let pivoted_ldl threshold m_A =
m_Lt, v_D, pi
let make_ldl ?(threshold=Constants.epsilon) m_A =
pivoted_ldl m_A

53
Utils/Coordinate.ml
View File

@ -1,19 +1,7 @@
type bohr
type angstrom
type bohr = Bohr.t
type angstrom = Angstrom.t
type xyz = {
x : float ;
y : float ;
z : float ;
}
type 'a point = xyz
type t = bohr point
external make : 'a point -> t = "%identity"
external make_angstrom : 'a point -> angstrom point = "%identity"
type t = bohr
type axis = X | Y | Z
@ -21,17 +9,19 @@ type axis = X | Y | Z
let a_to_b a = Constants.a0_inv *. a
let b_to_a b = Constants.a0 *. b
let bohr_to_angstrom { x ; y ; z } =
make
let bohr_to_angstrom { Bohr.x ; y ; z } =
Angstrom.make
{
Point.
x = b_to_a x ;
y = b_to_a y ;
z = b_to_a z ;
}
let angstrom_to_bohr { x ; y ; z } =
make
let angstrom_to_bohr { Angstrom.x ; y ; z } =
Bohr.make
{
Point.
x = a_to_b x ;
y = a_to_b y ;
z = a_to_b z ;
@ -39,27 +29,30 @@ let angstrom_to_bohr { x ; y ; z } =
let zero =
make { x = 0. ; y = 0. ; z = 0. }
Bohr.make { Point.x = 0. ; y = 0. ; z = 0. }
(** Linear algebra *)
let ( |. ) s { x ; y ; z } =
make {
let ( |. ) s { Bohr.x ; y ; z } =
Bohr.make {
Point.
x = s *. x ;
y = s *. y ;
z = s *. z ;
}
let ( |+ ) { x = x1 ; y = y1 ; z = z1 } { x = x2 ; y = y2 ; z = z2 } =
make {
let ( |+ ) { Bohr.x = x1 ; y = y1 ; z = z1 } { Bohr.x = x2 ; y = y2 ; z = z2 } =
Bohr.make {
Point.
x = x1 +. x2 ;
y = y1 +. y2 ;
z = z1 +. z2 ;
}
let ( |- ) { x = x1 ; y = y1 ; z = z1 } { x = x2 ; y = y2 ; z = z2 } =
make {
let ( |- ) { Bohr.x = x1 ; y = y1 ; z = z1 } { Bohr.x = x2 ; y = y2 ; z = z2 } =
Bohr.make {
Point.
x = x1 -. x2 ;
y = y1 -. y2 ;
z = z1 -. z2 ;
@ -69,7 +62,7 @@ let ( |- ) { x = x1 ; y = y1 ; z = z1 } { x = x2 ; y = y2 ; z = z2 } =
let neg a = -1. |. a
let dot { x = x1 ; y = y1 ; z = z1 } { x = x2 ; y = y2 ; z = z2 } =
let dot { Bohr.x = x1 ; y = y1 ; z = z1 } { Bohr.x = x2 ; y = y2 ; z = z2 } =
x1 *. x2 +. y1 *. y2 +. z1 *. z2
@ -77,7 +70,7 @@ let norm u =
sqrt ( dot u u )
let get axis { x ; y ; z } =
let get axis { Bohr.x ; y ; z } =
match axis with
| X -> x
| Y -> y
@ -86,12 +79,14 @@ let get axis { x ; y ; z } =
open Format
let pp ppf c =
let open Bohr in
fprintf ppf "@[@[%8.4f@] @[%8.4f@] @[%8.4f@]@]" c.x c.y c.z
let pp_bohr ppf c =
let open Bohr in
fprintf ppf "@[(@[%10f@], @[%10f@], @[%10f@] Bohr)@]" c.x c.y c.z
let pp_angstrom ppf c =
let c = bohr_to_angstrom c in
let open Angstrom in
fprintf ppf "@[(@[%10f@], @[%10f@], @[%10f@] Angs)@]" c.x c.y c.z

70
Utils/Coordinate.mli
View File

@ -1,108 +1,92 @@
(** Type for coordinates in 3D space.
All operations on points are done in atomic units. Therefore,
all coordinates are given in bohr and manipulated with this
all coordinates are given in {!Bohr.t} and manipulated with this
module.
*)
type bohr
type angstrom
type xyz = {
x : float ;
y : float ;
z : float ;
}
type 'a point = xyz
type t = bohr point
val make : 'a point -> t
(** Creates a point in atomic units *)
val make_angstrom : 'a point -> angstrom point
(** Creates a point in Angstrom *)
type t = Bohr.t
type axis = X | Y | Z
val bohr_to_angstrom : bohr point -> angstrom point
(** Converts a point in bohr to angstrom. *)
val bohr_to_angstrom : Bohr.t -> Angstrom.t
(** Converts a point in Bohr to Angstrom. *)
val angstrom_to_bohr : angstrom point -> bohr point
(** Converts a point in Angstrom to bohr. *)
val angstrom_to_bohr : Angstrom.t -> Bohr.t
(** Converts a point in Angstrom to Bohr. *)
val zero : bohr point
(** [zero = { x = 0. ; y=0. ; z=0. }] *)
val zero : Bohr.t
(** [zero = { Bohr.x = 0. ; y=0. ; z=0. }] *)
val get : axis -> bohr point -> float
val get : axis -> Bohr.t -> float
(** Extracts the projection of the coordinate on an axis.
Example:
[Coordinate.(get Y) { x=1. ; y=2. ; z=3. } -> 2.]
[Coordinate.(get Y) { Bohr.x=1. ; y=2. ; z=3. } -> 2.]
*)
(** {1 Vector operations} *)
val ( |. ) : float -> t -> t
val ( |. ) : float -> Bohr.t -> Bohr.t
(** Scale by a float.
Example:
[ 2. |. { x=1. ; y=2. ; z=3. } -> { x=2. ; y=4. ; z=6. } ]
[ 2. |. { Bohr.x=1. ; y=2. ; z=3. } -> { Bohr.x=2. ; y=4. ; z=6. } ]
*)
val ( |+ ) : t -> t -> t
val ( |+ ) : Bohr.t -> Bohr.t -> Bohr.t
(** Add two vectors.
Example:
{[{ x=1. ; y=2. ; z=3. } |+ { x=2. ; y=3. ; z=1. } ->
{ x=3. ; y=5. ; z=4. }]}
{[{ Bohr.x=1. ; y=2. ; z=3. } |+ { Bohr.x=2. ; y=3. ; z=1. } ->
{ Bohr.x=3. ; y=5. ; z=4. }]}
*)
val ( |- ) : t -> t -> t
val ( |- ) : Bohr.t -> Bohr.t -> Bohr.t
(** Subtract two vectors.
Example:
{[{ x=1. ; y=2. ; z=3. } |- { x=2. ; y=3. ; z=1. } ->
{ x=-1. ; y=-1. ; z=2. }]}
{[{ Bohr.x=1. ; y=2. ; z=3. } |- { Bohr.x=2. ; y=3. ; z=1. } ->
{ Bohr.x=-1. ; y=-1. ; z=2. }]}
*)
val neg : t -> t
val neg : Bohr.t -> Bohr.t
(** Example:
{[Coordinate.neg { x=1. ; y=2. ; z=-3. } ->
{ x=-1. ; y=-2. ; z=3. }]}
{[Coordinate.neg { Bohr.x=1. ; y=2. ; z=-3. } ->
{ Bohr.x=-1. ; y=-2. ; z=3. }]}
*)
val dot : t -> t -> float
val dot : Bohr.t -> Bohr.t -> float
(** Dot product. *)
val norm : t -> float
val norm : Bohr.t -> float
(** L{^2} norm of the vector. *)
(** {2 Printers} *)
val pp: Format.formatter -> t -> unit
val pp: Format.formatter -> Bohr.t -> unit
val pp_bohr: Format.formatter -> t -> unit
val pp_bohr: Format.formatter -> Bohr.t -> unit
val pp_angstrom : Format.formatter -> t -> unit
val pp_angstrom : Format.formatter -> Bohr.t -> unit