Dice Fortran Backend Documentation
driver.f90
Go to the documentation of this file.
1 ! ---------------------------------------------------------------------------------------
2 ! Dice Quantum Monte Carlo
3 ! ---------------------------------------------------------------------------------------
4 ! MAIN DRIVER CODE
5 !
6 ! DESCRIPTION:
9 !
12 ! ---------------------------------------------------------------------------------------
13 
14 program main
15 
16  use shared_data
17  use input_parser
18  use solvers
19  use vqmc
20  use write_netcdf
21  use testing
22  use restart_fns
23  use mpi
24 
25  implicit none
26 
27  call driver()
28 
29 contains
30 
31 ! ---------------------------------------------------------------------------------------
32 ! ROUTINE: driver
33 !
36 !
52 ! ---------------------------------------------------------------------------------------
53 subroutine driver()
54 
55 ! ---------------------------------------------------------------------------------------
56 ! SECTION: Variable Declarations
57 ! ---------------------------------------------------------------------------------------
58 
59  real(real64), dimension(:), allocatable :: p_arr
60  real(real64), dimension(:), allocatable :: bonds
61  real(real64), dimension(:), allocatable :: H2plus_param_c_array
62  real(real64), dimension(:), allocatable :: H2_param_a_array
63  real(real64), dimension(:), allocatable :: H2_param_beta_array
64  real(real64), dimension(:), allocatable :: energies
65  real(real64), dimension(:), allocatable :: variances
66  real(real64), dimension(:), allocatable :: accept_rates
68  real(real64) :: H2plus_param_c_guess
69  real(real64) :: H2plus_param_c_opt
70  real(real64) :: H2_param_beta_guess
71  real(real64) :: H2_param_beta_opt
72  real(real64) :: H2_param_a
73  real(real64) :: energy
74  real(real64) :: variance
75  real(real64) :: accept_rate
76 
77  integer :: b, i
78  integer :: logiostat
79  integer :: ierr
80  character(len=1) :: H2plus_c_converged
81  character(len=1) :: H2_beta_converged
82 
83 ! ---------------------------------------------------------------------------------------
84 ! SECTION: Initial Setup
85 ! ---------------------------------------------------------------------------------------
86  ! Initialise MPI
87  ! Initialize MPI
88  call mpi_init(mpierr)
89 
90  ! Find number of processors and current rank
91  call mpi_comm_size(mpi_comm_world,p,mpierr)
92  call mpi_comm_rank(mpi_comm_world,my_rank,mpierr)
93 
94  ! Read the input parameters.
95  call read_inputs()
96 
97  ! File management from rank 0
98  if (my_rank==0) then
99 
100  ! Notify user of the system
101  write(*, "('System: ', A)") p_system
102  write(*, "('')")
103 
104  ! Make sure the equilibration output directory exists if needed.
105  ! There is no simple way to check directory existence that works
106  ! with both ifort and gfortran, so mkdir is just run with stderr
107  ! suppressed instead.
108  if (run_equil) then
109  call system('mkdir ./equil_out 2>/dev/null')
110  write(*, "('Dice QMC started in Equilibration Mode')")
111  ! Check step length
112  if(run_steps>max_chain_len) then
113  write(*,"('')")
114  write(*, "(A, I7, A)") 'WARNING: The program will be slow at writing equilibration files, &
115  &as the requested number of steps exceeds the recommended &
116  &upper limit of ', max_chain_len, ' steps. Refer to the user manual for details.'
117  write(*,"('')")
118  endif
119  end if
120 
121  ! Same directory check as above for write_chains.
122  if (write_chains) then
123  if (.not. run_equil) then
124  call system('mkdir ./chains_out 2>/dev/null')
125  write(*, "('Dice QMC started in full VQMC mode with additional Markov chain outputting')")
126  ! Check step length
127  if(run_steps>max_chain_len) then
128  write(*,*) 'NOTE: For efficiency, only the first ', max_chain_len,' entries of the Markov chain &
129  &(after burning and thinning) will be written out to the chains_out folder.'
130  endif
131  else
132  error stop "Dice QMC cannot be started with both run_equil and write_chains set &
133  &to TRUE. Check your input params.txt and try again."
134  end if
135  end if
136 
137  ! Initialise the log file if required.
138  if (write_log) then
139 
140  ! Open the logfile
141  open(logid, file=logfile, status='unknown', access='sequential', &
142  form='formatted', iostat=logiostat)
143  if (logiostat .ne. 0) write(0, *) 'Error opening log file ', logfile
144 
145  ! Write program details to logfile
146  write(logid, "('*** Dice QMC log file ***')") ! include date/time?
147  write(logid, "('-------------------------------')")
148  write(logid, *) 'System: ', p_system
149  write(logid, "('')")
150 
151  ! Update log-writer tag for rank 0
152  i_log = .true.
153 
154  endif
155 
156  endif
157 
158  ! If p_system='test', run test module.
159  if (trim(p_system)=='test') then
160  ! Tests need to run on a single processor
161  if (my_rank==0) then
162  call unit_tests()
163  endif
164  ! Exit the program
165  call mpi_finalize(mpierr)
166  call exit()
167  endif
168 
169  ! If not a test job, setup MPI worksharing.
170  ! Equilibrations will be run serially on rank 0.
171  ! For non-equilibration runs:
172  if (.not. run_equil) then
173 
174  ! Split workload for VQMC runs
175  my_steps = run_steps / p ! Integer division
176  if (my_rank < (run_steps-(my_steps*p))) then
177  my_steps = my_steps+1
178  endif
179 
180  ! Pre-calculate frequently used factors.
181  my_size = my_steps/tstep ! Integer division
182  call mpi_reduce(my_size,tot_size,1,mpi_integer,mpi_sum,0,mpi_comm_world,mpierr)
183  if (my_rank==0) then
184  ! Factor representing division by total chain size after burning & thinning
185  per_step = 1.0_real64 / real(tot_size,kind=real64)
186  ! Factor representing division by number of processors
187  per_proc = 1.0_real64 / real(p,kind=real64)
188  ! Factor representing division by square of number of processors
189  per_proc2 = per_proc**2
190  endif
191  ! The factor 'per_proc' is required on all ranks for auto_grad calculations.
192  ! To ensure floating point agreement, broadcaste it from rank 0, rather than re-calculating on every rank.
193  call mpi_bcast(per_proc,1,mpi_double_precision,0,mpi_comm_world,mpierr)
194 
195  endif
196 
197  ! Initialise the random number generator
198  if (run_equil) then
199  ! Run calculations on one processor (rank 0)
200  p = 1
201  if(my_rank==0) then
202  ! Initialise the seed on rank 0
203  call init_ran1()
204  ! Update logfile if required
205  if (write_log) then
206  write(logid, "('Initial seed from /dev/urandom: ', I20)") seed
207  write(logid, "('')")
208  endif
209  ! Ensure auto_params is disabled
210  if (trim(p_system) .eq. 'H2plus' .and. h2plus%auto_params) then
211  h2plus%auto_params = .false.
212  write(*, "('WARNING: H2plus%auto_params should not be .TRUE. when in Equilibration Mode, disabling...')")
213  else if (trim(p_system) .eq. 'H2' .and. h2%auto_params) then
214  h2%auto_params = .false.
215  write(*, "('WARNING: H2plus%auto_params should not be .TRUE. when in Equilibration Mode, disabling...')")
216  end if
217  endif
218  else
219  ! For non-equil VQMC runs,
220  ! Initialise a seed on each MPI rank
221  call init_ran1()
222  ! Update logfile if required
223  if (write_log) then
224  if (my_rank>0) then
225  ! Send the seeds to rank 0
226  call mpi_send(seed,1,mpi_integer8,0,tag_seed,mpi_comm_world,mpierr)
227  else
228  ! Write seed of rank=0 to logfile
229  write(logid, "('Rank ',I2,'. Initial seed from /dev/urandom: ', I20)") my_rank, seed
230  do ip = 1, p-1
231  ! Receive seed value from rank=ip
232  call mpi_recv(ip_seed,1,mpi_integer8,ip,tag_seed,mpi_comm_world,status,mpierr)
233  ! Write seed of rank=ip to logfile
234  write(logid, "('Rank ',I2,'. Initial seed from /dev/urandom: ', I20)") ip, ip_seed
235  end do
236  write(logid, "('')")
237  endif
238  endif
239 
240  endif
241 
242 ! ---------------------------------------------------------------------------------------
243 ! SECTION: VQMC Loops
244 ! ---------------------------------------------------------------------------------------
245 
246  ! Execute problem-dependent routines.
247  select case(trim(p_system))
248 
249  ! -----------------------------------------------------------------------------------
250  ! SUBSECTION: Quantum Harmonic Oscillator
251  ! -----------------------------------------------------------------------------------
252 
253  case('QHO')
254 
255  ! Print out system info.
256  if (my_rank==0) then
257  write(*, "('Doing VQMC runs across ', I3, ' parameter(s).')") int(qho%alpha(3))
258  write(*, "('Steps: ', I10, ', Burn Steps: ', I7, ', Thin Every ', I7, ' Steps.')") &
259  qho%steps, qho%burn_step, qho%thin_step
260  write(*, "('')")
261  endif
262 
263  ! Set up alpha array.
264  p_arr = linspace(qho%alpha(1), qho%alpha(2), qho%alpha(3))
265 
266  if(my_rank==0) then
267  write(*, "('Grid search across linearly spaced parameters enabled.')")
268  write(*, "('')")
269  if (write_log) then
270  write(logid, "('Grid search across linearly spaced parameters enabled.')")
271  write(logid, "('Array of alpha values: ')")
272  do i = 1, int(qho%alpha(3))
273  write(logid, "(F12.8)") p_arr(i)
274  end do
275  write(logid, "('')")
276  end if
277  endif
278 
279  ! Main VQMC callers.
280  if (run_equil) then
281  ! For restarts, input_parser would have already set alpha_start
282  if (my_rank==0) then
283  call equil_qho(p_arr)
284  endif
285  else
286  ! Initialise arrays
287  allocate(energies(int(qho%alpha(3)))); energies(:) = 0.0_real64
288  allocate(variances(int(qho%alpha(3)))); variances(:) = 0.0_real64
289  allocate(accept_rates(int(qho%alpha(3)))); accept_rates(:) = 0.0_real64
290 
291  ! Sort out restarting
292  if (run_restart) then
293  ! Insert optimised parameters, energies, accept_rates, and variances into arrays
294  call read_energies(resfile_etot, energies, variances, accept_rates)
295  endif
296 
297  ! Run the grid evaluation
298  call grid_qho(p_arr, energies, accept_rates, variances)
299  end if
300 
301  ! -----------------------------------------------------------------------------------
302  ! SUBSECTION: Hydrogen Molecular Ion
303  ! -----------------------------------------------------------------------------------
304 
305  case('H2plus')
306 
307  ! Initialise arrays
308  allocate(energies(int(h2plus%bond(3))))
309  allocate(variances(int(h2plus%bond(3))))
310  allocate(accept_rates(int(h2plus%bond(3))))
311  allocate(h2plus_param_c_array(int(h2plus%bond(3))))
312 
313  ! Print out system info.
314  if (my_rank==0) then
315  write(*, "('Doing VQMC runs across ', I7, ' bond lengths.')") int(h2plus%bond(3))
316  write(*, "('Steps: ', I10, ', Burn Steps: ', I7, ', Thin Every ', &
317  & I7, ' Steps.')") h2plus%steps, h2plus%burn_step, h2plus%thin_step
318  write(*, "('')")
319  endif
320 
321  ! Set up bond length array
322  bonds = linspace(h2plus%bond(1), h2plus%bond(2), h2plus%bond(3))
323  if (i_log) then
324  write(logid, "('Array of bond lengths: ')")
325  do i = 1, int(h2plus%bond(3))
326  write(logid, "('Bond length ', I4, ': ', F7.4)") i, bonds(i)
327  end do
328  write(logid, "('')")
329  end if
330 
331  ! Set up parameter array, if doing a grid search
332  if (h2plus%auto_params .eqv. .false.) then
333  p_arr = linspace(h2plus%c_grid(1), h2plus%c_grid(2), h2plus%c_grid(3))
334  endif
335 
336  ! Update user from rank 0
337  if (my_rank==0) then
338  if (h2plus%auto_params .eqv. .false.) then
339  write(*, "('Grid search across linearly spaced parameters enabled.')")
340  write(*, "('')")
341  if (write_log) then
342  write(logid, "('Grid search across linearly spaced parameters enabled.')")
343  write(logid, "('Array of c values: ')")
344  do i = 1, int(h2plus%c_grid(3))
345  write(logid, "(F12.8)") p_arr(i)
346  end do
347  write(logid, "('')")
348  end if
349  else
350  write(*, "('Automatic search across parameters enabled.')")
351  write(*, "('')")
352  if (write_log) then
353  write(logid, "('Automatic search across parameters enabled.')")
354  write(logid, "('')")
355  end if
356  endif
357  endif
358 
359  ! Sort out restarting
360  if (run_restart .and. (.not. run_equil)) then
361  ! Insert optimised parameters, energies, accept_rates, and variances into arrays
362  call read_energies_h2plus(resfile_etot, energies, variances, accept_rates,h2plus_param_c_array)
363  endif
364 
365  ! Set starting guess value
366  if ((.not. run_equil) .and. h2plus%auto_params) then
367  if ( run_restart) then
368  if (.not. res%auto) error stop 'Error: Inconsistency in auto_params setting between params.txt and res_nml.txt.'
369  h2plus_param_c_guess = res%guess
370  else
371  h2plus_param_c_guess = h2plus%c
372  endif
373  endif
374 
375  ! Set up printed output for the loop.
376  if (my_rank==0) then
377  write(*, "('------------------ Begin VQMC Runs ------------------')")
378  write(*, "('Bond Length Converged? Optimised c Energy')")
379  if (write_log) then
380  write(logid, "('------------------ Begin VQMC Runs ------------------')")
381  end if
382  endif
383 
384  ! Loop over bond lengths
385  do b = bond_start, int(h2plus%bond(3))
386 
387  ! Store bond index for restart writing
388  i_bond = b
389 
390  if (i_log) then
391  write(logid, "(' ------------- Bond Length = ', F7.4, ' ------------- ')") &
392  bonds(b)
393  end if
394 
395  ! Main VQMC callers
396  if (.not. run_equil) then
397 
398  if (h2plus%auto_params) then
399 
400  ! Search for best param at this bond length
401  call auto_h_2_plus(h2plus_param_c_guess, bonds(b), h2plus_param_c_opt, &
402  energy, variance, accept_rate, h2plus_c_converged)
403 
404  ! Reset start counter
405  auto_start=1
406 
407  ! Use current optimised value of c as guess for next bond length's c,
408  ! provided calculations converged this loop
409  if (h2plus_c_converged .eq. 'Y') h2plus_param_c_guess = h2plus_param_c_opt
410 
411  else
412 
413  ! Find best param in grid at this bond length
414  call grid_h_2_plus(p_arr, bonds(b), h2plus_param_c_opt, energy, &
415  variance, accept_rate)
416 
417  ! Reset start counter
418  grid_start=1
419 
420  end if
421 
422  if(my_rank==0) then
423 
424  ! Append energy and parameters to arrays, on rank 0
425  h2plus_param_c_array(b) = h2plus_param_c_opt
426  energies(b) = energy
427  variances(b) = variance
428  accept_rates(b) = accept_rate
429 
430  ! Write restart files, if required
431  if (write_restart) then
432  ! Store the arrays
433  call write_h2plus_main(h2plus_param_c_array, bonds, energies, variances, &
434  resfile_etot, ierr, h2plus, accept_rates)
435  ! Write res_nml.txt
436  if (h2plus%auto_params) then
437  call write_res_nml(b+1,1,h2plus%auto_params,h2plus_param_c_guess)
438  else
439  call write_res_nml(b+1,1,h2plus%auto_params)
440  endif
441  endif
442 
443  endif
444 
445  else
446  ! Run equilibrations serially on rank 0
447  if (my_rank==0) then
448  call equil_h_2_plus(p_arr, bonds(b), b)
449  ! Write restart files, if required
450  if (write_restart) then
451  call write_res_nml(b+1,1,h2plus%auto_params)
452  endif
453 
454  endif
455  end if
456 
457  end do
458 
459  if (my_rank==0) then
460 
461  write(*, "('------------------- End VQMC Runs -------------------')")
462  write(*, "('')")
463 
464  ! If not in equilibration mode, write results
465  if (.not. run_equil) then
466  write(*, "(2A)") 'Calculations finished, writing results to NetCDF output file ', &
467  outfile
468  if (write_log) then
469  write(logid, "(2A)") 'Calculations finished, writing results to NetCDF &
470  &output file ', outfile
471  end if
472 
473  call write_h2plus_main(h2plus_param_c_array, bonds, energies, variances, &
474  outfile, ierr, h2plus, accept_rates)
475  end if
476 
477  endif
478 
479 
480  ! -----------------------------------------------------------------------------------
481  ! SUBSECTION: Hydrogen Molecule
482  ! -----------------------------------------------------------------------------------
483 
484  case('H2')
485 
486  ! Initialise arrays
487  allocate(energies(int(h2%bond(3))))
488  allocate(variances(int(h2%bond(3))))
489  allocate(accept_rates(int(h2%bond(3))))
490  allocate(h2_param_a_array(int(h2%bond(3))))
491  allocate(h2_param_beta_array(int(h2%bond(3))))
492 
493  ! Print out system info.
494  if (my_rank==0) then
495  write(*, "('Doing VQMC runs across ', I7, ' bond lengths.')") int(h2%bond(3))
496  write(*, "('Steps: ', I10, ', Burn Steps: ', I7, ', Thin Every ', I7, ' Steps.')") &
497  h2%steps, h2%burn_step, h2%thin_step
498  write(*, "('')")
499  endif
500 
501  ! Set up bond length array
502  bonds = linspace(h2%bond(1), h2%bond(2), h2%bond(3))
503  if (i_log) then
504  write(logid, "('Array of bond lengths: ')")
505  do i = 1, int(h2%bond(3))
506  write(logid, "('Bond length ', I4, ': ', F7.4)") i, bonds(i)
507  end do
508  write(logid, "('')")
509  end if
510 
511  ! Set up parameter array, if doing a grid search
512  if (h2%auto_params .eqv. .false.) then
513  p_arr = linspace(h2%beta_grid(1), h2%beta_grid(2), h2%beta_grid(3))
514  endif
515 
516  ! Update user from rank 0
517  if (my_rank==0) then
518  if (h2%auto_params .eqv. .false.) then
519  write(*, "('Grid search across linearly spaced parameters enabled.')")
520  write(*, "('')")
521  if (write_log) then
522  write(logid, "('Grid search across linearly spaced parameters enabled.')")
523  write(logid, "('Array of beta values: ')")
524  do i = 1, int(h2%beta_grid(3))
525  write(logid, "(F12.8)") p_arr(i)
526  end do
527  write(logid, "('')")
528  end if
529  else
530  write(*, "('Automatic search across parameters enabled.')")
531  write(*, "('')")
532  if (write_log) then
533  write(logid, "('Automatic search across parameters enabled.')")
534  write(logid, "('')")
535  end if
536  endif
537  end if
538 
539  ! Sort out restarting
540  if (run_restart .and. (.not. run_equil)) then
541  ! Insert optimised parameters, energies, accept_rates, and variances into arrays
542  call read_energies_h2(resfile_etot, energies, variances, accept_rates,h2_param_a_array,h2_param_beta_array)
543  endif
544 
545  ! Set starting guess value
546  if ((.not. run_equil) .and. h2%auto_params) then
547  if (run_restart) then
548  if (.not. res%auto) error stop 'Error: Inconsistency in auto_params setting between params.txt and res_nml.txt.'
549  h2_param_beta_guess = res%guess
550  else
551  h2_param_beta_guess = h2%beta
552  endif
553  endif
554 
555  ! Set up printed output for the loop.
556  if (my_rank==0) then
557  write(*, "('------------------ Begin VQMC Runs ------------------')")
558  write(*, "('Bond Length Converged? Optimised Beta Energy')")
559  if (i_log) then
560  write(logid, "('------------------ Begin VQMC Runs ------------------')")
561  end if
562  endif
563 
564  ! Loop over bond lengths
565  do b = bond_start, int(h2%bond(3))
566 
567  ! Store bond index for restart writing
568  i_bond = b
569 
570  if (i_log) then
571  write(logid, "(' ------------- Bond Length = ', F7.4, ' ------------- ')") &
572  bonds(b)
573  end if
574 
575  ! Calculate parameter a at this bond length by Newton-Raphson iteration.
576  h2_param_a = h_2_cusp(bonds(b))
577 
578  ! Main VQMC callers
579  if (.not. run_equil) then
580 
581  ! Either iterate through beta values, or automatically search for best beta.
582  if (h2%auto_params) then ! ie. if auto searching
583 
584  ! Solve value of beta at this bond length.
585  call auto_h_2(h2_param_a, h2_param_beta_guess, bonds(b), h2_param_beta_opt, &
586  energy, variance, accept_rate, h2_beta_converged)
587 
588  ! Use current optimised value of beta as guess for next bond length's beta,
589  ! provided calculations converged this loop
590  if (h2_beta_converged .eq. 'Y') h2_param_beta_guess = h2_param_beta_opt
591 
592  ! Reset start counter
593  auto_start=1
594 
595  else ! ie. if doing a grid search
596 
597  ! Get best value of beta at this bond length
598  call grid_h_2(h2_param_a, p_arr, bonds(b), h2_param_beta_opt, energy, &
599  variance, accept_rate)
600 
601  ! Reset start counter
602  grid_start=1
603 
604  ! Write restart files, if required
605  if (write_restart .and. my_rank==0) then
606  call write_h2_main(h2_param_a_array, h2_param_beta_array, bonds, &
607  energies, variances, resfile_etot, ierr, h2, accept_rates)
608  call write_res_nml(b+1,1,h2%auto_params)
609  endif
610 
611  end if
612 
613  if (my_rank==0) then
614 
615  ! Append energy and parameters to arrays
616  h2_param_a_array(b) = h2_param_a
617  h2_param_beta_array(b) = h2_param_beta_opt
618  energies(b) = energy
619  variances(b) = variance
620  accept_rates(b) = accept_rate
621 
622  ! Write restart files, if required
623  if (write_restart) then
624  ! Store the arrays
625  call write_h2_main(h2_param_a_array, h2_param_beta_array, bonds, &
626  energies, variances, resfile_etot, ierr, h2, accept_rates)
627  if (h2%auto_params) then
628  call write_res_nml(b+1,1,h2%auto_params,h2_param_beta_guess)
629  else
630  call write_res_nml(b+1,1,h2%auto_params)
631  endif
632  endif
633 
634  endif
635 
636  else
637 
638  ! Run equilibrations
639  if (my_rank==0) then
640  call equil_h_2(h2_param_a, p_arr, bonds(b), b)
641  ! Write restart files, if required
642  if (write_restart) then
643  call write_res_nml(b+1,1,h2%auto_params)
644  endif
645  endif
646 
647  end if
648 
649  end do
650 
651  if (my_rank==0) then
652 
653  write(*, "('------------------- End VQMC Runs -------------------')")
654  write(*, "('')")
655 
656  ! If not in equilibration mode, write results
657  if (.not. run_equil) then
658 
659  write(*, "(2A)") 'Calculations finished, writing results to NetCDF output file ', &
660  outfile
661 
662  if (write_log) then
663  write(logid, "(2A)") 'Calculations finished, writing results to NetCDF &
664  &output file ', outfile
665  end if
666 
667  call write_h2_main(h2_param_a_array, h2_param_beta_array, bonds, &
668  energies, variances, outfile, ierr, h2, accept_rates)
669 
670  end if
671 
672  endif
673 
674  ! Default case, for if none of the cases match. Should be invalidated by checks in input_parser.
675  case default
676  error stop 'Unknown case. How did you get here?'
677  end select
678 
679  ! Close logfile from rank 0
680  if (i_log) close(logid)
681 
682  ! Remove restart files if complete
683  call system('rm res_etot.nc 2>/dev/null')
684  call system('rm res_nml.txt 2>/dev/null')
685  call system('rm res_epar.nc 2>/dev/null')
686 
687  ! End MPI parallelism
688  call mpi_finalize(mpierr)
689 
690 end subroutine driver
691 
692 end program main
main
program main
Main program, contains the driver subroutine that runs all other code.
Definition: driver.f90:14
driver
subroutine driver()
Subroutine for running the main code (needed for Doxygen compatability).
Definition: driver.f90:54
input_parser
Contains routines that read the input file params.txt.
Definition: input_parser.f90:16
input_parser::read_inputs
subroutine read_inputs()
Main input parser function that reads the input file params.txt and calls the required namelist parse...
Definition: input_parser.f90:165
restart_fns
Contains routines to read and write restart files.
Definition: restart_fns.f90:15
restart_fns::write_res_nml
subroutine write_res_nml(i_qoi, i_par, auto, guess)
Writes the restart file res_nml.txt, using the current states of simulation variables.
Definition: restart_fns.f90:48
restart_fns::read_energies_h2plus
subroutine read_energies_h2plus(filename, energies, variances, accept_rates, c_array)
Reads restart files with energies, variances and acceptence arrays. Also reads one extra array of opt...
Definition: restart_fns.f90:239
restart_fns::read_energies
subroutine read_energies(filename, energies, variances, accept_rates)
Reads restart files with energies, variances and acceptence arrays.
Definition: restart_fns.f90:177
restart_fns::read_energies_h2
subroutine read_energies_h2(filename, energies, variances, accept_rates, a_array, b_array)
Reads restart files with energies, variances and acceptence arrays. Also reads two extra arrays of op...
Definition: restart_fns.f90:310
shared_data
Contains derived types, and global variables to store input values of simulation parameters.
Definition: shared_data.f90:49
shared_data::seed
integer(int64) seed
Seed for the random number generator.
Definition: shared_data.f90:152
shared_data::run_steps
integer(int32) run_steps
Number of MMC steps after burning (for MPI worksharing)
Definition: shared_data.f90:153
shared_data::h2
type(h2_type) h2
To store inputs for H2 molecule system.
Definition: shared_data.f90:140
shared_data::h2plus
type(h2plus_type) h2plus
To store inputs for H2 ion system.
Definition: shared_data.f90:139
shared_data::outfile
character(len= *), parameter outfile
Default NetCDF output file name.
Definition: shared_data.f90:134
shared_data::run_restart
logical run_restart
Tag for restart runs.
Definition: shared_data.f90:129
shared_data::res
type(restart_type) res
Instance of restart_type for reading restart files.
Definition: shared_data.f90:146
shared_data::tstep
integer(int32) tstep
Interval of steps for thinning in non-equil runs (for MPI worksharing)
Definition: shared_data.f90:154
shared_data::auto_start
integer(int32) auto_start
Start value for automatic search of parameter c, can be >1 for restarts.
Definition: shared_data.f90:158
shared_data::write_chains
logical write_chains
Tag to write final MMC chains.
Definition: shared_data.f90:132
shared_data::ip
integer(int32) ip
Processor rank for MPI communications.
Definition: shared_data.f90:169
shared_data::i_bond
integer(int32) i_bond
Index of bondlength grid value, for restart writing.
Definition: shared_data.f90:159
shared_data::per_step
real(real64) per_step
Reused multiplicative factor, representing division by total number of MMC steps (tot_size)
Definition: shared_data.f90:175
shared_data::per_proc2
real(real64) per_proc2
Reused multiplicative factor, representing division by square of number of processors (p^2)
Definition: shared_data.f90:177
shared_data::qho
type(qho_type) qho
To store inputs for QHO system.
Definition: shared_data.f90:138
shared_data::per_proc
real(real64) per_proc
Reused multiplicative factor, representing division by number of processors (p)
Definition: shared_data.f90:176
shared_data::i_log
logical i_log
Indicates whether the processor needs to write the logfile.
Definition: shared_data.f90:161
shared_data::tag_seed
integer(int32) tag_seed
Unique tag for communicating rank seeds.
Definition: shared_data.f90:171
shared_data::resfile_etot
character(len= *), parameter resfile_etot
Name of restart file storing total energies.
Definition: shared_data.f90:148
shared_data::run_equil
logical run_equil
Tag for equilibration runs.
Definition: shared_data.f90:128
shared_data::status
integer, dimension(mpi_status_size) status
Status of MPI Recv commands.
Definition: shared_data.f90:174
shared_data::grid_start
integer(int32) grid_start
Start value for grid search of parameter c, can be >1 for restarts.
Definition: shared_data.f90:157
shared_data::my_size
integer(int32) my_size
Length of chain on a processor, after burning & thinning.
Definition: shared_data.f90:172
shared_data::p
integer(int32) p
Number of processors.
Definition: shared_data.f90:166
shared_data::bond_start
integer(int32) bond_start
Start value for bondlength grid, can be >1 for restarts.
Definition: shared_data.f90:156
shared_data::p_system
character(len=10) p_system
System of interest: QHO, H2plus, H2.
Definition: shared_data.f90:127
shared_data::write_log
logical write_log
Tag to write logfile.
Definition: shared_data.f90:133
shared_data::logfile
character(len= *), parameter logfile
Default log name (if write_log is True)
Definition: shared_data.f90:135
shared_data::max_chain_len
integer max_chain_len
Maximum length of energy chains for writing to netCDF files.
Definition: shared_data.f90:162
shared_data::tot_size
integer(int32) tot_size
Total length of chains across processors, after burning and thinning.
Definition: shared_data.f90:173
shared_data::write_restart
logical write_restart
Tag for writing restart files.
Definition: shared_data.f90:130
shared_data::my_steps
integer(int32) my_steps
Number of steps per processor.
Definition: shared_data.f90:168
shared_data::ip_seed
integer(int64) ip_seed
Processor seed for MPI communications.
Definition: shared_data.f90:170
shared_data::my_rank
integer(int32) my_rank
Rank of current processor.
Definition: shared_data.f90:167
shared_data::mpierr
integer(int32) mpierr
MPI error flag.
Definition: shared_data.f90:165
shared_data::logid
integer logid
To store logfile unit.
Definition: shared_data.f90:160
solvers
Contains local energy and transition probability solvers for available problems.
Definition: quantum_solvers.f90:16
solvers::h_2_cusp
real(real64) function h_2_cusp(s)
Calculates a value for a given the bond length s by Newton-Raphson iteration.
Definition: quantum_solvers.f90:395
testing
Contains unit tests for all major functions in the program.
Definition: testing.f90:15
testing::unit_tests
subroutine unit_tests()
Runs all of the unit tests, reporting successes and failures.
Definition: testing.f90:794
vqmc
Contains routines for handline random number generation, the main VQMC algorithms and routines for de...
Definition: vqmc.f90:30
vqmc::grid_h_2
subroutine grid_h_2(a, beta_array, bond_length, beta_out, energy_out, variance_out, accept_out)
Performs a grid search to find the optimal value for the variational parameter in the problem.
Definition: vqmc.f90:2155
vqmc::auto_h_2_plus
subroutine auto_h_2_plus(c_in, bond_length, c_out, energy_out, variance_out, accept_out, converged)
Automatically finds the best variational parameter for the problem at given bond length.
Definition: vqmc.f90:1190
vqmc::equil_h_2_plus
subroutine equil_h_2_plus(c_array, bond_length, bonds_pos)
Performs VQMC on an array of values for the variational parameter , outputting energy traces for each...
Definition: vqmc.f90:1506
vqmc::equil_qho
subroutine equil_qho(alpha_array)
Performs VQMC on an array of values for the variational parameter , outputting energy traces for each...
Definition: vqmc.f90:743
vqmc::linspace
real(real64) function, dimension(:), allocatable linspace(start, end, num)
Adaptation of numpy's linspace function.
Definition: vqmc.f90:183
vqmc::grid_h_2_plus
subroutine grid_h_2_plus(c_array, bond_length, c_out, energy_out, variance_out, accept_out)
Performs a grid search to find the optimal value for the variational parameter in the problem.
Definition: vqmc.f90:1361
vqmc::auto_h_2
subroutine auto_h_2(a, beta_in, bond_length, beta_out, energy_out, variance_out, accept_out, converged)
Automatically finds the best variational parameter for the problem at a given bond length.
Definition: vqmc.f90:1982
vqmc::grid_qho
subroutine grid_qho(alpha_array, energies, acc_rates, variances)
Performs VQMC on an array of values for the variational parameter .
Definition: vqmc.f90:580
vqmc::equil_h_2
subroutine equil_h_2(a, beta_array, bond_length, bonds_pos)
Performs VQMC on an array of values for the variational parameter , outputting energy traces for each...
Definition: vqmc.f90:2300
vqmc::init_ran1
subroutine init_ran1()
Initialiser for the random number generator.
Definition: vqmc.f90:69
write_netcdf
Contains NetCDF output functions for different problem main & equilibration runs.
Definition: netcdf_out.f90:16
write_netcdf::write_h2plus_main
subroutine write_h2plus_main(param_array, bondlength_array, energy_array, uncertainty_array, filename, ierr, H2plus, accept_rate)
Output results of H2plus MCMC run.
Definition: netcdf_out.f90:439
write_netcdf::write_h2_main
subroutine write_h2_main(param_a_array, param_beta_array, bondlength_array, energy_array, uncertainty_array, filename, ierr, H2, accept_rate)
Output results of H2 MCMC run.
Definition: netcdf_out.f90:544