Working on f12 notebook

This commit is contained in:
Anthony Scemama 2019-12-12 18:41:14 +01:00
parent e66d790cb7
commit efe0be4bb2
1 changed files with 145 additions and 165 deletions

View File

@ -694,7 +694,7 @@
" List.map (fun particle ->\n",
" List.map (fun hole' -> \n",
" List.map (fun particle' ->\n",
" if hole' > hole && particle' < particle then\n",
" if hole' > hole && particle' < particle && hole <> particle && hole' <> particle' then\n",
" Some (Determinant.double_excitation Spin.Beta hole particle Spin.Beta hole' particle' ki)\n",
" else None\n",
" ) mos_abs\n",
@ -706,6 +706,7 @@
" List.map (fun particle ->\n",
" List.map (fun hole' -> \n",
" List.map (fun particle' ->\n",
" if hole = particle || hole' = particle' then None else\n",
" Some (Determinant.double_excitation Spin.Alfa hole particle Spin.Beta hole' particle' ki)\n",
" ) mos_abs\n",
" ) mos_b \n",
@ -836,7 +837,8 @@
" else\n",
" Printf.printf \"Failed: (%d, %d) | %e %e | %e\\n%!\" i j d v error\n",
"\n",
" "
" \n",
"let sum l f = List.fold_left (fun accu i -> accu +. f i) 0. l"
]
},
{
@ -898,27 +900,26 @@
" (* Alpha-Beta *)\n",
" let s, s' = Spin.(Alfa, Beta) in \n",
" let mos_a, mos_b = mos_a ki, mos_b ki in\n",
" List.fold_left (fun accu a -> accu +.\n",
" List.fold_left (fun accu i -> accu +. \n",
" \n",
" sum mos_cabs (fun a ->\n",
" sum mos_a (fun i -> \n",
" (h_one i a s +.\n",
" List.fold_left (fun accu j -> accu +. h_two i j a j s s ) 0. mos_a +. \n",
" List.fold_left (fun accu j -> accu +. h_two i j a j s s') 0. mos_b ) *. \n",
" sum mos_a (fun j -> h_two i j a j s s ) +.\n",
" sum mos_b (fun j -> h_two i j a j s s') ) *. \n",
" (f_one a i s +.\n",
" List.fold_left (fun accu j -> accu +. f_two a j i j s s ) 0. mos_a +.\n",
" List.fold_left (fun accu j -> accu +. f_two a j i j s s') 0. mos_b )\n",
" ) 0. mos_a\n",
" ) 0. mos_cabs\n",
" +.\n",
" List.fold_left (fun accu a -> accu +.\n",
" List.fold_left (fun accu i -> accu +. \n",
" sum mos_a (fun j -> f_two a j i j s s ) +.\n",
" sum mos_b (fun j -> f_two a j i j s s') )\n",
" )\n",
" +.\n",
" sum mos_b (fun i -> \n",
" (h_one i a s +.\n",
" List.fold_left (fun accu j -> accu +. h_two i j a j s' s') 0. mos_b +. \n",
" List.fold_left (fun accu j -> accu +. h_two i j a j s' s ) 0. mos_a ) *. \n",
" sum mos_b (fun j -> h_two i j a j s' s') +.\n",
" sum mos_a (fun j -> h_two i j a j s' s ) ) *. \n",
" (f_one a i s +.\n",
" List.fold_left (fun accu j -> accu +. f_two a j i j s' s ) 0. mos_a +.\n",
" List.fold_left (fun accu j -> accu +. f_two a j i j s' s') 0. mos_b )\n",
" ) 0. mos_b\n",
" ) 0. mos_cabs\n",
" sum mos_a (fun j -> f_two a j i j s' s ) +.\n",
" sum mos_b (fun j -> f_two a j i j s' s') )\n",
" ) \n",
" ) \n",
" \n",
"let _ =\n",
" check 0 integral_value 0 1 1 1\n",
@ -992,40 +993,29 @@
" \n",
" (* Alpha-Beta *)\n",
" let s, s' = Spin.(Alfa, Beta) in\n",
" List.fold_left (fun accu j -> accu +. \n",
" List.fold_left (fun accu i -> accu +. \n",
" List.fold_left (fun accu k -> accu +. \n",
" List.fold_left (fun accu a ->\n",
" accu +. h_two i j a k s s' *. f_two a k i j s s'\n",
" ) 0. mos_cabs\n",
" ) 0. mos_virt_b\n",
" ) 0. (mos_a ki)\n",
" ) 0. (mos_b ki)\n",
" sum mos_cabs (fun a ->\n",
" sum (mos_b ki) (fun j -> \n",
" sum (mos_a ki) (fun i -> \n",
" sum mos_virt_b (fun k -> \n",
" h_two i j a k s s' *. f_two a k i j s s'\n",
" )))\n",
" +.\n",
" List.fold_left (fun accu j -> accu +. \n",
" List.fold_left (fun accu i -> accu +. \n",
" List.fold_left (fun accu k -> accu +. \n",
" List.fold_left (fun accu a ->\n",
" accu +. h_two i j a k s s' *. f_two a k i j s s'\n",
" ) 0. mos_cabs\n",
" ) 0. mos_virt_a\n",
" ) 0. (mos_b ki)\n",
" ) 0. (mos_a ki)\n",
" sum (mos_a ki) (fun j -> \n",
" sum (mos_b ki) (fun i -> \n",
" sum mos_virt_a (fun k -> \n",
" h_two i j a k s s' *. f_two a k i j s s'\n",
" )))) \n",
" +.\n",
" (* Alpha-Alpha / Beta-Beta *)\n",
" List.fold_left (fun accu (mos_virt,mos,s) -> \n",
" let s' = s in accu +.\n",
" List.fold_left (fun accu j -> accu +. \n",
" List.fold_left (fun accu i -> accu +. \n",
" if i < j then accu else\n",
" List.fold_left (fun accu k -> accu +. \n",
" List.fold_left (fun accu a -> \n",
" accu +. h_two i j a k s s' *. f_two a k i j s s'\n",
" ) 0. mos_cabs\n",
" ) 0. mos_virt\n",
" ) 0. mos\n",
" ) 0. mos\n",
" ) 0. [ (mos_virt_a,mos_a ki,Spin.Alfa) ; (mos_virt_b,mos_b ki,Spin.Beta)]\n",
" sum mos (fun j -> \n",
" sum mos (fun i -> if i < j then 0. else\n",
" sum mos_virt (fun k -> \n",
" sum mos_cabs (fun a -> \n",
" h_two i j a k s s' *. f_two a k i j s s'\n",
" )))) \n",
" ) 0. [ (mos_virt_a,mos_a ki,Spin.Alfa) ; (mos_virt_b,mos_b ki,Spin.Beta)]\n",
" \n",
"let _ =\n",
" check 0 integral_value 0 1 2 2\n",
@ -1083,31 +1073,23 @@
"let integral_value ki kj = \n",
" (* Alpha-Beta *)\n",
" let s, s' = Spin.(Alfa, Beta) in\n",
" List.fold_left (fun accu j -> accu +. \n",
" List.fold_left (fun accu i -> accu +. \n",
" List.fold_left (fun accu b -> accu +. \n",
" List.fold_left (fun accu a ->\n",
" accu +. h_two i j a b s s' *. f_two a b i j s s'\n",
" ) 0. mos_cabs\n",
" ) 0. mos_cabs\n",
" ) 0. (mos_a ki)\n",
" ) 0. (mos_b ki)\n",
" sum mos_cabs (fun a -> \n",
" sum mos_cabs (fun b -> \n",
" sum (mos_a ki) (fun i -> \n",
" sum (mos_b ki) (fun j -> \n",
" h_two i j a b s s' *. f_two a b i j s s'\n",
" ))))\n",
" +.\n",
" (* Alpha-Alpha / Beta-Beta *)\n",
" List.fold_left (fun accu (mos,s) -> \n",
" let s' = s in accu +.\n",
" List.fold_left (fun accu j -> accu +. \n",
" List.fold_left (fun accu i -> accu +. \n",
" if i < j then accu else\n",
" List.fold_left (fun accu b -> accu +. \n",
" List.fold_left (fun accu a -> \n",
" if b > a then accu else\n",
" accu +. h_two i j a b s s' *. f_two a b i j s s'\n",
" ) 0. mos_cabs\n",
" ) 0. mos_cabs\n",
" ) 0. mos\n",
" ) 0. mos\n",
" ) 0. [ (mos_a ki,Spin.Alfa) ; (mos_b ki,Spin.Beta)]\n",
" sum mos_cabs (fun b -> \n",
" sum mos_cabs (fun a -> if b > a then 0. else\n",
" sum mos (fun j -> \n",
" sum mos (fun i -> if i < j then 0. else\n",
" h_two i j a b s s' *. f_two a b i j s s'\n",
" ))))\n",
" ) 0. [ (mos_a ki,Spin.Alfa) ; (mos_b ki,Spin.Beta)]\n",
" \n",
"let _ =\n",
" check 0 integral_value 0 2 2 2\n",
@ -1176,14 +1158,15 @@
"\n",
" let result = \n",
" let s' = Spin.other s in\n",
" List.fold_left (fun accu a -> accu +.\n",
" \n",
" sum mos_cabs (fun a -> \n",
" (h_one i a s +. \n",
" (List.fold_left (fun accu j -> accu +. h_two i j a j s s ) 0. mos_i ) +.\n",
" (List.fold_left (fun accu j -> accu +. h_two i j a j s s') 0. mos_i') ) *.\n",
" sum mos_i (fun j -> h_two i j a j s s ) +.\n",
" sum mos_i' (fun j -> h_two i j a j s s') ) *.\n",
" (f_one a k s +. \n",
" (List.fold_left (fun accu j -> accu +. f_two a j k j s s ) 0. mos_j ) +.\n",
" (List.fold_left (fun accu j -> accu +. f_two a j k j s s') 0. mos_j') ) \n",
" ) 0. mos_cabs\n",
" sum mos_j (fun j -> f_two a j k j s s ) +.\n",
" sum mos_j' (fun j -> f_two a j k j s s') ) \n",
" )\n",
" in\n",
" match phase with\n",
" | Phase.Pos -> result\n",
@ -1269,21 +1252,21 @@
"\n",
" let result = \n",
" let s' = Spin.other s in\n",
" List.fold_left (fun accu a -> accu +.\n",
" List.fold_left (fun accu j -> accu +.\n",
" sum mos_cabs (fun a -> \n",
" sum mos_j (fun j -> \n",
" (h_one j a s +. \n",
" (List.fold_left (fun accu m -> accu +. h_two j m a m s s ) 0. mos_i ) +.\n",
" (List.fold_left (fun accu m -> accu +. h_two j m a m s s') 0. mos_i') ) *.\n",
" sum mos_i (fun m -> h_two j m a m s s ) +.\n",
" sum mos_i' (fun m -> h_two j m a m s s') ) *.\n",
" (f_two a i j k s s )\n",
" ) 0. mos_j\n",
" ) \n",
" +.\n",
" List.fold_left (fun accu j -> accu +.\n",
" sum mos_j' (fun j -> \n",
" (h_one j a s +. \n",
" (List.fold_left (fun accu m -> accu +. h_two j m a m s' s ) 0. mos_i ) +.\n",
" (List.fold_left (fun accu m -> accu +. h_two j m a m s' s') 0. mos_i') ) *.\n",
" sum mos_i (fun m -> h_two j m a m s' s ) +.\n",
" sum mos_i' (fun m -> h_two j m a m s' s') ) *.\n",
" (f_two a i j k s' s )\n",
" ) 0. mos_j'\n",
" ) 0. mos_cabs\n",
" ) \n",
" ) \n",
" in\n",
" match phase with\n",
" | Phase.Pos -> result\n",
@ -1371,21 +1354,20 @@
"\n",
" let result = \n",
" let s' = Spin.other s in\n",
" List.fold_left (fun accu a -> accu +.\n",
" List.fold_left (fun accu j -> accu +.\n",
" sum mos_cabs (fun a -> \n",
" sum mos_j (fun j -> \n",
" (h_two i j k a s s ) *.\n",
" (f_one a j s +. \n",
" (List.fold_left (fun accu m -> accu +. f_two m a m j s s) 0. mos_i ) +.\n",
" (List.fold_left (fun accu m -> accu +. f_two m a m j s' s) 0. mos_i') ) \n",
" ) 0. mos_j\n",
" +.\n",
" List.fold_left (fun accu j -> accu +.\n",
" sum mos_i (fun m -> f_two m a m j s s) +.\n",
" sum mos_i' (fun m -> f_two m a m j s' s) ) \n",
" ) \n",
" +.\n",
" sum mos_j' (fun j -> \n",
" (h_two i j k a s s') *.\n",
" (f_one a j s' +. \n",
" (List.fold_left (fun accu m -> accu +. f_two m a m j s s') 0. mos_i ) +.\n",
" (List.fold_left (fun accu m -> accu +. f_two m a m j s' s') 0. mos_i') ) \n",
" ) 0. mos_j'\n",
" ) 0. mos_cabs\n",
" sum mos_i (fun m -> f_two m a m j s s') +.\n",
" sum mos_i' (fun m -> f_two m a m j s' s') ) \n",
" )) \n",
" in\n",
" match phase with\n",
" | Phase.Pos -> result\n",
@ -1417,26 +1399,32 @@
"source": [
"$|J\\rangle = \\hat{T}_{i}^{k}|I\\rangle$\n",
"\n",
"$\\sum_{j} \\sum_{a} \\sum_{m} \\hat{T}_{ij}^{ma} \\hat{T}_{ma}^{kj} + \n",
" \\hat{T}_{lj}^{ka} \\hat{T}_{ia}^{lj} +\n",
" \\hat{T}_{ij}^{am} \\hat{T}_{am}^{kj} \n",
" $\n",
"$\\sum_{j} \\sum_{a} \\sum_{m}\n",
" \\hat{T}_{ij}^{ma} \\hat{T}_{ma}^{kj} + \n",
" \\hat{T}_{lj}^{ka} \\hat{T}_{ia}^{lj} + \\\\\n",
"\\sum_{\\bar{j}} \\sum_{\\bar{a}} \\sum_{m}\n",
" \\hat{T}_{i\\bar{j}}^{m\\bar{a}} \\hat{T}_{m\\bar{a}}^{k\\bar{j}} + \n",
" \\hat{T}_{i\\bar{j}}^{a\\bar{m}} \\hat{T}_{a\\bar{m}}^{k\\bar{j}} +\n",
" \\hat{T}_{l\\bar{j}}^{k\\bar{a}} \\hat{T}_{i\\bar{a}}^{l\\bar{j}}\n",
"$\n",
"\n",
"$j,l$ : orbital indices of MOs occupied both in $|I\\rangle$ and $|J\\rangle$.\n",
"\n",
"$m$ : orbital indices of virtual MOs, unoccupied both in $|I\\rangle$ and $|J\\rangle$.\n",
"\n",
"$a$ : orbital indices of CABS MOs\n",
"\n",
"\\begin{align}\n",
" \\sum_{a} \\sum_{j} \\sum_{m}\n",
" \\langle i j || m a \\rangle \\left[ m a || k j \\right] + \n",
" \\sum_{\\bar{a}} \\sum_{\\bar{j}} \\sum_{m} \n",
" \\langle i \\bar{j} | m \\bar{a} \\rangle \\left[ m \\bar{a} | k \\bar{j} \\right] + \\\\\n",
" \\sum_{a} \\sum_{j} \\sum_{l}\n",
" \\langle l j || k a \\rangle \\left[ i a || l j \\right] + \n",
" \\sum_{\\bar{a}} \\sum_{\\bar{j}} \\sum_{l} \n",
" \\langle l \\bar{j} | k \\bar{a} \\rangle \\left[ i \\bar{a} | l \\bar{j} \\right] + \\\\\n",
" \\sum_{a} \\sum_{\\bar{j}} \\sum_{\\bar{m}} \n",
" \\langle i \\bar{j} | a \\bar{m} \\rangle \\left[ a \\bar{m} | k \\bar{j} \\right] \\\\\n",
"% \\sum_{a} \\sum_{j} \\sum_{m}\n",
"% \\langle i j || m a \\rangle \\left[ m a || k j \\right] + \n",
"% \\sum_{\\bar{a}} \\sum_{\\bar{j}} \\sum_{m} \n",
"% \\langle i \\bar{j} | m \\bar{a} \\rangle \\left[ m \\bar{a} | k \\bar{j} \\right] + \\\\\n",
"% \\sum_{a} \\sum_{j} \\sum_{l}\n",
"% \\langle l j || k a \\rangle \\left[ i a || l j \\right] + \n",
"% \\sum_{\\bar{a}} \\sum_{\\bar{j}} \\sum_{l} \n",
"% \\langle l \\bar{j} | k \\bar{a} \\rangle \\left[ i \\bar{a} | l \\bar{j} \\right] + \\\\\n",
"% \\sum_{a} \\sum_{\\bar{j}} \\sum_{\\bar{m}} \n",
"% \\langle i \\bar{j} | a \\bar{m} \\rangle \\left[ a \\bar{m} | k \\bar{j} \\right] \\\\\n",
"\\end{align}"
]
},
@ -1456,7 +1444,7 @@
" | _ -> assert false\n",
" in\n",
" \n",
" (* mos unoccupied in both I and J *)\n",
" (* MOs unoccupied in both I and J *)\n",
" let mos_virt_a, mos_virt_b = \n",
" Array.init mo_num (fun i -> Some (i+1)) , \n",
" Array.init mo_num (fun i -> Some (i+1)) \n",
@ -1471,10 +1459,10 @@
" Array.to_list mos_virt_b |> Util.list_some\n",
" in\n",
"\n",
" let mos_virt =\n",
" let mos_virt, mos_virt' =\n",
" match s with\n",
" | Spin.Alfa -> mos_virt_a\n",
" | Spin.Beta -> mos_virt_b\n",
" | Spin.Alfa -> mos_virt_a, mos_virt_b\n",
" | Spin.Beta -> mos_virt_b, mos_virt_a\n",
" in\n",
" \n",
" let mos_j, mos_j' =\n",
@ -1495,35 +1483,38 @@
"\n",
" let result = \n",
" let s' = Spin.other s in\n",
" List.fold_left (fun accu a -> accu +.\n",
" List.fold_left (fun accu m -> accu +.\n",
" List.fold_left (fun accu j -> accu +.\n",
" h_two i j m a s s *. f_two m a k j s s ) 0. mos_j\n",
" +.\n",
" List.fold_left (fun accu j -> accu +.\n",
" h_two i j m a s s' *. f_two m a k j s s' ) 0. mos_j'\n",
" ) 0. mos_virt\n",
" -.\n",
" List.fold_left (fun accu l -> accu +.\n",
" List.fold_left (fun accu j -> accu +.\n",
" h_two l j k a s s *. f_two i a l j s s ) 0. mos_j\n",
" +.\n",
" List.fold_left (fun accu j -> accu +.\n",
" h_two l j k a s s' *. f_two i a l j s s' ) 0. mos_j'\n",
" ) 0. mos_j\n",
" +.\n",
" List.fold_left (fun accu m -> accu +.\n",
" List.fold_left (fun accu j -> accu +.\n",
" h_two i j a m s s' *. f_two a m k j s s' ) 0. mos_j'\n",
" ) 0. mos_virt\n",
" ) 0. mos_cabs\n",
" sum mos_cabs (fun a -> \n",
" sum mos_j (fun j -> \n",
" sum mos_virt (fun m -> \n",
" h_two i j m a s s *. f_two m a k j s s \n",
" ))) +.\n",
" sum mos_cabs (fun a -> \n",
" sum mos_j' (fun j -> \n",
" sum mos_virt (fun m -> \n",
" h_two i j m a s s' *. f_two m a k j s s' \n",
" ))) +.\n",
" sum mos_cabs (fun a -> \n",
" sum mos_j' (fun j -> \n",
" sum mos_virt' (fun m -> \n",
" h_two i j a m s s' *. f_two a m k j s s'\n",
" ))) +.\n",
" sum mos_cabs (fun a -> \n",
" sum mos_j (fun j -> \n",
" sum mos_j (fun l -> if l = i then 0. else\n",
" h_two l j k a s s *. f_two i a l j s s \n",
" ))) +. \n",
" sum mos_cabs (fun a -> \n",
" sum mos_j' (fun j -> \n",
" sum mos_j (fun l -> if l = i then 0. else\n",
" h_two l j k a s s' *. f_two i a l j s s'\n",
" ))) \n",
" in\n",
" match phase with\n",
" | Phase.Pos -> result\n",
" | Phase.Neg -> -. result\n",
" \n",
"\n",
"let _ = check 20 integral_value 1 1 2 2\n",
"let _ = check 11 integral_value 1 1 2 2\n",
"\n",
"\n"
]
@ -1699,25 +1690,17 @@
"\n",
" let result = \n",
" (* Alpha-Beta *)\n",
" List.fold_left (fun accu j ->\n",
" let s' = Spin.other s in\n",
" accu +. List.fold_left (fun accu b ->\n",
" accu +. List.fold_left (fun accu a ->\n",
" accu +. h_two h j a b s s' *. f_two a b p j s s'\n",
" ) 0. mos_cabs\n",
" ) 0. mos_cabs\n",
" ) 0. mos' \n",
" let s' = Spin.other s in\n",
" sum mos_cabs (fun b -> \n",
" sum mos_cabs (fun a -> \n",
" sum mos' (fun j -> h_two h j a b s s' *. f_two a b p j s s')\n",
" ))\n",
" +.\n",
" (* Alpha-Alpha / Beta-Beta *)\n",
" List.fold_left (fun accu j ->\n",
" let s' = s in \n",
" accu +. List.fold_left (fun accu b ->\n",
" accu +. List.fold_left (fun accu a -> \n",
" if b > a then accu else\n",
" accu +. h_two h j a b s s' *. f_two a b p j s s'\n",
" ) 0. mos_cabs\n",
" ) 0. mos_cabs\n",
" ) 0. mos \n",
" sum mos_cabs (fun b -> \n",
" sum mos_cabs (fun a -> if b >= a then 0. else\n",
" sum mos (fun j -> h_two h j a b s s *. f_two a b p j s s)\n",
" )) \n",
" in\n",
" match phase with\n",
" | Phase.Pos -> result\n",
@ -1821,18 +1804,15 @@
"\n",
" let result = \n",
" if s <> s' then (* Alpha-Beta *)\n",
" List.fold_left (fun accu b ->\n",
" accu +. List.fold_left (fun accu a ->\n",
" accu +. h_two h h' a b s s' *. f_two a b p p' s s'\n",
" ) 0. mos_cabs\n",
" ) 0. mos_cabs\n",
" sum mos_cabs (fun b -> \n",
" sum mos_cabs (fun a -> \n",
" h_two h h' a b s s' *. f_two a b p p' s s'\n",
" )) \n",
" else (* Alpha-Alpha / Beta-Beta *)\n",
" List.fold_left (fun accu b ->\n",
" accu +. List.fold_left (fun accu a -> \n",
" if b > a then accu else\n",
" accu +. h_two h h' a b s s' *. f_two a b p p' s s'\n",
" ) 0. mos_cabs\n",
" ) 0. mos_cabs\n",
" sum mos_cabs (fun b -> \n",
" sum mos_cabs (fun a -> if b >= a then 0. else\n",
" h_two h h' a b s s' *. f_two a b p p' s s'\n",
" )) \n",
" in\n",
" match phase with\n",
" | Phase.Pos -> result\n",