Working on F12 notebook

This commit is contained in:
Anthony Scemama 2019-12-17 19:03:55 +01:00
parent efe0be4bb2
commit 3b22e014b4
1 changed files with 190 additions and 129 deletions

View File

@ -123,9 +123,9 @@
"metadata": {},
"outputs": [],
"source": [
"let basis_filename = \"/home/scemama/qp2/data/basis/6-31g\" \n",
"let aux_basis_filename = \"/home/scemama/qp2/data/basis/cc-pvdz\" \n",
"let nuclei = Nuclei.of_zmt_string \"be\" \n",
"let basis_filename = \"/home/scemama/qp2/data/basis/sto-3g\" \n",
"let aux_basis_filename = \"/home/scemama/qp2/data/basis/6-31g\" \n",
"let nuclei = Nuclei.of_zmt_string \"c\" \n",
"let frozen_core = false\n",
"let multiplicity = 1\n",
"let state = 1\n",
@ -1381,7 +1381,7 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"# 7. *1 1 22*"
"# 7. 1 1 22"
]
},
{
@ -1400,11 +1400,11 @@
"$|J\\rangle = \\hat{T}_{i}^{k}|I\\rangle$\n",
"\n",
"$\\sum_{j} \\sum_{a} \\sum_{m}\n",
" \\hat{T}_{ij}^{ma} \\hat{T}_{ma}^{kj} + \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}_{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",
@ -1418,9 +1418,9 @@
"% \\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",
"% \\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",
"% \\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",
@ -1497,12 +1497,12 @@
" 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",
" ))) -.\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",
" ))) -. \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",
@ -1514,126 +1514,11 @@
" | Phase.Neg -> -. result\n",
" \n",
"\n",
"let _ = check 11 integral_value 1 1 2 2\n",
"let _ = check 0 integral_value 1 1 2 2\n",
"\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"scrolled": true
},
"outputs": [],
"source": [
"let ki = det_I.(0)\n",
"let kj = det_I.(11)\n",
"\n",
"let alpha = \n",
" generate_alphas ki kj 1 1 2 2\n",
" |> List.filter (fun alpha -> abs_float (compute_HaaF ki [alpha] kj) < 2.0e-05)\n",
" |> List.fold_left (fun (a,mx) alpha ->\n",
" let x = abs_float (compute_HaaF ki [alpha] kj) in\n",
" if x > mx then (alpha, x) else (a,mx) ) (ki,-1.)\n",
" |> fst\n",
" \n",
"\n",
"let _ = Determinant.to_lists ki\n",
"let _ = Determinant.to_lists alpha\n",
"let _ = Determinant.to_lists kj\n",
"\n",
"let _ = \n",
" Format.printf \"Excitation: %a\\n@.\" Excitation.pp_exc (Excitation.of_det ki alpha);\n",
" Format.printf \"Excitation: %a\\n@.\" Excitation.pp_exc (Excitation.of_det alpha kj);\n",
" Format.printf \"Excitation: %a\\n@.\" Excitation.pp_exc (Excitation.of_det ki kj)\n",
"\n",
"let integral_value ki kj = \n",
" let i, k, s, phase =\n",
" match Excitation.of_det ki kj with\n",
" | Excitation.(Single (phase, { hole ; particle ; spin })) ->\n",
" hole, particle, spin, phase\n",
" | _ -> assert false\n",
" in\n",
" \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",
" in\n",
" List.iter (fun i -> mos_virt_a.(i-1) <- None) (mos_a ki);\n",
" List.iter (fun i -> mos_virt_a.(i-1) <- None) (mos_a kj);\n",
" List.iter (fun i -> mos_virt_b.(i-1) <- None) (mos_b ki);\n",
" List.iter (fun i -> mos_virt_b.(i-1) <- None) (mos_b kj);\n",
" \n",
" let mos_virt_a, mos_virt_b = \n",
" Array.to_list mos_virt_a |> Util.list_some,\n",
" Array.to_list mos_virt_b |> Util.list_some\n",
" in\n",
"\n",
" let mos_virt =\n",
" match s with\n",
" | Spin.Alfa -> mos_virt_a\n",
" | Spin.Beta -> mos_virt_b\n",
" in\n",
" \n",
" let mos_j, mos_j' =\n",
" let alfa = \n",
" let i = Spindeterminant.bitstring @@ Determinant.alfa ki in\n",
" let j = Spindeterminant.bitstring @@ Determinant.alfa kj in\n",
" Bitstring.to_list (Bitstring.logand i j)\n",
" in\n",
" let beta = \n",
" let i = Spindeterminant.bitstring @@ Determinant.beta ki in\n",
" let j = Spindeterminant.bitstring @@ Determinant.beta kj in\n",
" Bitstring.to_list (Bitstring.logand i j)\n",
" in\n",
" match s with\n",
" | Spin.Alfa -> alfa, beta\n",
" | Spin.Beta -> beta, alfa\n",
" in\n",
"\n",
" let result = \n",
" let s' = Spin.other s in\n",
" let a = 13 in\n",
" (*\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",
" *)\n",
" let a = 16\n",
" and j = 1\n",
" and m = 7\n",
" in\n",
" (*\n",
" h_two i j a m s s' *. f_two a m k j s s'\n",
" *)\n",
" h_two i j m a s s' *. f_two m a k j s s'\n",
" (*\n",
" let a = 18\n",
" and j = 1\n",
" and l = 2\n",
" in\n",
" -. h_two l j k a s s' *. f_two i a l j s s' \n",
" *)\n",
"\n",
" in\n",
" \n",
" match phase with\n",
" | Phase.Pos -> result\n",
" | Phase.Neg -> -. result\n",
" \n",
"\n",
"let _ = compute_HaaF ki [alpha] kj\n",
"let _ = integral_value ki kj\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
@ -1658,11 +1543,18 @@
"\n",
"$j$ : orbital indices of MOs occupied both in $|I\\rangle$ and $|J\\rangle$.\n",
"\n",
"$\\sum_{j} \\sum_{a} \\sum_{b}\n",
" \\hat{T}_{ij}^{ab} \\hat{T}_{ab}^{kj} +\n",
"\\sum_{\\bar{j}} \\sum_{a} \\sum_{\\bar{b}}\n",
" \\hat{T}_{i\\bar{j}}^{a\\bar{b}} \\hat{T}_{a\\bar{b}}^{k\\bar{j}} \n",
"$\n",
"\n",
"\n",
"$$\n",
"\\sum_{j} \\sum_{b}\\sum_{a<b}\n",
" \\langle i j || a b \\rangle \n",
" \\left[ ab || kj \\right]\n",
" \\langle i j || a b \\rangle \\left[ ab || kj \\right] + \n",
"\\sum_{\\bar{j}} \\sum_{\\bar{b}}\\sum_{a}\n",
" \\langle i \\bar{j} | a \\bar{b} \\rangle \\left[ a \\bar{b} | k \\bar{j} \\right] \n",
"$$"
]
},
@ -1714,7 +1606,7 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"# 9. *2 1 12*"
"# 9. 2 1 12"
]
},
{
@ -1726,6 +1618,175 @@
"png_image \"2_1_12.png\""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"$|J\\rangle = \\hat{T}_{ij}^{kl}|I\\rangle$\n",
"\n",
"$m$ : orbital indices of MOs unoccupied both in $|I\\rangle$ and $|J\\rangle$.\n",
"\n",
"$n$ : orbital indices of MOs occupied in $|I\\rangle$\n",
"\n",
"$\\sum_{a} \n",
" \\hat{T}_{i}^{a} \\hat{T}_{aj}^{kl} +\n",
" \\hat{T}_{j}^{a} \\hat{T}_{ai}^{lk} +\n",
" \\hat{T}_{i}^{a} \\hat{T}_{a\\bar{j}}^{k\\bar{l}} + \n",
" \\hat{T}_{j}^{a} \\hat{T}_{a\\bar{i}}^{l\\bar{k}} \n",
"$\n",
"\n",
"\n",
"$$\n",
"\\sum_{a}\n",
" \\left( h_{ia} + \\sum_{n} \\langle i n || a n \\rangle + \\sum_{\\bar{n}} \\langle i \\bar{n} | a \\bar{n} \\rangle \\right) \\left[ a j || k l \\right] + \\\\\n",
"\\sum_{a}\n",
" \\left( h_{ja} + \\sum_{n} \\langle j n || a n \\rangle + \\sum_{\\bar{n}} \\langle j \\bar{n} | a \\bar{n} \\rangle \\right) \\left[ a i || l k \\right] \n",
"$$"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"let integral_value ki kj = \n",
" let i, j, k, l, s, s', phase =\n",
" match Excitation.of_det ki kj with\n",
" | Excitation.(Double (phase, { hole ; particle ; spin }, {hole=hole' ; particle=particle' ; spin=spin' })) ->\n",
" hole, hole', particle, particle', spin, spin', phase\n",
" | _ -> assert false\n",
" in\n",
"\n",
" let mos, mos' =\n",
" match s with\n",
" | Spin.Alfa -> mos_a ki, mos_b ki\n",
" | Spin.Beta -> mos_b ki, mos_a ki\n",
" in\n",
"\n",
" let result = \n",
" let s'' = Spin.other s in\n",
" sum mos_cabs (fun a -> \n",
" ( h_one i a s +.\n",
" sum mos (fun n -> h_two n i n a s s) +.\n",
" sum mos' (fun n -> h_two n i n a s'' s) \n",
" ) *. f_two j a l k s' s\n",
" ) +.\n",
" sum mos_cabs (fun a -> \n",
" ( h_one j a s +.\n",
" sum mos (fun n -> h_two n j n a s s) +.\n",
" sum mos' (fun n -> h_two n j n a s'' s) \n",
" ) *. f_two i a k l s' s\n",
" ) \n",
" in\n",
"\n",
" match phase with\n",
" | Phase.Pos -> result\n",
" | Phase.Neg -> -. result\n",
"\n",
"\n",
"let _ = check 0 integral_value 2 1 1 2"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"scrolled": true
},
"outputs": [],
"source": [
"let ki = det_I.(0)\n",
"let kj = det_I.(11)\n",
"\n",
"let _ = \n",
"generate_alphas ki kj 2 1 1 2\n",
"|> Array.of_list\n",
"|> Array.mapi (fun kk alpha -> \n",
"\n",
"(*\n",
"let _ = Determinant.to_lists ki\n",
"let _ = Determinant.to_lists alpha\n",
"let _ = Determinant.to_lists kj\n",
"\n",
" \n",
"let _ = \n",
" Format.printf \"|I> -> |J> : %a\\n@.\" Excitation.pp (Excitation.of_det ki kj);\n",
" Format.printf \"|I> -> |a> : %a\\n@.\" Excitation.pp (Excitation.of_det ki alpha);\n",
" Format.printf \"|a> -> |J> : %a\\n@.\" Excitation.pp (Excitation.of_det alpha kj)\n",
"*)\n",
" let i, j, k, l, s, s', phase =\n",
" match Excitation.of_det ki kj with\n",
" | Excitation.(Double (phase, { hole ; particle ; spin }, {hole=hole' ; particle=particle' ; spin=spin' })) ->\n",
" hole, hole', particle, particle', spin, spin', phase\n",
" | _ -> assert false\n",
" in\n",
"\n",
"\n",
" let integral_value ki kj = \n",
" (*\n",
" let i, j, k, l, s, s', phase =\n",
" match Excitation.of_det ki kj with\n",
" | Excitation.(Double (phase, { hole ; particle ; spin }, {hole=hole' ; particle=particle' ; spin=spin' })) ->\n",
" hole, hole', particle, particle', spin, spin', phase\n",
" | _ -> assert false\n",
" in\n",
" *)\n",
" \n",
" let i, a, s, _ = \n",
" match Excitation.of_det ki alpha with\n",
" | Excitation.(Single (phase, { hole ; particle ; spin })) ->\n",
" hole, particle, spin, phase\n",
" | _ -> assert false\n",
" in\n",
"\n",
" let i', j', k', l', t, t', _ =\n",
" match Excitation.of_det alpha kj with\n",
" | Excitation.(Double (phase, { hole ; particle ; spin }, {hole=hole' ; particle=particle' ; spin=spin' })) ->\n",
" hole, hole', particle, particle', spin, spin', phase\n",
" | _ -> assert false\n",
" in\n",
" let j, k, l, s' = \n",
" if i' = a then\n",
" (assert (k' = k) ; j', k', l', t')\n",
" else\n",
" (assert (j'=a) ; assert (l'=k) ; i', l', k', t)\n",
" in\n",
" \n",
" Format.printf \"|I> -> |J> : %a | %d %d %d %d\\n@.\" Excitation.pp (Excitation.of_det ki kj) i j k l;\n",
" Format.printf \"|I> -> |a> : %a | %d %d\\n@.\" Excitation.pp (Excitation.of_det ki alpha) i a ;\n",
" Format.printf \"|a> -> |J> : %a | %d %d %d %d\\n@.\" Excitation.pp (Excitation.of_det alpha kj) j a l k;\n",
" \n",
" \n",
" \n",
" let mos, mos' =\n",
" match s with\n",
" | Spin.Alfa -> mos_a ki, mos_b ki\n",
" | Spin.Beta -> mos_b ki, mos_a ki\n",
" in\n",
"\n",
" let result = \n",
" let s'' = Spin.other s in\n",
" ( h_one i a s +.\n",
" sum mos (fun n -> h_two n i n a s s) +.\n",
" sum mos' (fun n -> h_two n i n a s'' s) \n",
" ) *. f_two j a l k s' s\n",
" in\n",
"\n",
" match phase with\n",
" | Phase.Pos -> result\n",
" | Phase.Neg -> -. result\n",
"\n",
"in\n",
"let a = (compute_HaaF ki [alpha] kj)\n",
"and b = (integral_value ki kj)\n",
"in\n",
"if a <> b then\n",
"Format.printf \"%6d %e %e@.\" (kk) a b\n",
")\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},