diff --git a/Basis/F12.ml b/Basis/F12.ml index 22ffaa2..4f6e6e5 100644 --- a/Basis/F12.ml +++ b/Basis/F12.ml @@ -31,12 +31,24 @@ let make_gaussian_corr_factor expo_s coef_g expo_sg = } -(* exp (-expo_s r) *) +(* -1/expo_s *. exp (-expo_s r) *) let gaussian_geminal expo_s = let coef_g = - [| 0.3144 ; 0.3037 ; 0.1681 ; 0.09811 ; 0.06024 ; 0.03726 |]; + [| 0.3144 ; 0.3037 ; 0.1681 ; 0.09811 ; 0.06024 ; 0.03726 |] + |> Array.map (fun x -> -. x /. expo_s) and expo_sg = - [| 0.2209 ; 1.004 ; 3.622 ; 12.16 ; 45.87 ; 254.4 |]; + [| 0.2209 ; 1.004 ; 3.622 ; 12.16 ; 45.87 ; 254.4 |] + in + make_gaussian_corr_factor expo_s coef_g expo_sg + + + +(* exp (-expo_s r) *) +let simple_gaussian_geminal expo_s = + let coef_g = + [| 0.3144 ; 0.3037 ; 0.1681 ; 0.09811 ; 0.06024 ; 0.03726 |] + and expo_sg = + [| 0.2209 ; 1.004 ; 3.622 ; 12.16 ; 45.87 ; 254.4 |] in make_gaussian_corr_factor expo_s coef_g expo_sg @@ -52,7 +64,7 @@ let gaussian_geminal_times_r12 expo_s = (* exp (-expo_s r) *) -let gaussian_geminal' expo_s = +let simple_gaussian_geminal' expo_s = let coef_g = [| -3.4793465193721626604883567779324948787689208984375 ; @@ -135,7 +147,7 @@ let filter_contracted_shell_pairs f12 ?(cutoff=integrals_cutoff) shell_pairs = *) -let store_class basis ?(cutoff=integrals_cutoff) data contracted_shell_pair_couple cls = +let store_class basis ?(cutoff=integrals_cutoff) data f12 contracted_shell_pair_couple cls = let to_powers x = let open Zkey in match to_powers x with @@ -148,11 +160,7 @@ let store_class basis ?(cutoff=integrals_cutoff) data contracted_shell_pair_coup in (* - let s = - Overlap.of_basis basis - |> Overlap.matrix - in - let lambda_inv = 1. /. f12.expo_s in + let lambda_inv = -. 1. /. f12.expo_s in *) Array.iteri (fun i_c powers_i -> let i_c = Cs.index (Csp.shell_a shell_p) + i_c + 1 in @@ -170,6 +178,7 @@ let store_class basis ?(cutoff=integrals_cutoff) data contracted_shell_pair_coup let value = Zmap.find cls key in (* lambda_inv *. (s.{i_c,j_c} *. s.{k_c,l_c} -. value) + lambda_inv *. value *) value |> set_chem data i_c j_c k_c l_c @@ -240,7 +249,7 @@ let of_basis_serial f12 basis = let cls = class_of_contracted_shell_pair_couple f12 cspc in - store_class basis ~cutoff eri_array cspc cls + store_class basis ~cutoff eri_array f12 cspc cls | None -> () ) shell_pairs with Exit -> () @@ -364,18 +373,15 @@ let of_basis_parallel f12 basis = Fis.create ~size:0 `Dense in - (* - let s = - Overlap.of_basis basis - |> Overlap.matrix - in - let lambda_inv = 1. /. f12.expo_s in + (* + let lambda_inv = -. 1. /. f12.expo_s in *) Farm.run ~ordered:false ~f input_stream |> Stream.iter (fun l -> Array.iter (fun (i_c,j_c,k_c,l_c,value) -> (* lambda_inv *. (s.{i_c,j_c} *. s.{k_c,l_c} -. value) + lambda_inv *. value *) value |> set_chem eri_array i_c j_c k_c l_c) l); diff --git a/CI/F12CI.ml b/CI/F12CI.ml index 772cb49..b4f6704 100644 --- a/CI/F12CI.ml +++ b/CI/F12CI.ml @@ -244,6 +244,11 @@ let make ~simulation ?(threshold=1.e-12) ~frozen_core ~mo_basis ~aux_basis_filen let f12_amplitudes = + fun c -> + let result = lacpy c in + Mat.scal (0.5) result ; + result + (* (* While in a sequential region, initiate the parallel 4-idx transformation to avoid nested parallel jobs *) @@ -262,8 +267,16 @@ let make ~simulation ?(threshold=1.e-12) ~frozen_core ~mo_basis ~aux_basis_filen |> Lazy.force in fun ci_coef -> - Matrix.ax_eq_b m_F (Matrix.dense_of_mat ci_coef) - |> Matrix.to_mat + let result = + Matrix.ax_eq_b m_F (Matrix.dense_of_mat ci_coef) + |> Matrix.to_mat + in + let norm = gemm ~transa:`T result result in + Printf.printf "Norm of |F> : %f\n%!" norm.{1,1}; + if abs_float norm.{1,1} > 1. then + failwith "Norm of |F> > 1"; + result + *) in let e_shift = diff --git a/Utils/Bitstring.ml b/Utils/Bitstring.ml index 3e6382c..ff86504 100644 --- a/Utils/Bitstring.ml +++ b/Utils/Bitstring.ml @@ -189,7 +189,7 @@ let permtutations m n = let t'' = shift_right (minus_one (logand (lognot t) t')) (trailing_zeros u + 1) in aux (k-1) (logor t' t'') (u :: rest) in - aux (Util.binom n m) (minus_one (shift_left_one m m)) [] + aux (Util.binom n m) (minus_one (shift_left_one n m)) [] (*-----------------------------------------------------------------------------------*)