From 56839c3332576d12add8795cedb38201a8c72312 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Tue, 12 Mar 2019 14:07:59 +0100 Subject: [PATCH] Added Zero_m_parameters --- Basis/ERI.ml | 20 +++++++++++++------- Basis/TwoElectronIntegrals.ml | 3 +-- Basis/TwoElectronIntegrals.mli | 3 +-- Basis/TwoElectronRR.ml | 9 +++++++-- Basis/TwoElectronRRVectorized.ml | 9 +++++++-- Basis/Zero_m_parameters.ml | 15 +++++++++++++++ 6 files changed, 44 insertions(+), 15 deletions(-) create mode 100644 Basis/Zero_m_parameters.ml diff --git a/Basis/ERI.ml b/Basis/ERI.ml index eff125f..1237a38 100644 --- a/Basis/ERI.ml +++ b/Basis/ERI.ml @@ -3,19 +3,24 @@ open Constants open Util + + module Zm = struct let name = "Electron repulsion integrals" - let zero_m ~maxm ~expo_p_inv ~expo_q_inv ~norm_pq_sq = - let expo_pq_inv = expo_p_inv +. expo_q_inv in + open Zero_m_parameters + + let zero_m z = + let expo_pq_inv = z.expo_p_inv +. z.expo_q_inv in assert (expo_pq_inv <> 0.); - let norm_pq_sq = - if norm_pq_sq > integrals_cutoff then norm_pq_sq else 0. - in let exp_pq = 1. /. expo_pq_inv in - let t = norm_pq_sq *. exp_pq in - let f = two_over_sq_pi *. (sqrt exp_pq) in + let t = + if z.norm_pq_sq > integrals_cutoff then + z.norm_pq_sq *. exp_pq + else 0. + in + let maxm = z.maxm in let result = boys_function ~maxm t in let rec aux accu k = function | 0 -> result.(k) <- result.(k) *. accu @@ -26,6 +31,7 @@ module Zm = struct aux new_accu (k+1) (l-1) end in + let f = two_over_sq_pi *. (sqrt exp_pq) in aux f 0 maxm; result diff --git a/Basis/TwoElectronIntegrals.ml b/Basis/TwoElectronIntegrals.ml index ade4e5f..786e1b8 100644 --- a/Basis/TwoElectronIntegrals.ml +++ b/Basis/TwoElectronIntegrals.ml @@ -17,8 +17,7 @@ module Fis = FourIdxStorage module type Zero_mType = sig val name : string - val zero_m : maxm:int -> expo_p_inv:float -> expo_q_inv:float -> - norm_pq_sq:float -> float array + val zero_m : Zero_m_parameters.t -> float array end diff --git a/Basis/TwoElectronIntegrals.mli b/Basis/TwoElectronIntegrals.mli index fc827b8..c95761b 100644 --- a/Basis/TwoElectronIntegrals.mli +++ b/Basis/TwoElectronIntegrals.mli @@ -11,8 +11,7 @@ module type Zero_mType = val name : string (** Name of the kind of integrals, for printing purposes. *) -val zero_m : maxm:int -> expo_p_inv:float -> expo_q_inv:float -> - norm_pq_sq:float -> float array +val zero_m : Zero_m_parameters.t -> float array (** The returned float array contains all the {% $(00|00)^m$ %} values, where [m] is the index of the array. diff --git a/Basis/TwoElectronRR.ml b/Basis/TwoElectronRR.ml index fb626d5..387e041 100644 --- a/Basis/TwoElectronRR.ml +++ b/Basis/TwoElectronRR.ml @@ -11,6 +11,7 @@ module Po = Powers module Psp = PrimitiveShellPair module Pspc = PrimitiveShellPairCouple module Ps = PrimitiveShell +module Zp = Zero_m_parameters let cutoff = Constants.integrals_cutoff let cutoff2 = cutoff *. cutoff @@ -340,7 +341,9 @@ let contracted_class_shell_pair_couple ~zero_m shell_pair_couple : float Zmap.t let expo_q_inv = Psp.exponent_inv sp_cd in let zero_m_array = - zero_m ~maxm ~expo_p_inv ~expo_q_inv ~norm_pq_sq + Zp.{ + maxm ; expo_p_inv ; expo_q_inv ; norm_pq_sq + } |> zero_m in begin @@ -451,7 +454,9 @@ let contracted_class_atomic_shell_pair_couple ~zero_m atomic_shell_pair_couple : let expo_q_inv = Psp.exponent_inv sp_cd in let zero_m_array = - zero_m ~maxm ~expo_p_inv ~expo_q_inv ~norm_pq_sq + Zp.{ + maxm ; expo_p_inv ; expo_q_inv ; norm_pq_sq + } |> zero_m in begin diff --git a/Basis/TwoElectronRRVectorized.ml b/Basis/TwoElectronRRVectorized.ml index f232cbc..a5856c5 100644 --- a/Basis/TwoElectronRRVectorized.ml +++ b/Basis/TwoElectronRRVectorized.ml @@ -10,6 +10,7 @@ module Cspc = ContractedShellPairCouple module Po = Powers module Psp = PrimitiveShellPair module Ps = PrimitiveShell +module Zp = Zero_m_parameters exception NullQuartet exception Found @@ -648,7 +649,9 @@ let contracted_class_shell_pairs ~zero_m ?schwartz_p ?schwartz_q shell_p shell_q in let zero_m_array = - zero_m ~maxm:0 ~expo_p_inv ~expo_q_inv ~norm_pq_sq + Zp.{ + maxm=0 ; expo_p_inv ; expo_q_inv ; norm_pq_sq + } |> zero_m in zero_m_array.(0) with NullQuartet -> 0. @@ -761,7 +764,9 @@ let contracted_class_shell_pairs ~zero_m ?schwartz_p ?schwartz_q shell_p shell_q x *. x +. y *. y +. z *. z in - zero_m ~maxm ~expo_p_inv ~expo_q_inv ~norm_pq_sq + Zp.{ + maxm ; expo_p_inv ; expo_q_inv ; norm_pq_sq + } |> zero_m ) sq in (* Transpose result *) diff --git a/Basis/Zero_m_parameters.ml b/Basis/Zero_m_parameters.ml new file mode 100644 index 0000000..efb6131 --- /dev/null +++ b/Basis/Zero_m_parameters.ml @@ -0,0 +1,15 @@ +type t = +{ + expo_p_inv : float ; + expo_q_inv : float ; + norm_pq_sq : float ; + maxm : int; +} + +let zero = +{ + maxm=0 ; + expo_p_inv = 0.; + expo_q_inv = 0.; + norm_pq_sq = 0.; +}