10
1
mirror of https://gitlab.com/scemama/QCaml.git synced 2024-06-19 19:52:06 +02:00
This commit is contained in:
Anthony Scemama 2018-01-19 23:31:10 +01:00
parent bba7b6e8e4
commit b92eaf9fe5

View File

@ -12,12 +12,8 @@ open Util
and $\phi_q$
*)
let zero_m ~maxm ~expo_pq_inv ~norm_pq_sq =
let exp_pq =
1. /. expo_pq_inv
in
let t =
norm_pq_sq *. exp_pq
in
let exp_pq = 1. /. expo_pq_inv in
let t = norm_pq_sq *. exp_pq in
boys_function ~maxm t
|> Array.mapi (fun m fm ->
two_over_sq_pi *. (if m mod 2 = 0 then fm else -.fm) *.
@ -34,8 +30,14 @@ type n_cls = { n : int ; cls : Z.t array }
(** Write all integrals to a file with the <ij|kl> convention *)
let to_file ~filename basis =
let to_int_tuple x =
let open Zkey in
match to_int_tuple Kind_3 x with
| Three x -> x
| _ -> assert false
in
let oc = open_out filename in
let key_array = Array.make 12 0 in
for i=0 to (Array.length basis) - 1 do
print_int basis.(i).Contracted_shell.indice ; print_newline ();
for j=0 to i do
@ -46,48 +48,36 @@ let to_file ~filename basis =
contracted_class basis.(i) basis.(j) basis.(k) basis.(l)
in
for i_c = 0 to (Array.length basis.(i).Contracted_shell.powers) - 1 do
let x = Zkey.(to_int_array Kind_3 basis.(i).Contracted_shell.powers.(i_c)) in
key_array.(0) <- x.(0);
key_array.(1) <- x.(1);
key_array.(2) <- x.(2);
for j_c = 0 to (Array.length basis.(j).Contracted_shell.powers) - 1 do
let x = Zkey.(to_int_array Kind_3 basis.(j).Contracted_shell.powers.(j_c)) in
key_array.(3) <- x.(0);
key_array.(4) <- x.(1);
key_array.(5) <- x.(2);
for k_c = 0 to (Array.length basis.(k).Contracted_shell.powers) - 1 do
let x = Zkey.(to_int_array Kind_3 basis.(k).Contracted_shell.powers.(k_c)) in
key_array.(6) <- x.(0);
key_array.(7) <- x.(1);
key_array.(8) <- x.(2);
for l_c = 0 to (Array.length basis.(l).Contracted_shell.powers) - 1 do
let x = Zkey.(to_int_array Kind_3 basis.(l).Contracted_shell.powers.(l_c)) in
key_array.( 9) <- x.(0);
key_array.(10) <- x.(1);
key_array.(11) <- x.(2);
(* Write the data in the output file *)
Array.iteri (fun i_c powers_i ->
let i_c = basis.(i).Contracted_shell.indice + i_c + 1 in
let xi = to_int_tuple powers_i in
Array.iteri (fun j_c powers_j ->
let j_c = basis.(j).Contracted_shell.indice + j_c + 1 in
let xj = to_int_tuple powers_j in
Array.iteri (fun k_c powers_k ->
let k_c = basis.(k).Contracted_shell.indice + k_c + 1 in
let xk = to_int_tuple powers_k in
Array.iteri (fun l_c powers_l ->
let l_c = basis.(l).Contracted_shell.indice + l_c + 1 in
let xl = to_int_tuple powers_l in
let key =
Zkey.(of_int_array Kind_12 key_array)
Zkey.of_int_tuple (Zkey.Twelve (xi,xj,xk,xl))
in
let value =
Zmap.find cls key
in
if (abs_float value > cutoff) then
Printf.fprintf oc "%4d %4d %4d %4d %20.12e\n"
(basis.(i).Contracted_shell.indice + i_c + 1)
(basis.(k).Contracted_shell.indice + k_c + 1)
(basis.(j).Contracted_shell.indice + j_c + 1)
(basis.(l).Contracted_shell.indice + l_c + 1)
value
done
done
done
i_c k_c j_c l_c value
) basis.(l).Contracted_shell.powers
) basis.(k).Contracted_shell.powers
) basis.(j).Contracted_shell.powers
) basis.(i).Contracted_shell.powers;
done;
done;
done;
done;
done
;
close_out oc
(*