From b5637661fa80c3277df196d52a8ceb4cc0a38106 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Wed, 15 Nov 2023 13:10:14 +0100 Subject: [PATCH 1/4] Fix shell_index when converting file --- bin/qp_convert_output_to_ezfio | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/bin/qp_convert_output_to_ezfio b/bin/qp_convert_output_to_ezfio index 091423e4..32721e1e 100755 --- a/bin/qp_convert_output_to_ezfio +++ b/bin/qp_convert_output_to_ezfio @@ -224,14 +224,18 @@ def write_ezfio(res, filename): exponent += [p.expo for p in b.prim] ang_mom.append(str.count(s, "z")) shell_prim_num.append(len(b.prim)) - shell_index += [nshell_tot+1] * len(b.prim) + shell_index += [nshell_tot] * len(b.prim) + + shell_num = len(ang_mom) + assert(shell_index[0] = 1) + assert(shell_index[-1] = shell_num) # ~#~#~#~#~ # # W r i t e # # ~#~#~#~#~ # ezfio.set_basis_basis("Read from ResultsFile") - ezfio.set_basis_shell_num(len(ang_mom)) + ezfio.set_basis_shell_num(shell_num) ezfio.set_basis_basis_nucleus_index(nucl_index) ezfio.set_basis_prim_num(len(coefficient)) From ce43b16fc0c7f76023b1744bdfd4dcd23a8aee50 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Wed, 28 Feb 2024 14:46:23 +0100 Subject: [PATCH 2/4] Fixed bug in PT2 with fast stochastic convergence --- src/cipsi/run_pt2_slave.irp.f | 8 +++++++- src/tools/diagonalize_h.irp.f | 1 + src/tools/print_energy.irp.f | 3 ++- 3 files changed, 10 insertions(+), 2 deletions(-) diff --git a/src/cipsi/run_pt2_slave.irp.f b/src/cipsi/run_pt2_slave.irp.f index debae596..cb1dd1f5 100644 --- a/src/cipsi/run_pt2_slave.irp.f +++ b/src/cipsi/run_pt2_slave.irp.f @@ -186,6 +186,7 @@ subroutine run_pt2_slave_large(thread,iproc,energy) type(pt2_type) :: pt2_data integer :: n_tasks, k, N integer :: i_generator, subset + integer :: ifirst integer :: bsize ! Size of selection buffers logical :: sending @@ -202,6 +203,7 @@ subroutine run_pt2_slave_large(thread,iproc,energy) zmq_socket_push = new_zmq_push_socket(thread) + ifirst = 0 b%N = 0 buffer_ready = .False. n_tasks = 1 @@ -250,7 +252,11 @@ subroutine run_pt2_slave_large(thread,iproc,energy) call omp_set_lock(global_selection_buffer_lock) global_selection_buffer%mini = b%mini call merge_selection_buffers(b,global_selection_buffer) - b%cur=0 + if (ifirst /= 0 ) then + b%cur=0 + else + ifirst = 1 + endif call omp_unset_lock(global_selection_buffer_lock) if ( iproc == 1 ) then call omp_set_lock(global_selection_buffer_lock) diff --git a/src/tools/diagonalize_h.irp.f b/src/tools/diagonalize_h.irp.f index c9ae2033..ffc53aa2 100644 --- a/src/tools/diagonalize_h.irp.f +++ b/src/tools/diagonalize_h.irp.f @@ -20,4 +20,5 @@ subroutine routine call diagonalize_CI print*,'N_det = ',N_det call save_wavefunction_general(N_det,N_states,psi_det_sorted,size(psi_coef_sorted,1),psi_coef_sorted) + call print_mol_properties end diff --git a/src/tools/print_energy.irp.f b/src/tools/print_energy.irp.f index 4fe1572c..0e67828e 100644 --- a/src/tools/print_energy.irp.f +++ b/src/tools/print_energy.irp.f @@ -14,5 +14,6 @@ end subroutine run implicit none - print *, psi_energy + nuclear_repulsion + call print_mol_properties + print *, psi_energy + nuclear_repulsion end From 30096e07ea33720d6104c13f43f9eb5788967ddc Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Wed, 28 Feb 2024 15:02:39 +0100 Subject: [PATCH 3/4] Merging TC with cipsi --- .../cipsi_tc_bi_ortho/run_pt2_slave.irp.f | 45 +++++++++---------- .../run_selection_slave.irp.f | 29 ++++++------ .../cipsi_tc_bi_ortho/zmq_selection.irp.f | 3 +- src/cipsi/run_selection_slave.irp.f | 2 +- 4 files changed, 39 insertions(+), 40 deletions(-) diff --git a/plugins/local/cipsi_tc_bi_ortho/run_pt2_slave.irp.f b/plugins/local/cipsi_tc_bi_ortho/run_pt2_slave.irp.f index aa6546e7..d4f45649 100644 --- a/plugins/local/cipsi_tc_bi_ortho/run_pt2_slave.irp.f +++ b/plugins/local/cipsi_tc_bi_ortho/run_pt2_slave.irp.f @@ -31,11 +31,12 @@ subroutine run_pt2_slave(thread,iproc,energy) double precision, intent(in) :: energy(N_states_diag) integer, intent(in) :: thread, iproc - if (N_det > 100000 ) then - call run_pt2_slave_large(thread,iproc,energy) - else - call run_pt2_slave_small(thread,iproc,energy) - endif + call run_pt2_slave_large(thread,iproc,energy) +! if (N_det > 100000 ) then +! call run_pt2_slave_large(thread,iproc,energy) +! else +! call run_pt2_slave_small(thread,iproc,energy) +! endif end subroutine run_pt2_slave_small(thread,iproc,energy) @@ -178,15 +179,12 @@ subroutine run_pt2_slave_large(thread,iproc,energy) type(pt2_type) :: pt2_data integer :: n_tasks, k, N integer :: i_generator, subset + integer :: ifirst integer :: bsize ! Size of selection buffers logical :: sending - double precision :: time_shift - PROVIDE global_selection_buffer global_selection_buffer_lock - call random_number(time_shift) - time_shift = time_shift*15.d0 zmq_to_qp_run_socket = new_zmq_to_qp_run_socket() @@ -198,15 +196,13 @@ subroutine run_pt2_slave_large(thread,iproc,energy) zmq_socket_push = new_zmq_push_socket(thread) + ifirst = 0 b%N = 0 buffer_ready = .False. n_tasks = 1 sending = .False. done = .False. - double precision :: time0, time1 - call wall_time(time0) - time0 = time0+time_shift do while (.not.done) integer, external :: get_tasks_from_taskserver @@ -233,28 +229,29 @@ subroutine run_pt2_slave_large(thread,iproc,energy) ASSERT (b%N == bsize) endif + double precision :: time0, time1 + call wall_time(time0) call pt2_alloc(pt2_data,N_states) b%cur = 0 call select_connected(i_generator,energy,pt2_data,b,subset,pt2_F(i_generator)) + call wall_time(time1) integer, external :: tasks_done_to_taskserver if (tasks_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id,n_tasks) == -1) then done = .true. endif call sort_selection_buffer(b) - - call wall_time(time1) -! if (time1-time0 > 15.d0) then - call omp_set_lock(global_selection_buffer_lock) - global_selection_buffer%mini = b%mini - call merge_selection_buffers(b,global_selection_buffer) - b%cur=0 - call omp_unset_lock(global_selection_buffer_lock) - call wall_time(time0) -! endif - call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) - if ( iproc == 1 .or. i_generator < 100 .or. done) then + call omp_set_lock(global_selection_buffer_lock) + global_selection_buffer%mini = b%mini + call merge_selection_buffers(b,global_selection_buffer) + if (ifirst /= 0 ) then + b%cur=0 + else + ifirst = 1 + endif + call omp_unset_lock(global_selection_buffer_lock) + if ( iproc == 1 ) then call omp_set_lock(global_selection_buffer_lock) call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), global_selection_buffer, (/task_id/), 1,sending) global_selection_buffer%cur = 0 diff --git a/plugins/local/cipsi_tc_bi_ortho/run_selection_slave.irp.f b/plugins/local/cipsi_tc_bi_ortho/run_selection_slave.irp.f index d351cc79..39c83c4b 100644 --- a/plugins/local/cipsi_tc_bi_ortho/run_selection_slave.irp.f +++ b/plugins/local/cipsi_tc_bi_ortho/run_selection_slave.irp.f @@ -5,19 +5,22 @@ subroutine run_selection_slave(thread, iproc, energy) implicit none - double precision, intent(in) :: energy(N_states) - integer, intent(in) :: thread, iproc + double precision, intent(in) :: energy(N_states) + integer, intent(in) :: thread, iproc + integer :: rc, i - integer :: rc, i - integer :: worker_id, task_id(1), ctask, ltask - character*(512) :: task - integer(ZMQ_PTR) :: zmq_to_qp_run_socket - integer(ZMQ_PTR) :: zmq_socket_push - integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket - integer(ZMQ_PTR), external :: new_zmq_push_socket - type(selection_buffer) :: buf, buf2 - type(pt2_type) :: pt2_data - logical :: done, buffer_ready + integer :: worker_id, task_id(1), ctask, ltask + character*(512) :: task + + integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket + integer(ZMQ_PTR) :: zmq_to_qp_run_socket + + integer(ZMQ_PTR), external :: new_zmq_push_socket + integer(ZMQ_PTR) :: zmq_socket_push + + type(selection_buffer) :: buf, buf2 + logical :: done, buffer_ready + type(pt2_type) :: pt2_data PROVIDE psi_bilinear_matrix_columns_loc psi_det_alpha_unique psi_det_beta_unique PROVIDE psi_bilinear_matrix_rows psi_det_sorted_tc_order psi_bilinear_matrix_order @@ -64,7 +67,7 @@ subroutine run_selection_slave(thread, iproc, energy) stop '-1' end if end if - call select_connected(i_generator, energy, pt2_data, buf,subset, pt2_F(i_generator)) + call select_connected(i_generator, energy, pt2_data, buf, subset, pt2_F(i_generator)) endif integer, external :: task_done_to_taskserver diff --git a/plugins/local/cipsi_tc_bi_ortho/zmq_selection.irp.f b/plugins/local/cipsi_tc_bi_ortho/zmq_selection.irp.f index dc3e0f27..22db643f 100644 --- a/plugins/local/cipsi_tc_bi_ortho/zmq_selection.irp.f +++ b/plugins/local/cipsi_tc_bi_ortho/zmq_selection.irp.f @@ -11,7 +11,7 @@ subroutine ZMQ_selection(N_in, pt2_data) integer, external :: omp_get_thread_num type(pt2_type), intent(inout) :: pt2_data - PROVIDE psi_det psi_coef N_det qp_max_mem N_states pt2_F s2_eig N_det_generators +! PROVIDE psi_det psi_coef N_det qp_max_mem N_states pt2_F s2_eig N_det_generators N = max(N_in,1) N = min(N, (elec_alpha_num * (mo_num-elec_alpha_num))**2) @@ -61,7 +61,6 @@ subroutine ZMQ_selection(N_in, pt2_data) ipos=1 task = ' ' - do i= 1, N_det_generators do j=1,pt2_F(i) write(task(ipos:ipos+30),'(I9,1X,I9,1X,I9,''|'')') j, i, N diff --git a/src/cipsi/run_selection_slave.irp.f b/src/cipsi/run_selection_slave.irp.f index 91bd3a38..87ebca40 100644 --- a/src/cipsi/run_selection_slave.irp.f +++ b/src/cipsi/run_selection_slave.irp.f @@ -65,7 +65,7 @@ subroutine run_selection_slave(thread,iproc,energy) stop '-1' end if end if - call select_connected(i_generator,energy,pt2_data,buf,subset,pt2_F(i_generator)) + call select_connected(i_generator, energy, pt2_data, buf, subset, pt2_F(i_generator)) endif integer, external :: task_done_to_taskserver From 992732813881397b9e854381f803f0056b6616ba Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Wed, 28 Feb 2024 18:15:25 +0100 Subject: [PATCH 4/4] Changed dummy into ghost --- ocaml/Angmom.ml | 3 +-- ocaml/Basis.ml | 2 +- ocaml/Element.ml | 38 ++++++++++++++++++------------------ ocaml/qp_create_ezfio.ml | 42 ++++++++++++++++++++-------------------- 4 files changed, 42 insertions(+), 43 deletions(-) diff --git a/ocaml/Angmom.ml b/ocaml/Angmom.ml index ed13e8dc..2da09340 100644 --- a/ocaml/Angmom.ml +++ b/ocaml/Angmom.ml @@ -26,8 +26,7 @@ let of_string = function | "J" | "j" -> J | "K" | "k" -> K | "L" | "l" -> L - | x -> raise (Failure ("Angmom should be S|P|D|F|G|H|I|J|K|L, -not "^x^".")) + | x -> raise (Failure ("Angmom should be S|P|D|F|G|H|I|J|K|L, not "^x^".")) let of_char = function | 'S' | 's' -> S diff --git a/ocaml/Basis.ml b/ocaml/Basis.ml index 9b0c6a38..f951a5f3 100644 --- a/ocaml/Basis.ml +++ b/ocaml/Basis.ml @@ -17,7 +17,7 @@ let read in_channel at_number = (** Find an element in the basis set file *) let find in_channel element = seek_in in_channel 0; - let element_read = ref Element.X in + let element_read = ref Element.Og in while !element_read <> element do let buffer = input_line in_channel in diff --git a/ocaml/Element.ml b/ocaml/Element.ml index f0d4455d..a794b2bb 100644 --- a/ocaml/Element.ml +++ b/ocaml/Element.ml @@ -4,7 +4,7 @@ open Qptypes exception ElementError of string type t = X - + |H |He |Li|Be |B |C |N |O |F |Ne |Na|Mg |Al|Si|P |S |Cl|Ar @@ -20,7 +20,7 @@ type t = X let of_string x = match (String.capitalize_ascii (String.lowercase_ascii x)) with -| "X" | "Dummy" -> X +| "X" | "Ghost" -> X | "H" | "Hydrogen" -> H | "He" | "Helium" -> He | "Li" | "Lithium" -> Li @@ -265,7 +265,7 @@ let to_string = function let to_long_string = function -| X -> "Dummy" +| X -> "Ghost" | H -> "Hydrogen" | He -> "Helium" | Li -> "Lithium" @@ -492,20 +492,20 @@ let to_charge c = | No -> 102 | Lr -> 103 | Rf -> 104 - | Db -> 105 - | Sg -> 106 - | Bh -> 107 - | Hs -> 108 - | Mt -> 109 - | Ds -> 110 - | Rg -> 111 - | Cn -> 112 - | Nh -> 113 - | Fl -> 114 - | Mc -> 115 - | Lv -> 116 - | Ts -> 117 - | Og -> 118 + | Db -> 105 + | Sg -> 106 + | Bh -> 107 + | Hs -> 108 + | Mt -> 109 + | Ds -> 110 + | Rg -> 111 + | Cn -> 112 + | Nh -> 113 + | Fl -> 114 + | Mc -> 115 + | Lv -> 116 + | Ts -> 117 + | Og -> 118 in Charge.of_int result @@ -565,7 +565,7 @@ let of_charge c = match (Charge.to_int c) with | 52 -> Te | 53 -> I | 54 -> Xe -| 55 -> Cs +| 55 -> Cs | 56 -> Ba | 57 -> La | 58 -> Ce @@ -880,7 +880,7 @@ let vdw_radius x = | Ts -> None | Og -> None in - match result x with + match result x with | Some y -> Some (Positive_float.of_float @@ Units.angstrom_to_bohr *. y ) | None -> None diff --git a/ocaml/qp_create_ezfio.ml b/ocaml/qp_create_ezfio.ml index 8e452762..4e17c0ad 100644 --- a/ocaml/qp_create_ezfio.ml +++ b/ocaml/qp_create_ezfio.ml @@ -6,8 +6,8 @@ type element = | Element of Element.t | Int_elem of (Nucl_number.t * Element.t) -(** Handle dummy atoms placed on bonds *) -let dummy_centers ~threshold ~molecule ~nuclei = +(** Handle ghost atoms placed on bonds *) +let ghost_centers ~threshold ~molecule ~nuclei = let d = Molecule.distance_matrix molecule in @@ -68,11 +68,11 @@ let run ?o b au c d m p cart xyz_file = (Molecule.of_file xyz_file ~charge:(Charge.of_int c) ~multiplicity:(Multiplicity.of_int m) ) in - let dummy = - dummy_centers ~threshold:d ~molecule ~nuclei:molecule.Molecule.nuclei + let ghost = + ghost_centers ~threshold:d ~molecule ~nuclei:molecule.Molecule.nuclei in let nuclei = - molecule.Molecule.nuclei @ dummy + molecule.Molecule.nuclei @ ghost in @@ -145,8 +145,6 @@ let run ?o b au c d m p cart xyz_file = | i :: k :: [] -> (Nucl_number.of_int @@ int_of_string i, Element.of_string k) | _ -> failwith "Expected format is int,Element:basis" in Int_elem result - and basis = - String.lowercase_ascii basis in let key = match elem with @@ -313,7 +311,7 @@ let run ?o b au c d m p cart xyz_file = } in let nuclei = - molecule.Molecule.nuclei @ dummy + molecule.Molecule.nuclei @ ghost in @@ -491,11 +489,7 @@ let run ?o b au c d m p cart xyz_file = |> List.rev |> list_map (fun (x,i) -> try - let e = - match x.Atom.element with - | Element.X -> Element.H - | e -> e - in + let e = x.Atom.element in let key = Int_elem (i,x.Atom.element) in @@ -507,9 +501,15 @@ let run ?o b au c d m p cart xyz_file = in try Basis.read_element (basis_channel key) i e - with Not_found -> - failwith (Printf.sprintf "Basis not found for atom %d (%s)" (Nucl_number.to_int i) - (Element.to_string x.Atom.element) ) + with _ -> + try + if e = Element.X then + Basis.read_element (basis_channel key) i (Element.H) + else + raise Not_found + with Not_found -> + failwith (Printf.sprintf "Basis not found for atom %d (%s)" (Nucl_number.to_int i) + (Element.to_string x.Atom.element) ) with | End_of_file -> failwith ("Element "^(Element.to_string x.Atom.element)^" not found in basis set.") @@ -710,9 +710,9 @@ If a file with the same name as the basis set exists, this file will be read. O arg=With_arg ""; doc="Total charge of the molecule. Default is 0. For negative values, use m instead of -, for ex m1"} ; - { opt=Optional ; short='d'; long="dummy"; + { opt=Optional ; short='g'; long="ghost"; arg=With_arg ""; - doc="Add dummy atoms. x * (covalent radii of the atoms)."} ; + doc="Add ghost atoms. x * (covalent radii of the atoms)."} ; { opt=Optional ; short='m'; long="multiplicity"; arg=With_arg ""; @@ -756,8 +756,8 @@ If a file with the same name as the basis set exists, this file will be read. O int_of_string x ) in - let dummy = - match Command_line.get "dummy" with + let ghost = + match Command_line.get "ghost" with | None -> 0. | Some x -> float_of_string x in @@ -782,7 +782,7 @@ If a file with the same name as the basis set exists, this file will be read. O | x::_ -> x in - run ?o:output basis au charge dummy multiplicity pseudo cart xyz_filename + run ?o:output basis au charge ghost multiplicity pseudo cart xyz_filename ) with (* | Failure txt -> Printf.eprintf "Fatal error: %s\n%!" txt *)