1!dalton_copyright_start
2!
3!
4!dalton_copyright_end
5
6#ifdef VAR_MPI
7module sync_coworkers
8
9! stefan: - this module provides all necessary functionality
10!           to synchronize the co-workers in parallel mcscf/ci
11!           calculations.
12!
13!           written by sknecht for DALTON, december 2010.
14  use dalton_mpi
15  use lucita_cfg
16  use lucita_mcscf_ci_cfg
17  use vector_xc_file_type
18  use file_type_module, only : file_type, file_info
19#ifdef USE_MPI_MOD_F90
20  use mpi
21  implicit none
22#else
23  implicit none
24#include "mpif.h"
25#endif
26
27#if defined (VAR_INT64)
28#define my_MPI_INTEGER MPI_INTEGER8
29#else
30#define my_MPI_INTEGER MPI_INTEGER4
31#endif
32
33  public sync_coworkers_cfg
34  public set_sync_default
35
36  private
37
38  save
39
40  integer(kind=MPI_INTEGER_KIND), private   :: istat(MPI_STATUS_SIZE)
41  integer(kind=MPI_INTEGER_KIND), private   :: ierr
42  integer                       , parameter :: my_comm_world = MPI_COMM_WORLD
43
44  integer, public, parameter :: dim_sync_ctrl_array                    = 6
45
46  logical, public            :: sync_ctrl_array(1:dim_sync_ctrl_array) = .false.
47
48#ifdef PRG_DIRAC
49#include "dgroup.h"
50#endif
51
52contains
53
54  subroutine sync_coworkers_cfg(xarray1,xarray2)
55!******************************************************************************
56!
57!    purpose:  synchronize (if necessary) co-workers for parallel CI/MCSCF
58!
59!*******************************************************************************
60    real(8), optional, intent(inout) :: xarray1(*)
61    real(8), optional, intent(inout) :: xarray2(*)
62!-------------------------------------------------------------------------------
63    integer                          :: select_sync
64!-------------------------------------------------------------------------------
65
66      select_sync = 0
67
68      do ! loop over synchronization processes
69
70        select_sync = select_sync + 1
71
72        if(select_sync > dim_sync_ctrl_array) exit
73
74        select case(sync_ctrl_array(select_sync))
75          case(.true.)
76!           select synchronization module
77#ifdef LUCI_DEBUG
78            print *,' *** select_sync = ***',select_sync
79#endif
80            select case(select_sync)
81              case(1) ! CI default settings and MCSCF/CI input variables
82                call sync_coworkers_ci_cfg()
83              case(2) ! MCSCF environment settings
84                call sync_coworkers_mc_cfg()
85              case(3) ! synchronize 1-el (i|j) and 2-el (ab|cd) integrals
86                if(present(xarray1) .and. present(xarray2))then
87                  call sync_coworkers_ij_abcd(xarray1,xarray2)
88                else
89                  call quit('*** sync_coworkers_cfg: sync of 1-/2-el integrals requires'//          &
90                            ' the output arrays in the argument list.***')
91                end if
92              case(4) ! distribute previous solution vector (i.e., we restart a CI run)
93              case(5) ! update information for vector exchange
94                call sync_coworkers_mc_vector_xc()
95              case(6) ! set current CI task
96                call sync_coworkers_citask()
97              case default ! no option found, quit.
98                print *,' *** sync_coworkers_cfg: synchronization module does not exist, program exits. ***'
99                call quit('*** sync_coworkers_cfg: unknown synchronization module.***')
100            end select
101          case default ! .false.
102#ifdef LUCI_DEBUG
103!           issue a warning
104            print *,' *** sync_coworkers_cfg: synchronization requested but not needed, program continues. ***'
105#endif
106        end select
107
108!       reset control array
109        sync_ctrl_array(select_sync) = .false.
110
111      end do
112
113  end subroutine sync_coworkers_cfg
114
115  subroutine sync_coworkers_ci_cfg()
116!*******************************************************************************
117!
118!    purpose:  provide co-workers with basic common block/orbital space
119!              knowledge in parallel CI/MCSCF runs.
120!
121!*******************************************************************************
122!-------------------------------------------------------------------------------
123
124!     part 1 - the dynamic variables
125!     ------------------------------
126
127!     logical block
128      call dalton_mpi_bcast(lucita_cfg_analyze_cvec,        0, my_comm_world)
129      call dalton_mpi_bcast(lucita_cfg_natural_orb_occ_nr,  0, my_comm_world)
130      call dalton_mpi_bcast(lucita_cfg_transition_densm,    0, my_comm_world)
131      call dalton_mpi_bcast(lucita_cfg_initialize_cb,       0, my_comm_world)
132
133!     real(8) block
134      call dalton_mpi_bcast(lucita_cfg_accepted_truncation, 0, my_comm_world)
135      call dalton_mpi_bcast(lucita_cfg_convergence_c      , 0, my_comm_world)
136      call dalton_mpi_bcast(lucita_cfg_auxiliary_conv_c   , 0, my_comm_world)
137
138!     integer block
139      call dalton_mpi_bcast(lucita_cfg_csym,                0, my_comm_world)
140      call dalton_mpi_bcast(lucita_cfg_hcsym,               0, my_comm_world)
141      call dalton_mpi_bcast(lucita_cfg_icstate,             0, my_comm_world)
142      call dalton_mpi_bcast(lucita_cfg_spin1_2el_op,        0, my_comm_world)
143      call dalton_mpi_bcast(lucita_cfg_spin2_2el_op,        0, my_comm_world)
144      call dalton_mpi_bcast(lucita_cfg_ptg_symmetry,        0, my_comm_world)
145      call dalton_mpi_bcast(lucita_cfg_el_operator_level,   0, my_comm_world)
146      call dalton_mpi_bcast(lucita_cfg_nr_roots,            0, my_comm_world)
147      call dalton_mpi_bcast(lucita_cfg_global_print_lvl,    0, my_comm_world)
148      call dalton_mpi_bcast(lucita_cfg_local_print_lvl,     0, my_comm_world)
149      call dalton_mpi_bcast(lucita_cfg_density_calc_lvl,    0, my_comm_world)
150      call dalton_mpi_bcast(lucita_cfg_spindensity_calc_lvl,0, my_comm_world)
151      call dalton_mpi_bcast(lucita_cfg_is_spin_multiplett,  0, my_comm_world)
152      call dalton_mpi_bcast(lucita_cfg_max_dav_subspace_dim,0, my_comm_world)
153      call dalton_mpi_bcast(lucita_cfg_max_nr_dav_ci_iter,  0, my_comm_world)
154
155!     part 2 - the static variables (update only if necessary, i.e. SETCI has been called)
156!     ------------------------------------------------------------------------------------
157
158      if(lucita_cfg_initialize_cb)then
159
160#ifdef PRG_DIRAC
161!       transfer information from Dirac cb to LUCITA
162        lucita_cfg_nr_ptg_irreps = nbsym
163#endif
164
165!       character
166        call dalton_mpi_bcast(lucita_cfg_run_title,           0, my_comm_world)
167        call dalton_mpi_bcast(lucita_cfg_ci_type,             0, my_comm_world)
168!       real(8)
169        call dalton_mpi_bcast(lucita_cfg_core_energy        , 0, my_comm_world)
170!       logical
171        call dalton_mpi_bcast(lucita_cfg_inactive_shell_set,  0, my_comm_world)
172        call dalton_mpi_bcast(lucita_cfg_minmax_occ_gas_set,  0, my_comm_world)
173        call dalton_mpi_bcast(lucita_cfg_timing_par,          0, my_comm_world)
174!       integer
175!       a. general definitions/settings
176        call dalton_mpi_bcast(lucita_cfg_nr_active_e,         0, my_comm_world)
177        call dalton_mpi_bcast(lucita_cfg_restart_ci,          0, my_comm_world)
178        call dalton_mpi_bcast(lucita_cfg_max_batch_size,      0, my_comm_world)
179        call dalton_mpi_bcast(lucipar_cfg_ttss_dist_strategy, 0, my_comm_world)
180        call dalton_mpi_bcast(lucipar_cfg_mem_reduction_multp,0, my_comm_world)
181!       call dalton_mpi_bcast(lucita_cfg_nr_calc_sequences,   0, my_comm_world)
182
183!       b. symmetry, spin and orbital spaces
184        call dalton_mpi_bcast(lucita_cfg_nr_ptg_irreps,       0, my_comm_world)
185        call dalton_mpi_bcast(lucita_cfg_nr_gas_spaces,       0, my_comm_world)
186        call dalton_mpi_bcast(lucita_cfg_init_wave_f_type,    0, my_comm_world)
187        call dalton_mpi_bcast(nish_lucita,                    0, my_comm_world)
188        call dalton_mpi_bcast(naos_lucita,                    0, my_comm_world)
189        call dalton_mpi_bcast(nmos_lucita,                    0, my_comm_world)
190        call dalton_mpi_bcast(nash_lucita,                    0, my_comm_world)
191        call dalton_mpi_bcast(nfro_lucita,                    0, my_comm_world)
192        call dalton_mpi_bcast(nocc_lucita,                    0, my_comm_world)
193
194        call dalton_mpi_bcast(ngsh_lucita,                    0, my_comm_world)
195        call dalton_mpi_bcast(ngso_lucita,                    0, my_comm_world)
196
197        select case(lucita_cfg_init_wave_f_type)
198          case(1) ! GAS settings
199!           nothing particular here
200          case(2) ! RAS settings
201!           logical
202!           call dalton_mpi_bcast(lucita_cfg_ras1_set,        0, my_comm_world)
203!           call dalton_mpi_bcast(lucita_cfg_ras2_set,        0, my_comm_world)
204!           call dalton_mpi_bcast(lucita_cfg_ras3_set,        0, my_comm_world)
205!           integer
206!           call dalton_mpi_bcast(nas1_lucita,                0, my_comm_world)
207!           call dalton_mpi_bcast(nas2_lucita,                0, my_comm_world)
208!           call dalton_mpi_bcast(nas3_lucita,                0, my_comm_world)
209            call dalton_mpi_bcast(lucita_cfg_min_e_ras1,      0, my_comm_world)
210            call dalton_mpi_bcast(lucita_cfg_max_e_ras1,      0, my_comm_world)
211            call dalton_mpi_bcast(lucita_cfg_min_e_ras3,      0, my_comm_world)
212            call dalton_mpi_bcast(lucita_cfg_max_e_ras3,      0, my_comm_world)
213        end select
214      end if
215
216  end subroutine sync_coworkers_ci_cfg
217!******************************************************************************
218
219  subroutine sync_coworkers_mc_cfg()
220!*******************************************************************************
221!
222!    purpose:  provide co-workers with the basic settings from the mcscf
223!              environment.
224!
225!*******************************************************************************
226#ifdef MOD_SRDFT
227      use lucita_mcscf_srdftci_cfg
228#endif
229!-------------------------------------------------------------------------------
230
231!     logical
232
233      call dalton_mpi_bcast(docisrdft_mc2lu,                  0,my_comm_world)
234      srdft_ci_with_lucita = docisrdft_mc2lu
235      call dalton_mpi_bcast(integrals_from_mcscf_env,         0,my_comm_world)
236      call dalton_mpi_bcast(mcscf_ci_update_ijkl,             0,my_comm_world)
237      call dalton_mpi_bcast(mcscf_orbital_trial_vector,       0,my_comm_world)
238      call dalton_mpi_bcast(mcscf_ci_trial_vector,            0,my_comm_world)
239      call dalton_mpi_bcast(io2io_vector_exchange_mc2lu_lu2mc,0,my_comm_world)
240      call dalton_mpi_bcast(vector_update_mc2lu_lu2mc,        0,my_comm_world)
241
242!     real(8)
243!     call dalton_mpi_bcast(einact_mc2lu,                     0, my_comm_world) ! no need to sync - we sync ecore/ecore_orig in lucita internally
244
245!     integer
246      call dalton_mpi_bcast(len_cref_mc2lu,                   0, my_comm_world)
247      call dalton_mpi_bcast(len_hc_mc2lu,                     0, my_comm_world)
248      call dalton_mpi_bcast(len_resolution_mat_mc2lu,         0, my_comm_world)
249      call dalton_mpi_bcast(len_int1_or_rho1_mc2lu,           0, my_comm_world)
250      call dalton_mpi_bcast(len_int2_or_rho2_mc2lu,           0, my_comm_world)
251      call dalton_mpi_bcast(vector_exchange_type1,            0, my_comm_world)
252      call dalton_mpi_bcast(vector_exchange_type2,            0, my_comm_world)
253
254      call dalton_mpi_bcast(file_info%current_file_nr_active1,0, my_comm_world)
255      call dalton_mpi_bcast(file_info%current_file_nr_active2,0, my_comm_world)
256
257  end subroutine sync_coworkers_mc_cfg
258!******************************************************************************
259
260  subroutine sync_coworkers_mc_vector_xc()
261!*******************************************************************************
262!
263!    purpose:  provide co-workers with the basic information to initiate
264!              the vector exchange process.
265!
266!*******************************************************************************
267
268!     logical
269      call dalton_mpi_bcast(exchange_f_info%exchange_file_io2io, 0,my_comm_world)
270
271!     integer
272      call dalton_mpi_bcast(exchange_f_info%present_sym_irrep,   0, my_comm_world)
273      call dalton_mpi_bcast(exchange_f_info%push_pull_switch,    0, my_comm_world)
274      call dalton_mpi_bcast(exchange_f_info%total_nr_vectors,    0, my_comm_world)
275
276      call dalton_mpi_bcast(file_info%current_file_nr_active1,   0, my_comm_world)
277      call dalton_mpi_bcast(file_info%current_file_nr_active2,   0, my_comm_world)
278
279  end subroutine sync_coworkers_mc_vector_xc
280!******************************************************************************
281
282  subroutine sync_coworkers_ij_abcd(xarray1,xarray2)
283!******************************************************************************
284!
285!    purpose:  provide co-workers with 1-/2-electron integrals in parallel
286!              CI/MCSCF runs.
287!
288!******************************************************************************
289    use lucita_energy_types
290    real(8), intent(inout) :: xarray1(*)
291    real(8), intent(inout) :: xarray2(*)
292!-------------------------------------------------------------------------------
293!   common block from lucita (note: the one below does not have its own include file)
294    integer    :: i12s,i34s,i1234s,nint1,nint2,nbint1,nbint2
295    COMMON/CINTFO/i12s,i34s,i1234s,nint1,nint2,nbint1,nbint2
296!#include "parluci.h"
297!-------------------------------------------------------------------------------
298
299!     total # 1-/2-electron integrals
300      call dalton_mpi_bcast(nint1,      0, my_comm_world)
301      call dalton_mpi_bcast(nint2,      0, my_comm_world)
302!     core energy + inactive energy and core energy
303      call dalton_mpi_bcast(ecore,      0, my_comm_world)
304      call dalton_mpi_bcast(ecore_orig, 0, my_comm_world)
305
306!     synchronize the 1-electron integrals (use the generic bcast because
307!     the 1-/2-electron integral arrays have too large dimensions in the calling
308!     lucita routines (for historic reasons/other purposes, i do not know. stefan dec 2010)
309
310!     call mpi_bcast(xarray1, nint1, mpi_real8, 0, my_comm_world, ierr)
311!     if(.not.mcscf_ci_trial_vector) call mpi_bcast(xarray2, nint2, mpi_real8, 0, my_comm_world, ierr)
312!     31-jan-2020 hjaaj: try again with dalton_mpi_bcast ...
313      call dalton_mpi_bcast(xarray1(1:nint1), 0, my_comm_world)
314      if(.not.mcscf_ci_trial_vector) call dalton_mpi_bcast(xarray2(1:nint2), 0, my_comm_world)
315
316!     set sync_ctrl option
317      sync_ctrl_array(3) = .false.
318
319  end subroutine sync_coworkers_ij_abcd
320!******************************************************************************
321
322  subroutine sync_coworkers_citask()
323!*******************************************************************************
324!
325!    purpose:  provide co-workers with knowledge about the parallel CI run task.
326!
327!*******************************************************************************
328!-------------------------------------------------------------------------------
329
330!     character
331      call dalton_mpi_bcast(lucita_ci_run_id, 0, my_comm_world)
332
333  end subroutine sync_coworkers_citask
334!******************************************************************************
335
336  subroutine set_sync_default(set_value_sync_ctrl_array,sync_ctrl_array_pos)
337!*******************************************************************************
338!
339!    purpose:  set/reset the co-workers synchronization control array.
340!
341!*******************************************************************************
342    logical, intent(in) :: set_value_sync_ctrl_array
343    integer, intent(in) :: sync_ctrl_array_pos
344!-------------------------------------------------------------------------------
345    logical             :: set_sync_default_val
346!-------------------------------------------------------------------------------
347
348      set_sync_default_val = set_value_sync_ctrl_array
349
350      call dalton_mpi_bcast(set_sync_default_val, 0, my_comm_world)
351
352!     insert value in ctrl-array
353      sync_ctrl_array(sync_ctrl_array_pos) = set_sync_default_val
354
355  end subroutine set_sync_default
356!******************************************************************************
357
358end module
359#else
360subroutine sync_coworkers
361! dummy routine for non-mpi compilation
362end
363#endif
364