diff --git a/Notebooks/F12_fast.ipynb b/Notebooks/F12_fast.ipynb index b46b371..ef18f66 100644 --- a/Notebooks/F12_fast.ipynb +++ b/Notebooks/F12_fast.ipynb @@ -1308,7 +1308,6 @@ " | _ -> assert false\n", " in\n", " \n", - " (* MOs unoccupied in both I and J *)\n", " let mos_novirt, mos_novirt' =\n", " let alfa = \n", " let i = Spindeterminant.bitstring @@ Determinant.alfa ki in\n", @@ -1410,103 +1409,424 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 129, "metadata": {}, - "outputs": [], + "outputs": [ + { + "data": { + "text/plain": [ + "val m_2112_1H_2Fa :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2112_1H_2Fb :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2112_2Ha_2Fa :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2112_2Hb_2Fa :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2112_2Ha_2Fb :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2121_2Ha_2Fa :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2121_2Hb_2Fa :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2121_2Ha_2Fb :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2122_2Ha_2Fa_ij :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2122_2Hb_2Fb_ij :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2122_2Hb_2Fb_ij2 :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2122_2Ha_2Fa_ij2 :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2122_2Ha_2Fa_nv :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2122_2Hb_2Fb_nv :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_2122_2Hb_2Fb_nv2 :\n", + " (float, Bigarray.float64_elt, Bigarray.fortran_layout) Bigarray.Genarray.t =\n", + " \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val integral_value : Determinant.t -> Determinant.t -> float = \n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Checking ... \n", + " - 10 %\n", + " - 20 %\n", + " - 30 %\n", + " - 40 %\n", + " - 50 %\n", + " - 60 %\n", + " - 70 %\n", + " - 80 %\n", + " - 90 %\n", + " - 100 %\n", + "OK: (2, 45) | 2.332967e-02 2.332967e-02 | 2.428613e-17 | cpu : 109.679883 0.039116\n" + ] + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val ki : Determinant.t =\n", + " phase:+1\n", + " a:+1 +--++-----------------------------------------------------------\n", + " b:+1 +++-------------------------------------------------------------\n", + "\n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val kj : Determinant.t =\n", + " phase:+1\n", + " a:+1 +-+-+-----------------------------------------------------------\n", + " b:+1 ++-+------------------------------------------------------------\n", + "\n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : float = 0.000896378780435275082\n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : float = 0.000896378780435274107\n" + ] + }, + "execution_count": 129, + "metadata": {}, + "output_type": "execute_result" + } + ], "source": [ "let m_2112_1H_2Fa =\n", " array_4_init mo_num mo_num mo_num mo_num (fun i j k l ->\n", " sum mos_cabs (fun a ->\n", " h_one i a Spin.Alfa *. f_two a j k l Spin.Alfa Spin.Alfa +.\n", - " h_one j a Spin.Alfa *. f_two i a k l Spin.Alfa Spin.Alfa )\n", + " h_one j a Spin.Alfa *. f_two i a k l Spin.Alfa Spin.Alfa +. \n", + " h_two i j a l Spin.Alfa Spin.Alfa *. f_one a k Spin.Alfa +.\n", + " h_two i j k a Spin.Alfa Spin.Alfa *. f_one a l Spin.Alfa +.\n", + " sum mos_in (fun m -> -. h_two i j a m Spin.Alfa Spin.Alfa *.\n", + " f_two m a k l Spin.Alfa Spin.Alfa) +.\n", + " sum mos_cabs (fun b -> if b >= a then 0. else\n", + " h_two i j a b Spin.Alfa Spin.Alfa *. f_two a b k l Spin.Alfa Spin.Alfa \n", + " )\n", + " )\n", " )\n", "\n", "let m_2112_1H_2Fb =\n", " array_4_init mo_num mo_num mo_num mo_num (fun i j k l ->\n", " sum mos_cabs (fun a ->\n", " h_one i a Spin.Alfa *. f_two a j k l Spin.Alfa Spin.Beta +.\n", - " h_one j a Spin.Alfa *. f_two i a k l Spin.Alfa Spin.Beta)\n", + " h_one j a Spin.Alfa *. f_two i a k l Spin.Alfa Spin.Beta +.\n", + " h_two i j a l Spin.Alfa Spin.Beta *. f_one a k Spin.Alfa +.\n", + " h_two i j k a Spin.Alfa Spin.Beta *. f_one a l Spin.Alfa +.\n", + " sum mos_in (fun m ->\n", + " h_two i j a m Spin.Alfa Spin.Beta *. f_two a m k l Spin.Alfa Spin.Beta +.\n", + " h_two i j m a Spin.Alfa Spin.Beta *. f_two m a k l Spin.Alfa Spin.Beta ) +.\n", + " sum mos_cabs (fun b -> \n", + " h_two i j a b Spin.Alfa Spin.Beta *. f_two a b k l Spin.Alfa Spin.Beta \n", + " )\n", + " )\n", " )\n", "\n", "let m_2112_2Ha_2Fa =\n", - " array_5_init mo_num mo_num mo_num mo_num mo_num (fun i j k l n ->\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", " sum mos_cabs (fun a ->\n", - " h_two i n a n Spin.Alfa Spin.Alfa *.\n", - " f_two a j k l Spin.Alfa Spin.Alfa )\n", + " h_two i n a n Spin.Alfa Spin.Alfa *. f_two a j k l Spin.Alfa Spin.Alfa +. \n", + " h_two j n a n Spin.Alfa Spin.Alfa *. f_two a i l k Spin.Alfa Spin.Alfa \n", + " )\n", " )\n", " \n", "let m_2112_2Hb_2Fa =\n", - " array_5_init mo_num mo_num mo_num mo_num mo_num (fun i j k l n ->\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", " sum mos_cabs (fun a ->\n", - " h_two i n a n Spin.Alfa Spin.Beta *.\n", - " f_two a j k l Spin.Alfa Spin.Alfa )\n", + " h_two i n a n Spin.Alfa Spin.Beta *. f_two a j k l Spin.Alfa Spin.Alfa +.\n", + " h_two j n a n Spin.Alfa Spin.Beta *. f_two a i l k Spin.Alfa Spin.Alfa \n", + " )\n", " )\n", " \n", "let m_2112_2Ha_2Fb =\n", - " array_5_init mo_num mo_num mo_num mo_num mo_num (fun i j k l n ->\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", " sum mos_cabs (fun a ->\n", - " h_two i n a n Spin.Alfa Spin.Alfa *.\n", - " f_two a j k l Spin.Alfa Spin.Beta )\n", + " h_two i n a n Spin.Alfa Spin.Alfa *. f_two a j k l Spin.Alfa Spin.Beta +.\n", + " h_two j n a n Spin.Alfa Spin.Beta *. f_two a i l k Spin.Alfa Spin.Beta )\n", " )\n", " \n", - "let m_2112_2Hb_2Fb =\n", - " array_5_init mo_num mo_num mo_num mo_num mo_num (fun i j k l n ->\n", - " sum mos_cabs (fun a ->\n", - " h_two i n a n Spin.Alfa Spin.Beta *.\n", - " f_two a j k l Spin.Alfa Spin.Beta )\n", - " )\n", " \n", - "\n", - "let m_2121_2Ha_1F =\n", - " array_4_init mo_num mo_num mo_num mo_num (fun i j k l ->\n", - " sum mos_cabs (fun a ->\n", - " h_two i j a l Spin.Alfa Spin.Alfa *. f_one a k Spin.Alfa +.\n", - " h_two i j k a Spin.Alfa Spin.Alfa *. f_one a l Spin.Alfa)\n", - " )\n", - "\n", - "let m_2121_2Hb_1F =\n", - " array_4_init mo_num mo_num mo_num mo_num (fun i j k l ->\n", - " sum mos_cabs (fun a ->\n", - " h_two i j a l Spin.Alfa Spin.Beta *. f_one a k Spin.Alfa +.\n", - " h_two i j k a Spin.Alfa Spin.Beta *. f_one a l Spin.Alfa)\n", - " )\n", - "\n", "let m_2121_2Ha_2Fa =\n", - " array_5_init mo_num mo_num mo_num mo_num mo_num (fun i j k l n ->\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", " sum mos_cabs (fun a ->\n", - " h_two i j a l Spin.Alfa Spin.Alfa *.\n", - " f_two a n k n Spin.Alfa Spin.Alfa )\n", + " h_two i j a l Spin.Alfa Spin.Alfa *. f_two a n k n Spin.Alfa Spin.Alfa +.\n", + " h_two j i a k Spin.Alfa Spin.Alfa *. f_two a n l n Spin.Alfa Spin.Alfa \n", + " )\n", " )\n", " \n", + "\n", "let m_2121_2Hb_2Fa =\n", - " array_5_init mo_num mo_num mo_num mo_num mo_num (fun i j k l n ->\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", " sum mos_cabs (fun a ->\n", - " h_two i j a l Spin.Alfa Spin.Beta *.\n", - " f_two a n k n Spin.Alfa Spin.Alfa )\n", + " h_two i j a l Spin.Alfa Spin.Beta *. f_two a n k n Spin.Alfa Spin.Alfa +.\n", + " h_two j i a k Spin.Alfa Spin.Beta *. f_two a n l n Spin.Alfa Spin.Beta \n", + " )\n", " )\n", " \n", "let m_2121_2Ha_2Fb =\n", - " array_5_init mo_num mo_num mo_num mo_num mo_num (fun i j k l n ->\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", " sum mos_cabs (fun a ->\n", - " h_two i j a l Spin.Alfa Spin.Alfa *.\n", - " f_two a n k n Spin.Alfa Spin.Beta )\n", + " h_two i j a l Spin.Alfa Spin.Alfa *. f_two a n k n Spin.Alfa Spin.Beta +.\n", + " h_two j i a k Spin.Alfa Spin.Alfa *. f_two a n l n Spin.Alfa Spin.Beta \n", + " )\n", " )\n", " \n", - "let m_2121_2Hb_2Fb =\n", - " array_5_init mo_num mo_num mo_num mo_num mo_num (fun i j k l n ->\n", + "\n", + "let m_2122_2Ha_2Fa_ij =\n", + " let s = Spin.Alfa in\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", " sum mos_cabs (fun a ->\n", - " h_two i j a l Spin.Alfa Spin.Beta *.\n", - " f_two a n k n Spin.Alfa Spin.Beta )\n", + " h_two i n a k s s *. f_two j a n l s s \n", + " +. h_two i n a l s s *. f_two j a k n s s \n", + " -. h_two j n a k s s *. f_two i a n l s s \n", + " -. h_two j n a l s s *. f_two i a k n s s\n", + " )\n", " )\n", - " " - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ + " \n", + "let m_2122_2Hb_2Fb_ij =\n", + " let s, s' = Spin.(Alfa, Beta) in\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", + " sum mos_cabs (fun a ->\n", + " h_two n i a k s s *. f_two a j n l s s'\n", + " +. h_two n j a l s s' *. f_two i a k n s s \n", + " -. h_two n j k a s s' *. f_two i a n l s s'\n", + " )\n", + " )\n", + " \n", + "let m_2122_2Hb_2Fb_ij2 =\n", + " let s, s' = Spin.(Alfa, Beta) in\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", + " sum mos_cabs (fun a ->\n", + " -. h_two i n a l s s' *. f_two a j k n s s' +.\n", + " (if n < j then\n", + " h_two i n k a s s' *. f_two j a l n s' s'\n", + " +. h_two n j a l s' s' *. f_two i a k n s s'\n", + " else \n", + " -. h_two i n k a s s' *. f_two j a n l s' s' \n", + " -. h_two j n a l s' s' *. f_two i a k n s s'\n", + " )\n", + " )\n", + " )\n", + " \n", + "let m_2122_2Ha_2Fa_ij2 =\n", + " let s, s' = Spin.(Alfa, Beta) in\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", + " sum mos_cabs (fun a ->\n", + " h_two i n k a s s' *. f_two j a l n s s'\n", + " +. h_two j n l a s s' *. f_two i a k n s s'\n", + " -. h_two i n l a s s' *. f_two j a k n s s'\n", + " -. h_two j n k a s s' *. f_two i a l n s s'\n", + " )\n", + " ) \n", + " \n", + "let m_2122_2Ha_2Fa_nv =\n", + " let s = Spin.Alfa in\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", + " sum mos_cabs (fun a -> h_two i j a n s s *. f_two n a k l s s ) ) \n", + " \n", + " \n", + " \n", + "let m_2122_2Hb_2Fb_nv =\n", + " let s, s' = Spin.(Alfa, Beta) in\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", + " sum mos_cabs (fun a -> -. h_two i j a n s s' *. f_two a n k l s s' ) ) \n", + " \n", + "let m_2122_2Hb_2Fb_nv2 =\n", + " let s, s' = Spin.(Alfa, Beta) in\n", + " array_5_init mo_num mo_num mo_num mo_num mo_num (fun n i j k l ->\n", + " sum mos_cabs (fun a -> -. h_two i j n a s s' *. f_two n a k l s s' ) ) \n", + " \n", + " \n", + " \n", "let integral_value ki kj = \n", " let i, j, k, l, s, s', phase =\n", " match Excitation.of_det ki kj with\n", @@ -1541,158 +1861,59 @@ " | Spin.Alfa -> alfa, beta\n", " | Spin.Beta -> beta, alfa\n", " in\n", - "\n", - " \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", + " \n", + " let mos_novirt, mos_novirt' =\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.logor 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.logor i j)\n", + " in\n", + " match s with\n", + " | Spin.Alfa -> alfa, beta\n", + " | Spin.Beta -> beta, alfa\n", " in\n", " \n", - " let mos_virt, mos_virt' =\n", - " match s with\n", - " | Spin.Alfa -> mos_virt_a, mos_virt_b\n", - " | Spin.Beta -> mos_virt_b, mos_virt_a\n", - " in\n", - " \n", - " \n", - " \n", - " let result = \n", + " let result = \n", " if s = s' then\n", - " let s'' = Spin.other s in\n", " m_2112_1H_2Fa.{i,j,k,l} +. \n", - " sum mos_i (fun n ->\n", - " m_2112_2Ha_2Fa.{i,j,k,l,n} +. m_2112_2Ha_2Fa.{j,i,l,k,n}\n", - " ) +.\n", - " sum mos_i' (fun n ->\n", - " m_2112_2Hb_2Fa.{i,j,k,l,n} +. m_2112_2Hb_2Fa.{j,i,l,k,n}\n", - " ) +.\n", - " m_2121_2Ha_1F.{i,j,k,l} +. \n", - " sum mos_j (fun n ->\n", - " m_2121_2Ha_2Fa.{i,j,k,l,n} +. m_2121_2Ha_2Fa.{j,i,l,k,n}\n", - " ) +.\n", - " sum mos_j' (fun n ->\n", - " m_2121_2Ha_2Fb.{i,j,k,l,n} +. m_2121_2Ha_2Fb.{j,i,l,k,n}\n", - " ) +.\n", - " sum mos_cabs (fun a -> \n", - " sum mos_ij (fun n ->\n", - " h_two i n a k s s *. f_two j a n l s s \n", - " +. h_two i n a l s s *. f_two j a k n s s \n", - " -. h_two j n a k s s *. f_two i a n l s s \n", - " -. h_two j n a l s s *. f_two i a k n s s\n", - " )\n", - " +. sum mos_virt (fun m ->\n", - " -. h_two i j a m s s *. f_two m a k l s s )\n", - " +. sum mos_ij' (fun n ->\n", - " h_two i n k a s s'' *. f_two j a l n s s''\n", - " +. h_two j n l a s s'' *. f_two i a k n s s''\n", - " -. h_two i n l a s s'' *. f_two j a k n s s''\n", - " -. h_two j n k a s s'' *. f_two i a l n s s''\n", - " )\n", - " )\n", + " sum mos_i (fun n -> m_2112_2Ha_2Fa.{n,i,j,k,l} ) +.\n", + " sum mos_i' (fun n -> m_2112_2Hb_2Fa.{n,i,j,k,l} ) +.\n", + " sum mos_j (fun n -> m_2121_2Ha_2Fa.{n,i,j,k,l} ) +.\n", + " sum mos_j' (fun n -> m_2121_2Ha_2Fb.{n,i,j,k,l} ) +.\n", + " sum mos_ij (fun n -> m_2122_2Ha_2Fa_ij.{n,i,j,k,l} ) +. \n", + " sum mos_ij' (fun n -> m_2122_2Ha_2Fa_ij2.{n,i,j,k,l} ) +. \n", + " sum mos_novirt (fun n -> m_2122_2Ha_2Fa_nv.{n,i,j,k,l} ) \n", " else\n", - " m_2112_1H_2Fb.{i,j,k,l} +.\n", - " sum mos_i (fun n ->\n", - " m_2112_2Ha_2Fb.{i,j,k,l,n} +. m_2112_2Hb_2Fb.{j,i,l,k,n}\n", - " ) +.\n", - " sum mos_i' (fun n ->\n", - " m_2112_2Hb_2Fb.{i,j,k,l,n} +. m_2112_2Ha_2Fb.{j,i,l,k,n}\n", - " ) +.\n", - " m_2121_2Hb_1F.{i,j,k,l} +.\n", - " sum mos_j (fun n ->\n", - " m_2121_2Hb_2Fa.{i,j,k,l,n} +. m_2121_2Hb_2Fb.{j,i,l,k,n}\n", - " ) +.\n", - " sum mos_j' (fun n ->\n", - " m_2121_2Hb_2Fb.{i,j,k,l,n} +. m_2121_2Hb_2Fa.{j,i,l,k,n}\n", - " ) +.\n", - " sum mos_cabs (fun a ->\n", - " sum mos_virt' (fun m ->\n", - " h_two i j a m s s' *. f_two a m k l s s' ) +.\n", - " sum mos_virt (fun m ->\n", - " h_two i j m a s s' *. f_two m a k l s s' ) +.\n", - " sum mos_ij (fun n ->\n", - " h_two n i a k s s *. f_two a j n l s s'\n", - " +. h_two n j a l s s' *. f_two i a k n s s \n", - " -. h_two n j k a s s' *. f_two i a n l s s'\n", - " ) +.\n", - " sum mos_ij' (fun n -> if n >= j then 0. else \n", - " h_two i n k a s s' *. f_two j a l n s' s'\n", - " +. h_two n j a l s' s' *. f_two i a k n s s' ) +. \n", - " sum mos_ij' (fun n -> if n <= j then 0. else \n", - " -. h_two i n k a s s' *. f_two j a n l s' s' \n", - " -. h_two j n a l s' s' *. f_two i a k n s s' ) +.\n", - " sum mos_ij' (fun n -> \n", - " -. h_two i n a l s s' *. f_two a j k n s s'\n", - " )\n", - " )\n", + " m_2112_1H_2Fb.{i,j,k,l} +. \n", + " sum mos_i (fun n -> m_2112_2Ha_2Fb.{n,i,j,k,l} ) +.\n", + " sum mos_i' (fun n -> m_2112_2Ha_2Fb.{n,j,i,l,k} ) +.\n", + " sum mos_j (fun n -> m_2121_2Hb_2Fa.{n,i,j,k,l} ) +.\n", + " sum mos_j' (fun n -> m_2121_2Hb_2Fa.{n,j,i,l,k} ) +.\n", + " sum mos_novirt'(fun n -> m_2122_2Hb_2Fb_nv.{n,i,j,k,l} ) +.\n", + " sum mos_novirt (fun n -> m_2122_2Hb_2Fb_nv2.{n,i,j,k,l} )+.\n", + " sum mos_ij (fun n -> m_2122_2Hb_2Fb_ij.{n,i,j,k,l} ) +. \n", + " sum mos_ij' (fun n -> m_2122_2Hb_2Fb_ij2.{n,i,j,k,l} ) \n", " in\n", "\n", " match phase with\n", " | Phase.Pos -> result\n", " | Phase.Neg -> -. result\n", " \n", + "let _ = check 100 integral_value 2 \n", "\n", - "let _ = check 100 integral_value 2 1 \n", - "(*\n", - "\n", - "let ki = det_I.(2)\n", - "let kj = det_I.(33)\n", + "let ki = det_I.(7)\n", + "let kj = det_I.(89)\n", "\n", "let _ = integral_value ki kj\n", "let _ = \n", - " let alphas = generate_alphas ki kj 2 1\n", + " let alphas = generate_alphas ki kj 2 \n", " in compute_HaaF ki alphas kj\n", - "*)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# 12. 2 2 " - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "let integral_value ki kj = \n", - " let h, h', p, p', s, s', phase =\n", - " match Excitation.of_det ki kj with\n", - " | Excitation.(Double (phase,\n", - " { hole=h ; particle=p ; spin=s },\n", - " { hole=h'; particle=p'; spin=s'}) ) -> h, h', p, p', s, s', phase\n", - " | _ -> assert false\n", - " in\n", - "\n", - " let result = \n", - " if s <> s' then (* Alpha-Beta *)\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", - " 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", - " | Phase.Neg -> -. result\n", - " \n", - "\n", - "let _ = check 100 integral_value 2 2 2 2" + "\n" ] }, { @@ -1704,9 +1925,48 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 130, "metadata": {}, - "outputs": [], + "outputs": [ + { + "data": { + "text/plain": [ + "val integral_value : Determinant.t -> Determinant.t -> float = \n" + ] + }, + "execution_count": 130, + "metadata": {}, + "output_type": "execute_result" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Checking ... \n", + " - 10 %\n", + " - 20 %\n", + " - 30 %\n", + " - 40 %\n", + " - 50 %\n", + " - 60 %\n", + " - 70 %\n", + " - 80 %\n", + " - 90 %\n", + " - 100 %\n", + "OK: (13, 52) | 6.122906e-04 6.122906e-04 | 3.252607e-19 | cpu : 77.238524 0.119682\n" + ] + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 130, + "metadata": {}, + "output_type": "execute_result" + } + ], "source": [ "let integral_value ki kj = \n", " let i, j, m, k, l, n, s1, s2, s3, phase =\n", @@ -1765,7 +2025,14 @@ " | Phase.Neg -> -. result\n", " \n", "\n", - "let _ = check 200 integral_value 3 " + "let _ = check 100 integral_value 3 " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "---------------------------" ] }, { @@ -1909,13 +2176,6 @@ " let x = (integral_value ki kj) in\n", " Printf.printf \"%20.8e %20.8e %20.8e\\n%!\" x v (v-. x)\n" ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] } ], "metadata": {