MPI-SCATCI  2.0
An MPI version of SCATCI
Options_module.f90
Go to the documentation of this file.
1 ! Copyright 2019
2 !
3 ! For a comprehensive list of the developers that contributed to these codes
4 ! see the UK-AMOR website.
5 !
6 ! This file is part of UKRmol-in (UKRmol+ suite).
7 !
8 ! UKRmol-in is free software: you can redistribute it and/or modify
9 ! it under the terms of the GNU General Public License as published by
10 ! the Free Software Foundation, either version 3 of the License, or
11 ! (at your option) any later version.
12 !
13 ! UKRmol-in is distributed in the hope that it will be useful,
14 ! but WITHOUT ANY WARRANTY; without even the implied warranty of
15 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 ! GNU General Public License for more details.
17 !
18 ! You should have received a copy of the GNU General Public License
19 ! along with UKRmol-in (in source/COPYING). Alternatively, you can also visit
20 ! <https://www.gnu.org/licenses/>.
21 
41 module options_module
42 
43  use const_gbl, only: line_len, stdin, stdout
44  use consts_mpi_ci, only: mxint, no_ci_target, normal_csf, generate_target_wf, diagonalize, default_archer_memory, &
45  scatci_decide, ecp_type_null, save_all_coeffs, no_diagonalization, diagonalize_no_restart, &
46  normal_phase, name_len_max, do_ci_target_contract, prototype_csf, symtype_dinfh
47  use mpi_gbl, only: master, mpiint, mpi_xermsg
48  use precisn, only: longint, wp
49  use scatci_routines, only: rdnfto, movep, mkpt, ham_pars_io
50  use memorymanager_module, only: master_memory
51  use parallelization_module, only: grid => process_grid
52 
53  implicit none
54 
55  public optionsset, options, phase, ham_pars_io_wrapper
56 
61  type phase
62  integer :: size_phaze
63  integer, allocatable :: phase_index(:)
64  end type phase
65 
66 
74  type :: options
75  !--------------SCATCI-INPUTS----------------------------!
76  integer :: megul
77  integer :: ci_target_flag
78  integer :: csf_type
79  integer :: ci_target_switch
80  real(wp) :: exotic_mass
81  integer :: diagonalization_flag
82  integer :: integral_unit
83  integer :: hamiltonian_unit
84  integer :: ci_output_unit
85  integer :: num_matrix_elements_per_rec
86  integer :: phase_correction_flag
87  integer :: num_eigenpairs
88  logical :: use_ukrmol_integrals
89  logical :: quantamoln
90  integer :: integral_ordering
91  integer :: print_dataset_heading
92  integer :: output_ci_set_number
93  integer(longint) :: total_memory
94  integer :: diagonalizer_choice
95  integer :: use_scf
96  integer :: force_serial
97  integer :: print_flags(2)
98 
100  character(len=NAME_LEN_MAX) :: name
101  character(len=NAME_LEN_MAX) :: diag_name
102 
103  !------------Header information on CSF------------------------------!
104  integer :: lambda = 0
105  real(wp) :: spin = 0.0_wp
106  real(wp) :: spin_z = 0.0_wp
107  real(wp) :: refl_quanta = 0.0_wp
108  real(wp) :: pin = 0.0_wp
109  integer :: tot_num_orbitals = 0
110  integer :: tot_num_spin_orbitals = 0
111  integer :: num_csfs = 0
112  integer :: num_electrons = 0
113  integer :: length_coeff = 0
114  integer :: matrix_eval = -1
115  integer :: num_syms = 0
116  integer :: sym_group_flag = 0
117  integer :: length_dtrs = 0
118  integer :: positron_flag = 0
119  integer :: output_flag(6)
120  real(wp) :: threshold = 0.0_wp
121  integer :: num_continuum = 0
122  integer :: maximum_num_target_states = 0
123  integer :: total_num_ci_target = 0
124  integer :: total_num_target_states = 0
125  integer :: num_target_sym = 0
126  integer :: last_continuum = 0
127  integer :: num_expanded_continuum_csf = 0
128  integer :: seq_num_last_continuum_csf = 0
129  integer :: continuum_block_size = 0
130  integer :: num_l2_csf = 0
131  integer :: actual_num_csfs = 0
132  integer :: contracted_mat_size = 0
133  integer :: total_num_tgt_states = 0
134 
135  !---------------allocatable header information----------------!
136  integer, allocatable :: num_orbitals_sym(:)
137  integer, allocatable :: num_electronic_orbitals_sym(:)
138  integer, allocatable :: num_positron_orbitals_sym(:)
139  integer, allocatable :: num_unused_orbitals_sym(:)
140  integer, allocatable :: reference_dtrs(:)
141  integer, allocatable :: num_dtr_csf_out(:)
142  type(phase), allocatable :: ci_phase(:)
143  integer, allocatable :: phase_index(:)
144  integer, allocatable :: num_orbitals_sym_dinf(:)
145  integer, allocatable :: num_target_orbitals_sym_congen(:)
146  integer, allocatable :: num_target_orbitals_sym_dinf_congen(:)
147  integer, allocatable :: num_target_orbitals_sym(:)
148  integer, allocatable :: num_ci_target_sym(:)
149  integer, allocatable :: num_continuum_orbitals_target(:)
150  integer, allocatable :: num_target_state_sym(:)
151  integer, allocatable :: lambda_continuum_orbitals_target(:)
152  integer, allocatable :: gerade_sym_continuum_orbital_target(:)
153  integer, allocatable :: orbital_sequence_number(:)
154  integer, allocatable :: ci_set_number(:)
155  integer, allocatable :: ci_sequence_number(:)
156  real(wp), allocatable :: energy_shifts(:)
157 
158  !ECP variables
159  integer, allocatable :: target_multiplicity(:)
160  integer, allocatable :: target_spatial(:)
161  logical :: multiplicity_defined = .false.
162  logical :: spatial_defined = .false.
163  logical :: all_ecp_defined = .false.
164  integer :: ecp_type = ecp_type_null
165  character(len=line_len*5) :: ecp_filename
166  integer :: exclude_rowcolumn = -1
167 
168  ! Eigenvector output write order
169  integer(mpiint) :: preceding_writer = -1
170  integer(mpiint) :: following_writer = -1
171 
172  !----------------------miscellaneous variables
173  integer :: temp_orbs(20)
174 
175  real(wp) :: eigenvec_conv = 1d-10
176  real(wp) :: residual_conv = 1d-8
177  real(wp) :: orthog_trigger = 1d-7
178  real(wp) :: max_tolerance = 1d-12
179  integer :: max_iterations = -1
180 
181  integer :: vector_storage_method = save_all_coeffs
182  contains
183  procedure, public :: read => read_input
184  procedure, private :: read_congen
185  procedure, private :: compute_variables
186  procedure, private :: compute_expansions
187  procedure, private :: check_sanity_congen_header
188  procedure, private :: check_sanity_complete
189  procedure, private :: arrange_phase
190  procedure, public :: do_ci_contraction
191  end type options
192 
193 
202  type :: optionsset
203  type(options), allocatable :: opts(:)
204  integer, allocatable :: idtarg(:) ! denprop-to-congen state order permutation
205  logical :: write_amp = .false. ! write boundary amplitudes (swinterf "fort.21")
206  logical :: write_dip = .false. ! write transition dipole moments (cdenprop_target "fort.624")
207  logical :: write_rmt = .false. ! write RMT input file (rmt_interface "molecular_data")
208  logical :: all_props = .true. ! whether to calculate IRR-with-itself properties and energy-sort states (.true. for RMT)
209  integer :: luprop = -1 ! orbital properties file unit to read if different than NFTI
210  integer :: nfdm ! number of extra R-matrix sub-spheres for RMT (excluding the main R-matrix sphere)
211  real(wp) :: rmatr ! radius of the R-matrix sphere
212  real(wp) :: delta_r ! spacing between the R-matrix sub-spheres
213  character(len=1) :: cform = 'F' ! formatted or unformatted channel description file
214  character(len=1) :: rform = 'U' ! formatted or unformatted R-matrix file
215  contains
216  procedure, public :: read => read_all_namelists
217  procedure, private :: read_all_input
218  procedure, private :: read_all_cinorn
219  procedure, private :: read_outer_interface
220  procedure, public :: setup_write_order
221  end type optionsset
222 
223 contains
224 
233  subroutine read_all_namelists (this)
234 
235  class(optionsset), intent(inout) :: this
236 
237  integer :: num_arguments, io, err, num_input, num_cinorn, i
238  character(len=800) :: input_filename
239 
240  num_arguments = command_argument_count()
241 
242  if (num_arguments > 1) then
243  write (stdout, "('Only one argument allowed')")
244  write (stdout, "('Command line must be:')")
245  write (stdout, "('./mpi-scatci <input filename>')")
246  call mpi_xermsg('Options_module', 'read_all_namelists', 'Too many command-line arguments', 1, 1)
247  end if
248 
249  if (num_arguments == 1) then
250  call get_command_argument(1, input_filename)
251  input_filename = trim(input_filename)
252 
253  open (newunit = io, action = 'read', status = 'old', file = input_filename, iostat = err)
254  write (stdout, "('Opened input file in unit = ',i4)") io
255 
256  if (err /= 0) then
257  call mpi_xermsg('Options_module', 'read_all_namelists', 'Could not open input file ' // trim(input_filename), 2, 1)
258  end if
259  else
260  io = stdin
261  end if
262 
263  if (allocated(this % opts)) deallocate(this % opts)
264 
265  num_input = this % read_all_input(io)
266  num_cinorn = this % read_all_cinorn(io)
267 
268  call this % read_outer_interface(io)
269 
270  if (io /= stdin) close (io)
271 
272  ! we require at least one valid input
273  if (num_input == 0) then
274  call mpi_xermsg('Options_module', 'read_all_namelists', 'Missing &INPUT namelist in the input file', 3, 1)
275  end if
276 
277  ! we require equal number of &INPUT and &CINORN namelists
278  if (num_input /= num_cinorn) then
279  call mpi_xermsg('Options_module', 'read_all_namelists', 'Unequal number of &INPUT and &CINORN namelists', 4, 1)
280  end if
281 
282  ! let each Options instance finalize itself
283  do i = 1, num_input
284  call this % opts(i) % read
285  end do
286 
287  end subroutine read_all_namelists
288 
289 
302  integer function read_all_input (this, io) result (i)
303  class(optionsset), intent(inout) :: this
304  integer, intent(in) :: io
305 
306  type(options), allocatable :: tmp_opts(:)
307  character(NAME_LEN_MAX) :: name
308 
309  logical :: qmoln, ukrmolp_ints
310  integer :: icidg, icitg, iexpc, idiag, nfti, nfte, lembf, megul, nftg, ntgsym, ncont, ncorb, iord, scfuse, lusme
311  integer :: io_stat, ierr, nset, nrec, nnuc, nocsf, nstat, mgvn, nelt, k, j, iposit
312 
313  character(1000) :: io_msg
314  real(wp) :: scalem, thrhm, spin, spinz, e0
315  integer, allocatable :: numtgt(:), notgt(:), nctgt(:), ntgtf(:), ntgts(:), mcont(:), gucont(:), nobc(:)
316 
317  namelist /input/ icidg, icitg, iexpc, idiag, nfti, nfte, lembf, megul, nftg, ntgsym, numtgt, notgt, nctgt, ntgtf, ntgts, &
318  mcont, gucont, ncont, ncorb, nobc, iord, name, scalem, scfuse, thrhm, qmoln, ukrmolp_ints, lusme, iposit
319 
320  i = 0
321 
322  allocate (numtgt(mxint), notgt(mxint), nctgt(mxint), ntgtf(mxint), ntgts(mxint), mcont(mxint), gucont(mxint), &
323  nobc(mxint), stat = ierr)
324 
325  if (ierr /= 0) then
326  call mpi_xermsg('Options_module', 'read_all_input', 'Input arrays allocation failure', 1, 1)
327  end if
328 
329  rewind(io, iostat = io_stat)
330 
331  namelist_loop: do
332 
333  ! set namelist entries to their default values
334  icidg = diagonalize
335  icitg = no_ci_target
336  iexpc = normal_csf
337  idiag = -1
338  nfti = 16
339  nfte = 26
340  lembf = 5000
341  megul = 13
342  nftg = generate_target_wf
343  ntgsym = 0
344  numtgt(:) = 0
345  notgt(:) = 0
346  nctgt(:) = 1
347  ntgtf(:) = 0
348  ntgts(:) = 0
349  mcont(:) = 0
350  ncorb = normal_phase
351  iord = 1
352  name = ''
353  scalem = 0
354  scfuse = 0
355  iposit = 0
356  qmoln = .false.
357  ukrmolp_ints = .true.
358 
359  ! read the next &INPUT namelist
360  read (io, nml = input, iostat = io_stat, iomsg = io_msg, end = 991)
361  if (io_stat /= 0) then
362  call mpi_xermsg('Options_module', 'read_all_input', &
363  'Error when reading &INPUT namelist: ' // trim(io_msg), i + 1, 1)
364  end if
365 
366  ! seems to be a valid namelist
367  i = i + 1
368 
369  ! re-allocate array of Options structures
370  if (.not. allocated(this % opts) .or. size(this % opts) < i) then
371  if (allocated(this % opts)) call move_alloc(this % opts, tmp_opts)
372  allocate (this % opts(i), stat = ierr)
373  if (ierr /= 0) then
374  call mpi_xermsg('Options_module', 'read_all_input', 'Input structures allocation failure', 2, 1)
375  end if
376  if (allocated(tmp_opts)) this % opts(1:size(tmp_opts)) = tmp_opts(:)
377  end if
378 
379  ! prevent overflows in below assigments (but should be assured by length check of the namelist data already)
380  if (ntgsym > mxint) then
381  call mpi_xermsg('Options_module', 'read_all_input', 'NTGSYM exceeds built-in limit MXINT', 3, 1)
382  end if
383  this % opts(i) % total_num_tgt_states = sum(numtgt(1:ntgsym))
384  if (this % opts(i) % total_num_tgt_states > mxint) then
385  call mpi_xermsg('Options_module', 'read_all_input', 'Number of target states exceeds built-in limit MXINT', 4, 1)
386  end if
387 
388  ! copy data from the namelist to the structure
389  this % opts(i) % diagonalization_flag = icidg
390  this % opts(i) % ci_target_flag = icitg
391  this % opts(i) % csf_type = iexpc
392  this % opts(i) % matrix_eval = idiag
393  this % opts(i) % integral_unit = nfti
394  this % opts(i) % hamiltonian_unit = nfte
395  this % opts(i) % num_matrix_elements_per_rec = lembf
396  this % opts(i) % megul = megul
397  this % opts(i) % ci_target_switch = nftg
398  this % opts(i) % num_target_sym = ntgsym
399 
400  if (ntgsym > 0) then
401  this % opts(i) % num_target_state_sym = numtgt(1:ntgsym)
402  this % opts(i) % num_continuum_orbitals_target = notgt(1:ntgsym)
403  this % opts(i) % ci_set_number = ntgtf(1:this % opts(i) % total_num_tgt_states)
404  this % opts(i) % ci_sequence_number = ntgts(1:this % opts(i) % total_num_tgt_states)
405  this % opts(i) % lambda_continuum_orbitals_target = mcont(1:ntgsym)
406  end if
407 
408  this % opts(i) % phase_correction_flag = ncorb
409  this % opts(i) % integral_ordering = iord
410  this % opts(i) % name = name
411  this % opts(i) % exotic_mass = scalem
412  this % opts(i) % use_SCF = scfuse
413  this % opts(i) % positron_flag = iposit
414  this % opts(i) % QuantaMolN = qmoln
415  this % opts(i) % use_UKRMOL_integrals = ukrmolp_ints
416 
417  if ( this % opts(i) % positron_flag /= 0 ) then
418  call mpi_xermsg('Options_module', 'read_all_input', &
419  'Positrons are not yet implemented for mpi-scatci.', i, 1)
420  end if
421 
422  ! also read spin-symmetry information about the target states from the target CI expansion coefficients file (if given)
423  if (nftg > 0) then
424  allocate (this % opts(i) % target_spatial(this % opts(i) % total_num_tgt_states), &
425  this % opts(i) % target_multiplicity(this % opts(i) % total_num_tgt_states))
426  do j = 1, this % opts(i) % total_num_tgt_states
427  if (this % opts(i) % ci_set_number(j) <= 0) then
428  call mpi_xermsg('Options_module', 'read_all_input', &
429  'Non-zero NFTG requires setting also non-zero NTGTF/NTGTS.', i, 1)
430  end if
431  call movep(nftg, this % opts(i) % ci_set_number(j), ierr, 0, 0)
432  if (ierr /= 0) then
433  call mpi_xermsg('Options_module', 'read_all_input', &
434  'CI data not found in unit', nftg, 1)
435  end if
436  read (nftg) ! dummy read to skip header
437  read (nftg) nset, nrec, name, nnuc, nocsf, nstat, mgvn, spin, spinz, nelt, e0
438  this % opts(i) % target_spatial(j) = mgvn
439  this % opts(i) % target_multiplicity(j) = nint(2 * spin + 1)
440  end do
441  end if
442 
443  end do namelist_loop
444 
445  991 return ! no more namelists in the file...
446 
447  end function read_all_input
448 
449 
462  integer function read_all_cinorn (this, io) result (i)
463  class(optionsset), intent(inout) :: this
464  integer, intent(in) :: io
465 
466  type(options), allocatable :: tmp_opts(:)
467 
468  integer :: nfta, nftw, nciset, npflg(2), npcvc, itgt, nopvec, ntgt, notgt, ncipfg, nkey, large, thrprt
469  integer :: keycsf, nstat, igh, maxiter, forse, exrc, vecstore, ecp, targmul(mxint), targspace(mxint)
470  integer :: ierr, io_stat
471  real(wp) :: critc, critr, ortho, crite, memp, eshift(mxint)
472 
473  character(256) :: io_msg
474  character(NAME_LEN_MAX) :: name
475 
476  namelist /cinorn/ nfta, nftw, nciset, npflg, npcvc, name, itgt, nopvec, ntgt, notgt, eshift, ncipfg, nkey, large, thrprt, &
477  keycsf, nstat, igh, critc, critr, ortho, crite, maxiter, memp, forse, exrc, vecstore, ecp, &
478  targmul, targspace
479 
480  i = 0
481 
482  rewind(io, iostat = io_stat)
483 
484  namelist_loop: do
485 
486  ! set namelist entries to their default values
487  nftw = 25
488  nciset = 1
489  npcvc = 1
490  name = ''
491  nstat = 0
492  igh = scatci_decide
493  critc = 1e-10_wp
494  critr = 1e-08_wp
495  ortho = 1e-07_wp
496  crite = 1e-12_wp
497  maxiter = -1
498  memp = default_archer_memory
499  forse = 0
500  exrc = -1
501  vecstore = save_all_coeffs
502  ecp = ecp_type_null
503  targmul = 0
504  targspace = 0
505  eshift = 0.0_wp
506  npflg = 0
507 
508  ! read the next &CINORN namelist
509  read (io, nml = cinorn, iostat = io_stat, iomsg = io_msg, end = 991)
510  if (io_stat /= 0) then
511  call mpi_xermsg('Options_module', 'read_all_input', &
512  'Error when reading &CINORN namelist: ' // trim(io_msg), i + 1, 1)
513  end if
514 
515  ! seems to be a valid namelist
516  i = i + 1
517 
518  ! re-allocate array of Options structures as needed
519  if (.not. allocated(this % opts) .or. size(this % opts) < i) then
520  if (allocated(this % opts)) call move_alloc(this % opts, tmp_opts)
521  allocate (this % opts(i), stat = ierr)
522  if (ierr /= 0) then
523  call mpi_xermsg('Options_module', 'read_all_cinorn', 'Input structures allocation failure', 2, 1)
524  end if
525  if (allocated(tmp_opts)) this % opts(1:size(tmp_opts)) = tmp_opts(:)
526  end if
527 
528  ! copy data from the namelist to the structure
529  this % opts(i) % CI_output_unit = nftw
530  this % opts(i) % output_ci_set_number = nciset
531  this % opts(i) % print_dataset_heading = npcvc
532  this % opts(i) % diag_name = name
533  this % opts(i) % num_eigenpairs = nstat
534  this % opts(i) % diagonalizer_choice = igh
535  this % opts(i) % eigenvec_conv = critc
536  this % opts(i) % residual_conv = critr
537  this % opts(i) % orthog_trigger = ortho
538  this % opts(i) % max_tolerance = crite
539  this % opts(i) % max_iterations = maxiter
540  this % opts(i) % energy_shifts = eshift
541  this % opts(i) % print_flags = npflg
542 
543  ! MPI-SCATCI specific parameters
544  this % opts(i) % total_memory = int(memp * 1024_longint**3, longint)
545  this % opts(i) % force_serial = forse
546  this % opts(i) % exclude_rowcolumn = exrc
547  this % opts(i) % vector_storage_method = vecstore
548  this % opts(i) % ecp_type = ecp
549 
550  end do namelist_loop
551 
552  991 return ! no more namelists in the file...
553 
554  end function read_all_cinorn
555 
556 
566  subroutine read_outer_interface (this, io)
567 
568  class(optionsset), intent(inout) :: this
569  integer, intent(in) :: io
570 
571  logical :: write_amp, write_dip, write_rmt, all_props
572  integer :: luprop, nfdm, io_stat, ntarg, idtarg(mxint)
573  real(wp) :: delta_r, rmatr
574 
575  character(1) :: cform, rform
576  character(256) :: io_msg
577 
578  namelist /outer_interface/ write_amp, write_dip, write_rmt, all_props, luprop, nfdm, delta_r, rmatr, ntarg, idtarg, &
579  cform, rform
580 
581  write_amp = .true.
582  write_dip = .false.
583  write_rmt = .false.
584  all_props = .true.
585  luprop = -1
586  nfdm = 18
587  delta_r = 0.08_wp
588  rmatr = -1.0_wp
589  ntarg = 0
590  idtarg = 0
591  cform = 'F'
592  rform = 'U'
593 
594  rewind(io, iostat = io_stat)
595  read (io, nml = outer_interface, iostat = io_stat, iomsg = io_msg, end = 991)
596  if (io_stat /= 0) then
597  call mpi_xermsg('Options_module', 'read_outer_interface', &
598  'Error when reading &RMT_INTERFACE namelist: ' // trim(io_msg), 1, 1)
599  end if
600 
601  if (rmatr < 0 .and. (write_amp .or. write_rmt)) then
602  call mpi_xermsg('Options_module', 'read_outer_interface', 'Missing "rmatr" in &outer_interface namelist.', 1, 1)
603  end if
604 
605  this % write_amp = write_amp
606  this % write_dip = write_dip
607  this % write_rmt = write_rmt
608  this % all_props = all_props
609  this % luprop = luprop
610  this % nfdm = nfdm
611  this % delta_r = delta_r
612  this % rmatr = rmatr
613  this % cform = cform
614  this % rform = rform
615 
616  if (ntarg > 0) then
617  this % idtarg = idtarg(1:ntarg)
618  else if (this % write_amp .or. this % write_rmt) then
619  if (any(this % opts(:) % ci_target_switch == 0)) then
620  call mpi_xermsg('Options_module', 'read_outer_interface', &
621  'When NFTG = 0, outer interface requires setting NTARG and IDTARG.', 1, 1)
622  end if
623  end if
624 
625  991 return ! no more namelists in the file...
626 
627  end subroutine read_outer_interface
628 
629 
634  subroutine read_input (this)
635  class(options), intent(inout) :: this
636 
637  call this % read_congen
638 
639  call this % compute_variables
640  call this % compute_expansions
641  call this % check_sanity_complete
642 
643  call master_memory % construct(this % total_memory)
644  call master_memory % print_memory_report
645 
646  end subroutine read_input
647 
648 
649  subroutine check_sanity_congen_header (this)
650  class(options), intent(inout) :: this
651 
652  ! write(stdout,"(' ',A6,': ',A)") 'SCATCI',this%name
653  ! write(stdout,16) this%lambda,this%spin,this%num_electrons,this%spin_z,this%matrix_eval,this%refl_quanta,this%threshold,&
654  ! & this%num_syms,this%MEGUL,this%num_csfs,this%output_flag(1:6)
655  ! 16 FORMAT(' MGVN =',I10,5X,'S =',F5.1,/,' NELT =',I10,5X,'SZ =' &
656  ! & ,F5.1,/,' IDIAGT=',I10,5X,'R =',F5.1,/, &
657  ! & ' THRES=',D10.1,5X,'NSYM =',I5,/,' MEGUL=',I3,', NOCSF=' &
658  ! & ,I10,/,' NPFLG',I10,5I4)
659 
660  if (this % num_csfs <= 0) THEN
661  write (stdout, "('Options::check_sanity_congen_header - No CSFS detected...stopping')")
662  stop 'No of CSFS = 0'
663  end if
664 
665  end subroutine check_sanity_congen_header
666 
667 
674  subroutine check_sanity_complete (this)
675  class(options), intent(inout) :: this
676  integer :: I
677 
678  write (stdout, *)
679  write (stdout, "('---------------SCATCI run options------------------')")
680  write (stdout, "('Name of run is ', a)") trim(this % name)
681  write (stdout, "('CONGEN unit ', i4)") this % megul
682  write (stdout, "('Hamiltonian unit ', i4)") this % hamiltonian_unit
683  write (stdout, "('Number of hamiltonian records to write ', i4)") this % num_matrix_elements_per_rec
684  write (stdout, "('Integral unit ', i4)") this % integral_unit
685 
686  write (stdout, "('Integral Type - ')", advance = 'no')
687  if (this % sym_group_flag <= 1) then
688  write (stdout, "('ALCHEMY')")
689  else if (this % sym_group_flag == 2 .and. this % use_UKRMOL_integrals) then
690  write (stdout, "('UKRMOL+')")
691  else
692  write (stdout, "('SWEDEN')")
693  end if
694 
695  write (stdout, "('CI target Contraction? ')", advance = 'no')
696  if (this % ci_target_flag == no_ci_target) then
697  write (stdout, "('No')")
698  else
699  write (stdout, "('Yes')")
700  write (stdout, "('Generate CI vectors? ')", advance = 'no')
701  if (this % ci_target_switch == generate_target_wf) then
702  write (stdout, "('Yes')")
703  else
704  write (stdout, "('No, We are reading from unit ',i4)") this % ci_target_switch
705  end if
706  end if
707 
708  write (stdout, "('Do we diagonalize? ')", advance = 'no')
709  if (this % diagonalization_flag == no_diagonalization) then
710  write (stdout, "('No')")
711  else if (this % diagonalization_flag == diagonalize) then
712  write (stdout, "('Yes with restart for ARPACK')")
713  else if (this % diagonalization_flag == diagonalize_no_restart) then
714  write (stdout, "('Yes with NO restart for ARPACK')")
715  endif
716 
717  write (stdout, "('Num eigenpairs ,',i12)") this % num_eigenpairs
718  write (stdout, *)
719  write (stdout, "(' -----------Molecule options------------ ')")
720  write (stdout, "('Number of electrons = ', i10)") this % num_electrons
721  write (stdout, "('Number of molecular orbitals = ', i10)") this % tot_num_orbitals
722  write (stdout, "('Number of spin orbitals = ', i10)") this % tot_num_spin_orbitals
723  write (stdout, *)
724  write (stdout, "(' -----------CSF options------------ ')")
725  write (stdout, "('Number of configuration state functions = ', i10)") this % num_csfs
726  write (stdout, "('Number of continuum functions = ', i10)") this % last_continuum
727  this % num_L2_CSF = this % num_csfs - this % last_continuum
728  write (stdout, "('Number of L2 functions = ', i10)") this % num_L2_CSF
729  write (stdout, *)
730  write (stdout, "(' -----------Orbital/Symmetry options------------ ')")
731  write (stdout, "('Number of symmetries = ', i4)") this % num_syms
732  write (stdout, "('Number of orbitals per symmetry = ', 20i4)") this % num_orbitals_sym(:)
733  write (stdout, "('Number of electronic orbitals per symmetry = ', 20i4)") this % num_electronic_orbitals_sym(:)
734  write (stdout, "('Number of exotic orbitals per symmetry = ', 20i4)") this % num_positron_orbitals_sym(:)
735 
736  write (stdout, "('Use ECP? ')", advance = 'no')
737  if (this % ecp_type == ecp_type_null) then
738  write (stdout, "('No')")
739  else
740  write (stdout, "('Yes')")
741  end if
742 
743  end subroutine check_sanity_complete
744 
745 
754  subroutine read_congen (this)
755 
756  class(options), intent(inout) :: this
757 
758  !Dummy variables for RDNFTO
759  integer :: dum_nobe, dum_nobp, dum_nobv, iposit, ido, nobep, ntgcon, idiagt
760  character(len=NAME_LEN_MAX) :: congen_name
761 
762  !Begin reading first round of header information
763  rewind(this % megul)
764  read (this % megul) congen_name, this % lambda, this % spin, this % spin_z, this % refl_quanta, this % pin, &
765  this % tot_num_orbitals, this % tot_num_spin_orbitals, this % num_csfs, this % num_electrons, &
766  this % length_coeff,idiagt, this % num_syms, this % sym_group_flag, this % length_dtrs, &
767  this % output_flag, this % threshold, this % num_continuum, ntgcon
768 
769  if (trim(this % name) == '') this % name = congen_name
770  if (this % matrix_eval < 0) this % matrix_eval = idiagt
771  if (ntgcon > 0) this % num_target_sym = ntgcon
772 
773  call this % check_sanity_congen_header()
774 
775  !allocate header arrays
776  call allocate_integer_array(this % num_orbitals_sym, this % num_syms)
777  call allocate_integer_array(this % num_electronic_orbitals_sym, this % num_syms)
778  call allocate_integer_array(this % num_positron_orbitals_sym, this % num_syms)
779  call allocate_integer_array(this % num_unused_orbitals_sym, this % num_syms)
780  call allocate_integer_array(this % reference_dtrs, this % num_electrons)
781  call allocate_integer_array(this % num_dtr_csf_out, this % num_csfs)
782  call allocate_integer_array(this % phase_index, this % num_continuum)
783  call allocate_integer_array(this % num_orbitals_sym_dinf, this % num_syms * 2)
784  call allocate_integer_array(this % num_target_orbitals_sym_congen, this % num_syms)
785  call allocate_integer_array(this % num_target_orbitals_sym, this % num_syms)
786  call allocate_integer_array(this % num_target_orbitals_sym_dinf_congen, this % num_syms * 2)
787  call allocate_integer_array(this % num_ci_target_sym, this % num_target_sym)
788  call allocate_integer_array(this % num_continuum_orbitals_target, this % num_target_sym)
789  call allocate_integer_array(this % lambda_continuum_orbitals_target, this % num_target_sym)
790  call allocate_integer_array(this % gerade_sym_continuum_orbital_target, this % num_target_sym)
791  call allocate_integer_array(this % num_target_state_sym, this % num_target_sym)
792  call allocate_integer_array(this % orbital_sequence_number, this % num_csfs)
793 
794  call rdnfto(this % megul, this % num_orbitals_sym, this % num_target_orbitals_sym_congen, this % num_orbitals_sym_dinf, &
795  this % num_target_orbitals_sym_dinf_congen, this % num_syms, this % reference_dtrs, this % num_electrons, &
796  this % num_dtr_csf_out, this % num_csfs, this % phase_index, this % num_continuum, this % num_ci_target_sym, &
797  this % temp_orbs, this % lambda_continuum_orbitals_target, this % gerade_sym_continuum_orbital_target, &
798  this % num_target_sym, this % num_target_sym, this % num_electronic_orbitals_sym, &
799  this % num_positron_orbitals_sym, this % num_unused_orbitals_sym, this % positron_flag)
800 
801  do ido = 1, this % num_syms
802  if (this % num_positron_orbitals_sym(ido) == 0) then
803  this % num_electronic_orbitals_sym(ido) = this % num_orbitals_sym(ido)
804  write (stdout, *) 'NOB(i)=', this % num_orbitals_sym(ido)
805  write (stdout, *) 'setting NOBE(i)=', this % num_electronic_orbitals_sym(ido)
806  end if
807  nobep = this % num_electronic_orbitals_sym(ido) + this % num_positron_orbitals_sym(ido)
808  if (nobep /= this % num_orbitals_sym(ido)) then
809  write (stdout, *) 'ERROR on input:'
810  write (stdout, *) 'not: NOB(i)=NOBE(i)+NOBP(i)'
811  write (stdout, *) 'i=', ido
812  write (stdout, *) 'NOB(i)=', this % num_orbitals_sym(ido)
813  write (stdout, *) 'NOBE(i)=', this % num_electronic_orbitals_sym(ido)
814  write (stdout, *) 'NOBP(i)=', this % num_positron_orbitals_sym(ido)
815  end if
816  if (this % num_unused_orbitals_sym(ido) == 0) then
817  this % num_unused_orbitals_sym(ido) = this % num_target_orbitals_sym_congen(ido)
818  end if
819  if (this % num_target_orbitals_sym_congen(ido) > this % num_orbitals_sym(ido)) then
820  write (stdout, *) 'ERROR on input:'
821  write (stdout, *) 'NOB0(i) > NOB(i)'
822  write (stdout, *) 'i=', ido
823  write (stdout, *) 'NOB(i)=', this % num_orbitals_sym(ido)
824  write (stdout, *) 'NOB0(i)=', this % num_target_orbitals_sym_congen(ido)
825  end if
826  end do
827 
828  ! Allocate the phases for each target symmetry (must be allocated as it is used in main program)
829  allocate(this % ci_phase(this % num_target_sym))
830 
831  ! Set up the phases
832  !call this % arrange_phase(this % phase_index)
833 
834  end subroutine read_congen
835 
836 
837  subroutine arrange_phase (this, temp_phase)
838  class(options), intent(inout) :: this
839  integer, intent(in) :: temp_phase(:)
840  integer :: ido, ci_num
841 
842  ci_num = 0
843 
844  if (sum(this % num_ci_target_sym) == this % num_continuum) then
845  do ido = 1, this % num_target_sym
846  this % ci_phase(ido) % size_phaze = this % num_ci_target_sym(ido)
847  allocate(this % ci_phase(ido) % phase_index(this % ci_phase(ido) % size_phaze))
848  this % ci_phase(ido) % phase_index(1:this % ci_phase(ido) % size_phaze) = &
849  temp_phase(ci_num + 1 : ci_num + this % ci_phase(ido) % size_phaze)
850  ci_num = ci_num + this % ci_phase(ido) % size_phaze
851  end do
852  else if (this % num_continuum == 0) then
853  return
854  else
855  stop "[Options] Problem with phase not matching sum of NCTGT"
856  end if
857 
858  end subroutine arrange_phase
859 
860 
869  subroutine setup_write_order (this)
870 
871  class(optionsset), intent(inout) :: this
872 
873  integer :: i, j, k, n
874 
875  ! set up preceding/following relations for writing datasets in the same files
876  do i = 1, size(this % opts)
877 
878  n = 0 ! number of datasets written to the same file before this process
879  k = 0 ! which setup index corresponds to the last write to the same file
880 
881  do j = 1, i - 1
882  if (this % opts(j) % CI_output_unit == this % opts(i) % CI_output_unit) then
883  n = n + 1
884  k = j
885  end if
886  end do
887 
888  ! if automatic numbering of datasets is used, get the true unit number for this setup
889  if (this % opts(i) % output_ci_set_number == 0) then
890  this % opts(i) % output_ci_set_number = n + 1
891  end if
892 
893  ! if there is a previous dataset in the same file, link them to avoid lock/overwriting
894  if (k > 0) then
895  if (this % opts(i) % output_ci_set_number <= this % opts(k) % output_ci_set_number) then
896  call mpi_xermsg('Options_module', 'read_all_namelists', &
897  'The entries "nciset" must form ascending sequence.', i, 1)
898  end if
899  this % opts(k) % following_writer = grid % group_master_world_rank(grid % which_group_is_work(i))
900  this % opts(i) % preceding_writer = grid % group_master_world_rank(grid % which_group_is_work(k))
901  end if
902 
903  end do
904 
905  end subroutine setup_write_order
906 
907 
912  subroutine compute_variables (this)
913  class(options), intent(inout) :: this
914  integer :: itgtsym, isym, I, num_orb_temp
915 
916  if (this % csf_type == prototype_csf) then
917  this % num_target_orbitals_sym = 0
918  write (stdout, "('Computing start of contiuum expansion...')", advance = 'no')
919  do itgtsym = 1, this % num_target_sym
920 
921  !If we're working with D_inf_h then we need to 'double' the symmtery
922  if (this % sym_group_flag /= symtype_dinfh) then
923  isym = this % lambda_continuum_orbitals_target(itgtsym) + 1
924  else
925  isym = (this % lambda_continuum_orbitals_target(itgtsym) * 2) + 1
926  if (mod(this % lambda_continuum_orbitals_target(itgtsym), 2) == 0 .and. &
927  this % gerade_sym_continuum_orbital_target(itgtsym) == -1) isym = isym + 1
928  if (mod(this % lambda_continuum_orbitals_target(itgtsym), 2) == 1 .and. &
929  this % gerade_sym_continuum_orbital_target(itgtsym) == 1) isym = isym + 1
930  end if
931 
932  !If there is already a value here then don't bother (impossible since the array must be allocated)
933  if (this % num_target_orbitals_sym(isym) > 0) cycle
934 
935  num_orb_temp = this % num_orbitals_sym(isym) - this % num_continuum_orbitals_target(itgtsym)
936  if (num_orb_temp < this % num_target_orbitals_sym_dinf_congen(isym)) then
937  call mpi_xermsg('Options_module', 'compute_variables', &
938  'Incompatible number of target orbitals. Wrong CONGEN file or NOTGT?', itgtsym, 1)
939  end if
940 
941  this % num_target_orbitals_sym(isym) = num_orb_temp
942 
943  write (stdout, "('Num eigenpairs = ',i8)") this % num_eigenpairs
944  write (stdout, 1010) this % num_target_sym, (this % num_continuum_orbitals_target(i), i = 1, this % num_target_sym)
945  1010 format(/,' Number of target symmetries in expansion, NTGSYM =',i5,/, &
946  ' Number of continuum orbitals for each state, NOTGT =',20i5,/,(' ',20i5))
947  write (stdout, 1020) (this % num_ci_target_sym(i), i = 1, this % num_target_sym)
948  1020 format( ' Number of CI components for each state, NCTGT =',20i10,/,(' ',20i5))
949  write (stdout, 1025) (this % num_target_state_sym(i), i = 1, this % num_target_sym)
950  1025 format( ' Number of target states of each symmetry, NUMTGT =',20i5,/,(' ',20i5))
951  write (stdout, 1030) (this % lambda_continuum_orbitals_target(i), i = 1, this % num_target_sym)
952  1030 format( ' Continuum M projection for each state, MCONT =',20i5,/,(' ',20i5))
953  if (this % gerade_sym_continuum_orbital_target(1) /= 0) then
954  write (stdout, 1040) (this % gerade_sym_continuum_orbital_target(i), i = 1, this % num_target_sym)
955  end if
956  1040 format( ' Continuum G/U symmetry for each state, GUCONT =',20i5,/,(' ',20i5))
957  write (stdout, *)
958 
959  end do
960  else if (this % num_continuum > 0) then
961  this % num_continuum_orbitals_target(1:this%num_target_sym) = this % temp_orbs(1:this%num_target_sym)
962  end if
963 
964  write (stdout, "('done')")
965 
966  this % all_ecp_defined = this % multiplicity_defined .and. this % spatial_defined
967 
968  write (stdout, "(' NOBC =',I10,20I4)") (this % num_target_orbitals_sym(i), i = 1, this % num_syms)
969  write (stdout, "(' NOB =',I10,20I4)") (this % num_orbitals_sym_dinf(i), i = 1, this % num_syms)
970  write (stdout, "(' NOB0 =',I10,20I4)") (this % num_target_orbitals_sym_dinf_congen(i), i = 1 , this % num_syms)
971 
972  end subroutine compute_variables
973 
974 
975  subroutine compute_expansions (this)
976  class(options), intent(inout) :: this
977  integer :: lusme = 0, i
978 
979  if (this % csf_type /= normal_csf .or. this % ci_target_flag /= no_ci_target) then
980 
981  call mkpt(this % orbital_sequence_number, & ! kpt
982  this % num_orbitals_sym_dinf, & ! nob
983  this % num_target_orbitals_sym, & ! nobc
984  this % sym_group_flag, & ! symtype
985  this % num_continuum_orbitals_target, & ! notgt
986  this % num_ci_target_sym, & ! nctgt
987  this % lambda_continuum_orbitals_target, & ! mcont
988  this % gerade_sym_continuum_orbital_target, & ! gucont
989  this % contracted_mat_size, & ! mocsf
990  this % num_csfs, & ! nocsf0
991  this % actual_num_csfs, & ! mxcsf
992  this % last_continuum, & ! ncont
993  this % num_expanded_continuum_CSF, & ! ncont2
994  this % num_target_sym, & ! ntgsym
995  this % total_num_target_states, & ! ntgt
996  this % num_target_state_sym, & ! numtgt
997  this % total_num_ci_target, & ! ncitot
998  this % maximum_num_target_states, & ! nummx
999  this % csf_type, & ! iexpc
1000  this % ci_target_flag, & ! icitg
1001  lusme)
1002 
1003  1050 format(' SCATCI will expand NOCSF =',i7,' prototype CSFs',/15x, &
1004  'into MXCSF =',i7,' actual configurations',/15x, &
1005  'and MOCSF =',i7,' dimension final Hamiltonian')
1006  write (stdout, 1050) this % num_csfs, this % actual_num_csfs, this % contracted_mat_size
1007 
1008  else
1009 
1010  this % contracted_mat_size = this % num_csfs
1011  this % actual_num_csfs = this % num_csfs
1012  this % total_num_target_states = this % num_target_sym
1013 
1014  end if
1015 
1016  write (stdout, "(/' Number of last continuum CSF, NCONT =',i7)") this % last_continuum
1017 
1018  this % seq_num_last_continuum_csf = this % contracted_mat_size - (this % num_csfs - this % last_continuum)
1019  this % continuum_block_size = this % seq_num_last_continuum_csf
1020  if (this % num_eigenpairs == 0) this % num_eigenpairs = this % contracted_mat_size
1021 
1022  end subroutine compute_expansions
1023 
1024 
1025  logical function do_ci_contraction (this)
1026  class(options), intent(inout) :: this
1027 
1028  do_ci_contraction = (this % ci_target_flag == do_ci_target_contract)
1029 
1030  end function do_ci_contraction
1031 
1032 
1033  subroutine allocate_integer_array (array, num_elms)
1034  integer, allocatable :: array(:)
1035  integer :: num_elms, err
1036 
1037  if (.not. allocated(array)) then
1038  allocate(array(num_elms), stat = err)
1039  call master_memory % track_memory(kind(array), size(array), err, 'OPTION::INT_ARRAY')
1040  if (err /= 0) then
1041  write (stdout, "('[Option] Error in allocating one of the arrays')")
1042  stop "[Option] Array allocation error"
1043  end if
1044  array = 0
1045  end if
1046 
1047  end subroutine allocate_integer_array
1048 
1049 
1062  subroutine ham_pars_io_wrapper (option, write_ham_pars, nelms)
1063 
1064  class(options), intent(inout) :: option
1065  logical, intent(in) :: write_ham_pars
1066  integer, intent(inout) :: nelms
1067  integer :: size_phase
1068 
1069  size_phase = size(option % phase_index)
1070 
1071  if (grid % grank /= master) return
1072 
1073  call ham_pars_io(write_ham_pars, &
1074  option % seq_num_last_continuum_csf, &
1075  option % num_target_sym, &
1076  option % num_target_state_sym, &
1077  option % num_continuum_orbitals_target, &
1078  option % lambda_continuum_orbitals_target, &
1079  option % hamiltonian_unit, &
1080  nelms, &
1081  option % sym_group_flag, &
1082  option % spin, &
1083  option % spin_z, &
1084  option % num_electrons, &
1085  option % lambda, &
1086  option % integral_unit, &
1087  option % phase_index, &
1088  size_phase, &
1089  option % use_UKRMOL_integrals)
1090 
1091  end subroutine ham_pars_io_wrapper
1092 
1093 end module options_module
consts_mpi_ci::mxint
integer, parameter mxint
Maximal length of integer arrays in input namelists.
Definition: consts_mpi_ci.f90:40
memorymanager_module
Memory manager module.
Definition: MemoryManager_module.f90:31
consts_mpi_ci::normal_csf
integer, parameter normal_csf
Use configuration state functions as is.
Definition: consts_mpi_ci.f90:67
options_module::ham_pars_io_wrapper
subroutine, public ham_pars_io_wrapper(option, write_ham_pars, nelms)
Read/write information about Hamiltonian.
Definition: Options_module.f90:1063
options_module::read_all_cinorn
integer function read_all_cinorn(this, io)
Read all &CINORN namelists from file unit.
Definition: Options_module.f90:463
options_module::read_all_input
integer function read_all_input(this, io)
Read all &INPUT namelists from file unit.
Definition: Options_module.f90:303
options_module::allocate_integer_array
subroutine allocate_integer_array(array, num_elms)
Definition: Options_module.f90:1034
options_module::check_sanity_complete
subroutine check_sanity_complete(this)
Print a subset of options read from the input.
Definition: Options_module.f90:675
options_module::options
Calculation setup for a single Hamiltonian diagonalization.
Definition: Options_module.f90:74
consts_mpi_ci::no_ci_target
integer, parameter no_ci_target
No Ci target contraction (target only run)
Definition: consts_mpi_ci.f90:61
options_module::read_outer_interface
subroutine read_outer_interface(this, io)
Read contents of the OUTER interface namelist.
Definition: Options_module.f90:567
options_module::compute_expansions
subroutine compute_expansions(this)
Definition: Options_module.f90:976
consts_mpi_ci::scatci_decide
integer, parameter scatci_decide
SCATCI chooses the diagonalizer.
Definition: consts_mpi_ci.f90:99
parallelization_module::process_grid
type(processgrid) process_grid
Definition: Parallelization_module.F90:80
consts_mpi_ci
MPI SCATCI Constants module.
Definition: consts_mpi_ci.f90:31
consts_mpi_ci::symtype_dinfh
integer, parameter symtype_dinfh
This describes D_inf_h symmetries.
Definition: consts_mpi_ci.f90:55
consts_mpi_ci::name_len_max
integer, parameter name_len_max
The maximum length of a name.
Definition: consts_mpi_ci.f90:49
options_module::optionsset
MPI-SCATCI input.
Definition: Options_module.f90:202
parallelization_module
Distribution of processes into a grid.
Definition: Parallelization_module.F90:29
options_module::phase
?
Definition: Options_module.f90:61
consts_mpi_ci::diagonalize_no_restart
integer, parameter diagonalize_no_restart
Diagonalize but with no restart.
Definition: consts_mpi_ci.f90:89
options_module::read_congen
subroutine read_congen(this)
Reads all input information.
Definition: Options_module.f90:755
options_module::read_all_namelists
subroutine read_all_namelists(this)
Read all relevant namelists from the input.
Definition: Options_module.f90:234
options_module::compute_variables
subroutine compute_variables(this)
Computes additional variables required by scatci.
Definition: Options_module.f90:913
options_module::check_sanity_congen_header
subroutine check_sanity_congen_header(this)
Definition: Options_module.f90:650
consts_mpi_ci::default_archer_memory
real(wp), parameter default_archer_memory
A default memory value (in GiB), we use the default given for per proc in archer.
Definition: consts_mpi_ci.f90:46
options_module::read_input
subroutine read_input(this)
Read and check input files.
Definition: Options_module.f90:635
consts_mpi_ci::do_ci_target_contract
integer, parameter do_ci_target_contract
Perform CI target contraction (target+scattering run)
Definition: consts_mpi_ci.f90:63
options_module::arrange_phase
subroutine arrange_phase(this, temp_phase)
Definition: Options_module.f90:838
consts_mpi_ci::generate_target_wf
integer, parameter generate_target_wf
Generate target wavefunction parameter.
Definition: consts_mpi_ci.f90:71
consts_mpi_ci::ecp_type_null
integer, parameter ecp_type_null
Definition: consts_mpi_ci.f90:145
consts_mpi_ci::diagonalize
integer, parameter diagonalize
Diagonalize.
Definition: consts_mpi_ci.f90:87
options_module::setup_write_order
subroutine setup_write_order(this)
Set up write order flags.
Definition: Options_module.f90:870
consts_mpi_ci::save_all_coeffs
integer, parameter save_all_coeffs
Do not discard any coefficients.
Definition: consts_mpi_ci.f90:114
options_module
Options module.
Definition: Options_module.f90:41
options_module::do_ci_contraction
logical function do_ci_contraction(this)
Definition: Options_module.f90:1026
consts_mpi_ci::prototype_csf
integer, parameter prototype_csf
The configuration state functions are prototypes and therefore require expansion.
Definition: consts_mpi_ci.f90:69
consts_mpi_ci::no_diagonalization
integer, parameter no_diagonalization
No diagonalization.
Definition: consts_mpi_ci.f90:85
memorymanager_module::master_memory
type(memorymanager), public master_memory
This is the global memory manager.
Definition: MemoryManager_module.f90:69
consts_mpi_ci::normal_phase
integer, parameter normal_phase
No phase correction.
Definition: consts_mpi_ci.f90:93