From 1896d9c45f3c51e9733e36ab7e9289b67b82b9e4 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Fri, 29 Mar 2019 19:01:51 +0100 Subject: [PATCH] Reduced RAM in F12 --- CI/DeterminantSpace.ml | 5 -- CI/F12CI.ml | 161 +++++++---------------------------------- 2 files changed, 27 insertions(+), 139 deletions(-) diff --git a/CI/DeterminantSpace.ml b/CI/DeterminantSpace.ml index 48b1d02..024e32b 100644 --- a/CI/DeterminantSpace.ml +++ b/CI/DeterminantSpace.ml @@ -270,11 +270,6 @@ let arbitrary_of_mo_basis mo_basis f = let n_det_beta = Array.length det_beta in let n_det_alfa = Array.length det_alfa in - let ndet = n_det_alfa * n_det_beta in - if Parallel.master then - Format.printf "Number of determinants : %d %d %d\n%!" - n_det_alfa n_det_beta ndet; - let det = Array.make n_det_alfa (Array.init n_det_beta (fun i -> i)) in diff --git a/CI/F12CI.ml b/CI/F12CI.ml index 5027222..69384ee 100644 --- a/CI/F12CI.ml +++ b/CI/F12CI.ml @@ -86,76 +86,41 @@ let is_internal det_space = | 1 | 2 -> false | _ -> true - (* - if not (Bitstring.logand aux_mask alfa |> Bitstring.is_zero ) then - false - else - let beta = - Determinant.beta a - |> Spindeterminant.bitstring - in - Bitstring.logand aux_mask beta |> Bitstring.is_zero - *) - let dressing_vector ~frozen_core aux_basis f12_amplitudes ci = -(* - let i_o1_alfa = h_ij aux_basis in - let alfa_o2_i = f_ij aux_basis in - - let w_alfa _ _ = 1. in - - let mo_class = CI.mo_class ci in - - let list_holes = List.concat - [ MOClass.inactive_mos mo_class ; MOClass.active_mos mo_class ] - and list_particles2 = MOClass.auxiliary_mos mo_class - and list_particles1 = List.concat - [ MOClass.active_mos mo_class ; MOClass.virtual_mos mo_class ; MOClass.auxiliary_mos mo_class ] - in - (* Single state here *) - let result = - CI.second_order_sum ci list_holes list_particles1 list_holes list_particles2 - (is_internal ci.det_space) i_o1_alfa alfa_o2_i w_alfa f12_amplitudes ~unique:false - |> Vec.of_list - in - Matrix.sparse_of_vector_array [| Vector.sparse_of_vec result |] - -*) - - Printf.printf "Building matrix\n%!"; + if Parallel.master then + Printf.printf "Building matrix\n%!"; (* Determinants of the FCI space as a list *) let in_dets = - DeterminantSpace.determinants_array ci.CI.det_space - |> Array.to_list + DeterminantSpace.determinant_stream ci.CI.det_space + |> Util.stream_to_list in (* Stream that generates only singly and doubly excited determinants wrt FCI space *) let out_dets_stream = + (* Stream that generates all determinants of FCI space *) + let s = + DeterminantSpace.fci_of_mo_basis ~frozen_core aux_basis + |> DeterminantSpace.determinant_stream - (* Stream that generates all determinants of FCI space *) - let s = - DeterminantSpace.fci_of_mo_basis ~frozen_core aux_basis - |> DeterminantSpace.determinant_stream - - in - (* Select only singly and doubly excited determinants - wrt FCI space *) - Stream.from (fun _ -> - try - let rec result () = - let ki = Stream.next s in - if is_internal ci.CI.det_space ki then - result () - else - Some ki - in - result () - with Stream.Failure -> None - ) + in + (* Select only singly and doubly excited determinants + wrt FCI space *) + Stream.from (fun _ -> + try + let rec result () = + let ki = Stream.next s in + if is_internal ci.CI.det_space ki then + result () + else + Some ki + in + result () + with Stream.Failure -> None + ) in let make_h_and_f alpha_list = @@ -187,9 +152,8 @@ let dressing_vector ~frozen_core aux_basis f12_amplitudes ci = in - Printf.printf "Matrix product\n%!"; let m_HF = - let batch_size = 1 + 10_000_000 / (Mat.dim1 f12_amplitudes) in + let batch_size = 1 + 1_000_000 / (Mat.dim1 f12_amplitudes) in let input_stream = Stream.from (fun i -> let rec make_batch accu = function @@ -213,7 +177,7 @@ let dressing_vector ~frozen_core aux_basis f12_amplitudes ci = in let iteration input = - Printf.printf "gemm\n%!"; + Printf.printf ".%!"; let m_H_aux, m_F_aux = make_h_and_f input in let m_HF = gemm m_H_aux m_F_aux ~transb:`T @@ -225,44 +189,10 @@ let dressing_vector ~frozen_core aux_basis f12_amplitudes ci = |> Farm.run ~ordered:false ~f:iteration |> Stream.iter (fun hf -> ignore @@ Mat.add result hf ~c:result ); - result + Printf.printf "\n"; + Parallel.broadcast (lazy result) in - -(* - let m_HF = - - let in_dets = - DeterminantSpace.determinants_array ci.CI.det_space - in - - let fci_space = DeterminantSpace.fci_of_mo_basis ~frozen_core aux_basis in - Array.mapi (fun i ki -> - Printf.printf "%d / %d\r%!" i (Array.length in_dets); - Array.map (fun kj -> - match Determinant.degrees ki kj with - | (0,0) | (0,1) | (0,2) | (0,3) | (0,4) - | (1,0) | (2,0) | (3,0) | (4,0) - | (1,1) | (2,1) | (3,1) - | (1,2) | (1,3) - | (2,2) -> - ( let x = ref 0. in - DeterminantSpace.determinant_stream fci_space - |> Stream.iter (fun k_alfa -> - if not (is_internal ci.CI.det_space ki) then - let f = f_ij aux_basis k_alfa kj in - if f <> 0. then - let h = h_ij aux_basis ki k_alfa in - x := !x +. f *. h - ); - !x ) - | _ -> 0. - ) in_dets - ) in_dets - |> Mat.of_array - in - *) - Printf.printf "Done\n%!"; Matrix.dense_of_mat m_HF @@ -316,43 +246,6 @@ Printf.printf "det space\n%!"; Parallel.broadcast (lazy x) in - - (* - let f12_amplitudes = - fun c -> - let result = lacpy c in - Mat.scal (0.5) result ; - result - (* While in a sequential region, initiate the parallel - 4-idx transformation to avoid nested parallel jobs - *) - ignore @@ MOBasis.f12_ints aux_basis; - - let two_gamma_inv = -2. *. f12.F12.expo_s in - - let f = fun ki kj -> - if ki <> kj then - (f_ij aux_basis ki kj) - else - two_gamma_inv +. (f_ij aux_basis ki kj) - in - let m_F = - CI.create_matrix_spin f det_space - |> Lazy.force - in - fun ci_coef -> - let result = - Matrix.ax_eq_b m_F (Matrix.dense_of_mat ci_coef) - |> Matrix.to_mat - in - let norm = gemm ~transa:`T result result in - Printf.printf "Norm of |F> : %f\n%!" norm.{1,1}; - if abs_float norm.{1,1} > 1. then - failwith "Norm of |F> > 1"; - result - in - *) - let e_shift = let det = DeterminantSpace.determinant_stream det_space