diff --git a/Basis/TwoElectronRR.ml b/Basis/TwoElectronRR.ml index a0c76c6..06fef53 100644 --- a/Basis/TwoElectronRR.ml +++ b/Basis/TwoElectronRR.ml @@ -82,7 +82,7 @@ let hvrr_two_e (angMom_a, angMom_b, angMom_c, angMom_d) try Zmap.find map_2d key with | Not_found -> let result = - let am, cm, cmm, axyz, cxyz, xyz = + let am, cm, cmm, axyz, cmxyz, xyz = let angMom_ax, angMom_ay, angMom_az = angMom_a and angMom_cx, angMom_cy, angMom_cz = angMom_c in match angMom_c with @@ -90,19 +90,19 @@ let hvrr_two_e (angMom_a, angMom_b, angMom_c, angMom_d) (angMom_ax-1, angMom_ay, angMom_az), (angMom_cx-1, angMom_cy, angMom_cz), (angMom_cx-2, angMom_cy, angMom_cz), - angMom_ax,angMom_cx, 0 + angMom_ax,angMom_cx-1, 0 | (_,_,0) -> (* 612_002 *) (angMom_ax, angMom_ay-1, angMom_az), (angMom_cx, angMom_cy-1, angMom_cz), (angMom_cx, angMom_cy-2, angMom_cz), - angMom_ay,angMom_cy, 1 + angMom_ay,angMom_cy-1, 1 | _ -> (* 1_067_324 *) (angMom_ax, angMom_ay, angMom_az-1), (angMom_cx, angMom_cy, angMom_cz-1), (angMom_cx, angMom_cy, angMom_cz-2), - angMom_az,angMom_cz, 2 + angMom_az,angMom_cz-1, 2 in - if cxyz < 1 then empty else + if cmxyz < 0 then empty else let f1 = -. expo_d *. expo_inv_q *. (Coordinate.coord center_cd xyz) in @@ -118,9 +118,9 @@ let hvrr_two_e (angMom_a, angMom_b, angMom_c, angMom_d) (* 9_232_029 *) f1 *. v1.(m) -. (if m = maxm then 0. else f2 *. v1.(m+1)) ) in let result = - if cxyz < 2 then result else + if cmxyz < 1 then result else let f3 = - (float_of_int (cxyz-1)) *. expo_inv_q *. 0.5 + (float_of_int cmxyz) *. expo_inv_q *. 0.5 in if (abs_float f3 < cutoff) && (abs_float (f3 *. expo_inv_q) < cutoff) then result else ( @@ -132,7 +132,7 @@ let hvrr_two_e (angMom_a, angMom_b, angMom_c, angMom_d) ) in let result = - if axyz < 1 || cxyz < 1 then result else + if (axyz < 1) || (cmxyz < 0) then result else let f5 = (float_of_int axyz) *. expo_inv_p *. expo_inv_q *. 0.5 in @@ -158,34 +158,44 @@ let hvrr_two_e (angMom_a, angMom_b, angMom_c, angMom_d) (* 8_448_486 *) match totAngMom_b with - | 0 -> (* 0 *) (vrr (angMom_a.(0) ,angMom_a.(1) ,angMom_a.(2)) (angMom_c.(0) ,angMom_c.(1) ,angMom_c.(2)) totAngMom_a totAngMom_c).(0) + | 0 -> (* 0 *) (vrr angMom_a angMom_c totAngMom_a totAngMom_c).(0) | 1 -> - (* 5_045_008 *) let xyz = if angMom_b.(0) = 1 then 0 else if angMom_b.(1) = 1 then 1 else 2 in - let ap = [| angMom_a.(0) ; angMom_a.(1) ; angMom_a.(2) |] in - ap.(xyz) <- ap.(xyz) + 1; + let angMom_ax, angMom_ay, angMom_az = angMom_a in + (* 5_045_008 *) + let ap, xyz = + match angMom_b with + | (1,_,_) -> (angMom_ax+1,angMom_ay,angMom_az), 0 + | (_,1,_) -> (angMom_ax,angMom_ay+1,angMom_az), 1 + | (_,_,_) -> (angMom_ax,angMom_ay,angMom_az+1), 2 + in let v1 = - vrr (ap.(0),ap.(1),ap.(2)) (angMom_c.(0) ,angMom_c.(1) ,angMom_c.(2)) (totAngMom_a+1) totAngMom_c + vrr ap angMom_c (totAngMom_a+1) totAngMom_c in let f2 = (Coordinate.coord center_ab xyz) in if (abs_float f2 < cutoff) then v1.(0) else let v2 = - vrr (angMom_a.(0) ,angMom_a.(1) ,angMom_a.(2)) (angMom_c.(0) ,angMom_c.(1) ,angMom_c.(2)) totAngMom_a totAngMom_c + vrr angMom_a angMom_c totAngMom_a totAngMom_c in v1.(0) +. f2 *. v2.(0) | _ -> - (* 3_403_478 *) let ap = [| angMom_a.(0) ; angMom_a.(1) ; angMom_a.(2) |] - and bm = [| angMom_b.(0) ; angMom_b.(1) ; angMom_b.(2) |] - and xyz = - match angMom_b with - | [|_;0;0|] -> (* 677_315 *) 0 - | [|_;_;0|] -> (* 1_136_646 *) 1 - | _ -> (* 1_589_517 *) 2 + let angMom_ax, angMom_ay, angMom_az = angMom_a + and angMom_bx, angMom_by, angMom_bz = angMom_b in + (* 3_403_478 *) + let bxyz, xyz = + match angMom_b with + | (_,0,0) -> angMom_bx, 0 + | (_,_,0) -> angMom_by, 1 + | (_,_,_) -> angMom_bz, 2 in - ap.(xyz) <- ap.(xyz) + 1; - bm.(xyz) <- bm.(xyz) - 1; - if (bm.(xyz) < 0) then 0. else + if (bxyz < 1) then 0. else + let ap, bm = + match xyz with + | 0 -> (angMom_ax+1,angMom_ay,angMom_az),(angMom_bx-1,angMom_by,angMom_bz) + | 1 -> (angMom_ax,angMom_ay+1,angMom_az),(angMom_bx,angMom_by-1,angMom_bz) + | _ -> (angMom_ax,angMom_ay,angMom_az+1),(angMom_bx,angMom_by,angMom_bz-1) + in let h1 = hrr0 ap bm angMom_c (totAngMom_a+1) (totAngMom_b-1) totAngMom_c in @@ -205,20 +215,19 @@ let hvrr_two_e (angMom_a, angMom_b, angMom_c, angMom_d) match (totAngMom_b, totAngMom_d) with | (_,0) -> (* 3_608_781 *) if (totAngMom_b = 0) then - (vrr (angMom_a.(0) ,angMom_a.(1) ,angMom_a.(2)) (angMom_c.(0) ,angMom_c.(1) ,angMom_c.(2)) totAngMom_a totAngMom_c).(0) + (vrr angMom_a angMom_c totAngMom_a totAngMom_c).(0) else hrr0 angMom_a angMom_b angMom_c totAngMom_a totAngMom_b totAngMom_c | (_,_) -> - (* 4_130_325 *) let cp = [| angMom_c.(0) ; angMom_c.(1) ; angMom_c.(2) |] - and dm = [| angMom_d.(0) ; angMom_d.(1) ; angMom_d.(2) |] - and xyz = + (* 4_130_325 *) + let (angMom_cx, angMom_cy, angMom_cz) = angMom_c + and (angMom_dx, angMom_dy, angMom_dz) = angMom_d in + let cp, dm, xyz = match angMom_d with - | [|_;0;0|] -> (* 1_524_451 *) 0 - | [|_;_;0|] -> (* 1_302_937 *) 1 - | _ -> (* 1_302_937 *) 2 + | (_,0,0) -> (* 1_524_451 *) (angMom_cx+1, angMom_cy, angMom_cz), (angMom_dx-1, angMom_dy, angMom_dz), 0 + | (_,_,0) -> (* 1_302_937 *) (angMom_cx, angMom_cy+1, angMom_cz), (angMom_dx, angMom_dy-1, angMom_dz), 1 + | _ -> (* 1_302_937 *) (angMom_cx, angMom_cy, angMom_cz+1), (angMom_dx, angMom_dy, angMom_dz-1), 2 in - cp.(xyz) <- cp.(xyz) + 1; - dm.(xyz) <- dm.(xyz) - 1; let h1 = hrr angMom_a angMom_b cp dm totAngMom_a totAngMom_b (totAngMom_c+1) (totAngMom_d-1) in @@ -230,7 +239,12 @@ let hvrr_two_e (angMom_a, angMom_b, angMom_c, angMom_d) h1 +. f2 *. h2 in - hrr angMom_a angMom_b angMom_c angMom_d totAngMom_a totAngMom_b totAngMom_c totAngMom_d + hrr + (angMom_a.(0),angMom_a.(1),angMom_a.(2)) + (angMom_b.(0),angMom_b.(1),angMom_b.(2)) + (angMom_c.(0),angMom_c.(1),angMom_c.(2)) + (angMom_d.(0),angMom_d.(1),angMom_d.(2)) + totAngMom_a totAngMom_b totAngMom_c totAngMom_d diff --git a/Basis/TwoElectronRRVectorized.ml b/Basis/TwoElectronRRVectorized.ml index f17ecae..307370e 100644 --- a/Basis/TwoElectronRRVectorized.ml +++ b/Basis/TwoElectronRRVectorized.ml @@ -53,10 +53,10 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d) let am = [| angMom_a.(0) ; angMom_a.(1) ; angMom_a.(2) |] and amm = [| angMom_a.(0) ; angMom_a.(1) ; angMom_a.(2) |] and xyz = - match angMom_a with - | [|0;0;_|] -> 2 - | [|0;_;_|] -> 1 - | _ -> 0 + match (angMom_a.(0),angMom_a.(1),angMom_a.(2)) with + | (0,0,_) -> 2 + | (0,_,_) -> 1 + | _ -> 0 in am.(xyz) <- am.(xyz) - 1; amm.(xyz) <- amm.(xyz) - 2; @@ -112,10 +112,10 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d) and cm = [| angMom_c.(0) ; angMom_c.(1) ; angMom_c.(2) |] and cmm = [| angMom_c.(0) ; angMom_c.(1) ; angMom_c.(2) |] and xyz = - match angMom_c with - | [|0;0;_|] -> 2 - | [|0;_;_|] -> 1 - | _ -> 0 + match (angMom_c.(0),angMom_c.(1),angMom_c.(2)) with + | (0,0,_) -> 2 + | (0,_,_) -> 1 + | _ -> 0 in am.(xyz) <- am.(xyz) - 1; cm.(xyz) <- cm.(xyz) - 1; @@ -217,10 +217,10 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d) let ap = [| angMom_a.(0) ; angMom_a.(1) ; angMom_a.(2) |] and bm = [| angMom_b.(0) ; angMom_b.(1) ; angMom_b.(2) |] and xyz = - match angMom_b with - | [|0;0;_|] -> 2 - | [|0;_;_|] -> 1 - | _ -> 0 + match (angMom_b.(0),angMom_b.(1),angMom_b.(2)) with + | (0,0,_) -> 2 + | (0,_,_) -> 1 + | _ -> 0 in ap.(xyz) <- ap.(xyz) + 1; bm.(xyz) <- bm.(xyz) - 1; @@ -244,10 +244,10 @@ let hvrr_two_e_vector (angMom_a, angMom_b, angMom_c, angMom_d) let cp = [| angMom_c.(0) ; angMom_c.(1) ; angMom_c.(2) |] and dm = [| angMom_d.(0) ; angMom_d.(1) ; angMom_d.(2) |] and xyz = - match angMom_d with - | [|0;0;_|] -> 2 - | [|0;_;_|] -> 1 - | _ -> 0 + match (angMom_d.(0),angMom_d.(1),angMom_d.(2)) with + | (0,0,_) -> 2 + | (0,_,_) -> 1 + | _ -> 0 in cp.(xyz) <- cp.(xyz) + 1; dm.(xyz) <- dm.(xyz) - 1;