10
1
mirror of https://gitlab.com/scemama/QCaml.git synced 2024-07-25 04:07:24 +02:00

1d arrays

This commit is contained in:
Anthony Scemama 2018-02-14 20:00:12 +01:00
parent 5161cd8226
commit 7b8be87bfc
2 changed files with 84 additions and 40 deletions

View File

@ -78,21 +78,25 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
match v1_top with
| None -> ()
| Some v0 ->
for l=0 to np-1 do
let f0 =
-. expo_b.(l) *. expo_inv_p.(l) *. cab
Array.iteri (fun l result_l ->
let f0 = -. expo_b.(l) *. expo_inv_p.(l) *. cab
and v0_l = v0.(l)
and result_l = result.(l)
in
for k=0 to nq-1 do
result.(l).(k) <- v0.(l).(k) *. f0
done
done
Array.iteri (fun k v0_lk ->
result_l.(k) <- v0_lk *. f0) v0_l
) result
end;
for l=0 to np-1 do
for k=0 to nq-1 do
result.(l).(k) <- result.(l).(k)
+. expo_inv_p.(l) *. (center_pq xyz).(l).(k) *. p0.(l).(k)
done
done
Array.iteri (fun l result_l ->
let expo_inv_p_l = expo_inv_p.(l)
and center_pq_xyz_l = (center_pq xyz).(l)
and result_l = result.(l)
and p0_l = p0.(l)
in
Array.iteri (fun k p0_lk ->
result_l.(k) <- result_l.(k)
+. expo_inv_p_l *. center_pq_xyz_l.(k) *. p0_lk
) p0_l ) result
end
else
begin
@ -105,12 +109,12 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
match v1_top with
| None -> ()
| Some v0 ->
for l=0 to np-1 do
let f0 = -. expo_b.(l) *. expo_inv_p.(l) *. cab in
for k=0 to nq-1 do
result.(l).(k) <- v0.(l).(k) *. f0
done
done
Array.iteri (fun l result_l ->
let f0 = -. expo_b.(l) *. expo_inv_p.(l) *. cab
and v0_l = v0.(l)
in
Array.iteri (fun k v0_lk ->
result_l.(k) <- v0_lk *. f0) v0_l ) result
end;
let v1 =
match v1_top2 with
@ -122,15 +126,20 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
| Some p1_top2 -> p1_top2
| None -> assert false
in
for l=0 to np-1 do
let f = (float_of_int amxyz) *. expo_inv_p.(l) *. 0.5 in
for k=0 to nq-1
do
result.(l).(k) <- result.(l).(k) +.
expo_inv_p.(l) *. (center_pq xyz).(l).(k) *. p0.(l).(k) +.
f *. (v1.(l).(k) +. v2.(l).(k) *. expo_inv_p.(l))
done
done
Array.iteri (fun l result_l ->
let f = float_of_int amxyz *. expo_inv_p.(l) *. 0.5
and expo_inv_p_l = expo_inv_p.(l)
and center_pq_xyz_l = (center_pq xyz).(l)
and v1_l = v1.(l)
and v2_l = v2.(l)
and result_l = result.(l)
in
Array.iteri (fun k p0_lk ->
result_l.(k) <- result_l.(k) +.
expo_inv_p_l *. center_pq_xyz_l.(k) *. p0_lk +.
f *. (v1_l.(k) +. v2_l.(k) *. expo_inv_p_l)
) p0.(l)
) result
end;
Some result
end
@ -195,7 +204,11 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
result.(l).(k) <- v1.(l).(k) *. f1.(k)
done
done;
Some result
Some (
Array.init np (fun l ->
let v1_l = v1.(l) in
Array.init nq (fun k -> v1_l.(k) *. f1.(k))
))
end
else None
in
@ -203,8 +216,9 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
let v2 =
let f2 =
Array.init np (fun l ->
let cpq_l = (center_pq xyz).(l) in
Array.init nq (fun k ->
let x = expo_inv_q.(k) *. (center_pq xyz).(l).(k) in
let x = expo_inv_q.(k) *. cpq_l.(k) in
if (!do_compute) then x
else (if abs_float x > cutoff then do_compute := true ; x)
) )
@ -215,8 +229,11 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
| Some v2 ->
begin
for l=0 to np-1 do
let f2_l = f2.(l)
and v2_l = v2.(l)
in
for k=0 to nq-1 do
f2.(l).(k) <- -. v2.(l).(k) *. f2.(l).(k)
f2_l.(k) <- -. v2_l.(k) *. f2_l.(k)
done
done;
Some f2
@ -233,8 +250,11 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
| Some v1, Some v2 ->
begin
for l=0 to np-1 do
let v1_l = v1.(l)
and v2_l = v2.(l)
in
for k=0 to nq-1 do
v2.(l).(k) <- v2.(l).(k) +. v1.(l).(k)
v2_l.(k) <- v2_l.(k) +. v1_l.(k)
done
done;
Some v2
@ -259,8 +279,11 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
begin
let result = Array.make_matrix np nq 0. in
for l=0 to np-1 do
let v1_l = v1.(l)
and result_l = result.(l)
in
for k=0 to nq-1 do
result.(l).(k) <- v1.(l).(k) *. f1.(k)
result_l.(k) <- v1_l.(k) *. f1.(k)
done;
done;
Some result
@ -283,8 +306,11 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
begin
let result = Array.make_matrix np nq 0. in
for l=0 to np-1 do
let v3_l = v3.(l)
and result_l = result.(l)
in
for k=0 to nq-1 do
result.(l).(k) <- v3.(l).(k) *. f2.(k)
result_l.(k) <- v3_l.(k) *. f2.(k)
done
done;
Some result
@ -299,8 +325,12 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
| Some p1, Some v1, Some v3 ->
begin
for l=0 to np-1 do
let v3_l = v3.(l)
and v1_l = v1.(l)
and p1_l = p1.(l)
in
for k=0 to nq-1 do
v3.(l).(k) <- p1.(l).(k) +. v1.(l).(k) +. v3.(l).(k)
v3_l.(k) <- p1_l.(k) +. v1_l.(k) +. v3_l.(k)
done
done;
Some v3
@ -308,8 +338,11 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
| Some p1, Some v1, None ->
begin
for l=0 to np-1 do
let v1_l = v1.(l)
and p1_l = p1.(l)
in
for k=0 to nq-1 do
p1.(l).(k) <- v1.(l).(k) +. p1.(l).(k)
p1_l.(k) <- v1_l.(k) +. p1_l.(k)
done
done;
Some p1
@ -317,8 +350,11 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
| Some p1, None, Some v3 ->
begin
for l=0 to np-1 do
let v3_l = v3.(l)
and p1_l = p1.(l)
in
for k=0 to nq-1 do
p1.(l).(k) <- p1.(l).(k) +. v3.(l).(k)
p1_l.(k) <- p1_l.(k) +. v3_l.(k)
done
done;
Some p1
@ -326,8 +362,11 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
| None , Some v1, Some v3 ->
begin
for l=0 to np-1 do
let v3_l = v3.(l)
and v1_l = v1.(l)
in
for k=0 to nq-1 do
v3.(l).(k) <- v1.(l).(k) +. v3.(l).(k)
v3_l.(k) <- v1_l.(k) +. v3_l.(k)
done
done;
Some v3
@ -349,8 +388,11 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d)
in
for l=0 to np-1 do
let fa = (float_of_int axyz) *. expo_inv_p.(l) *. 0.5 in
let p2_l = p2.(l)
and v_l = v.(l)
in
for k=0 to nq-1 do
p2.(l).(k) <- p2.(l).(k) -. fa *. expo_inv_q.(k) *. v.(l).(k)
p2_l.(k) <- p2_l.(k) -. fa *. expo_inv_q.(k) *. v_l.(k)
done
done;
Some p2

View File

@ -1,4 +1,4 @@
(** Key for hastables that contain tuples of integers encoded in a Zarith integer *)
(** Key for hastables that contain tuples of integers encoded in small integers *)
type kind_array =
| Kind_3
| Kind_6
@ -218,6 +218,7 @@ let to_string ~kind { left ; right } =
|> String.concat ", "
) ^ " >"
(*
let debug () =
let k2 = of_int_array Kind_2 [| 1 ; 2 |]
and k3 = of_int_array Kind_3 [| 1 ; 2 ; 3 |]
@ -240,3 +241,4 @@ let debug () =
*)
let () = debug ()
*)