Re: [QE-users] Random seed for MD calculations

[Paolo Giannozzi]

Thank you for noticing again a problem. It will be fixed in next release:
  https://gitlab.com/QEF/q-e/-/merge_requests/2137
Paolo

On 8/14/23 19:24, Han Hsu (徐翰) wrote:


You don't often get email from han...@ncu.edu.tw. Learn why this is 
important 



Dear Professor Giannozzi, QE developers, and QE users:

Previously, I reported that while the patch file provided by Prof. 
Giannozzi works (namely, not all restart MD calculations got the same 
sequence of the random numbers), a particular sequence of random numbers 
is generated much more frequently than the others: Nearly half of my 
restart calculations got that sequence of random numbers. Now I think I 
have found out the reason why. By modifying Modules/random_numbers.f90, 
in particular, its subroutine "set_random_seed" (lines 74–87), this 
problem can be resolved, as described below.


In the subroutine "set_random_seed" , the parameter "iseed" is generated 
based on the system time (lines 81–84), as shown below. By calling 
"date_and_time", the system time is extracted in terms of eight numbers 
(year, month, ..., milliseonds). Four of these numbers are adopted to 
generate "iseed", including the time difference with respect to UTC (in 
minutes) [itime(4)], and minutes, seconds, and milliseconds 
[itime(6–8)]. Among them, I find the inclusion of itime(4) 
incomprehensible, as this term is a constant based on the location of 
the computer. Remarkably, in the UK, itime(4)= 0, meaning this term has 
ZERO effects on "iseed", so why should it be included anyway? By 
contrast, in my country Taiwan (GMT+8), itime(4)= 480, which results in 
a much larger "iseed" than in the UK, by up to hundreds of times...


CALL date_and_time ( values = itime )
! itime contains: year, month, day, time difference in minutes, hours,
!                 minutes, seconds and milliseconds.
iseed = ( itime(8) + itime(6) ) * ( itime(7) + itime(4) )
irand = randy ( iseed )

Given the above analysis, I modified this subroutine by replacing 
itime(4) with itime(5) (namely, hours), making "iseed" fully 
time-dependent. With this modification, all my test calculations 
(from_scratch and all the restart rounds) got different sequences of 
random numbers. If this modification is correct, perhaps it can be 
included in the next release of Quantum Espresso (7.3)?


Thank you for your suggestions and advice.

Best regards,

Han Hsu

On Sat, Aug 12, 2023 at 1:25 AM Han Hsu (徐翰) > wrote:


Dear Professor Giannozzi

Thank you very much for the patch file. It works!! Now, the random
numbers generated in all rounds of calculations (including
from_scratch and all the restart calculations) are not always
repetitive.

After more tests, however, I noticed that some "random numbers" are
generated much much more frequently than the others. Not
surprisingly, this occurs in both 'andersen' and 'svr' thermostats.
What surprises me is that this occurs on different machines at
different institutes as well. I therefore believe this should be
reported, as further described below.

Again, all the test runs were performed using a 8-atom cell of fcc
Si, with dt = 20, tempw = 600. Each run lasts for 12000 steps,
divided into a from_scratch (steps 1–1000) and 11 restart
calculations (steps 1001–2000, 2001–3000, ..., 11001–12000). Both
'andersen' and 'svr' are tested (nraise = 100). A few representative
results obtained on two different machines are shown below (JobID on
these two machines are very different). In order to print out the
random numbers, additional WRITE statements are added to the patched
dynamics_module.f90 file.


(1) svr (random numbers R1, R2, ..., R_Nf are generated in each step)

In run_4323, same sequence of random numbers were generated in
rounds 2, 4, 6, 8, 10, and 12. In run_6632321 (on a different
machine), same sequence of random numbers were generated in rounds
3, 5, 6, 7, 8, and 12. Furthermore, all these 12 rounds on two
different machines got the same sequence of random numbers...

===
run_4323
===
In pwscf.md.1.out, ...
      Entering Dynamics:    iteration =     1
      rr**2 = R1**2 =  0.5153
      sum_of_gaussians2 = R2**2+...+R_Nf**2 =     26.0457
      Entering Dynamics:    iteration =     2
      rr**2 = R1**2 =  0.0411
      sum_of_gaussians2 = R2**2+...+R_Nf**2 =     14.9090
      Entering Dynamics:    iteration =     3
      rr**2 = R1**2 =  0.2940
      sum_of_gaussians2 = R2**2+...+R_Nf**2 =     20.9707
      Entering Dynamics:    iteration =     4
      rr**2 = R1**2 =  1.0311
      sum_of_gaussians2 = R2**2+...+R_Nf**2 =     23.6863
      Entering Dynamics:    iteration =     5
      rr**2 = R1**2 =  0.0971
      sum_of_gaussians2 = R2**2+...+R_Nf**2 =      8.4677
      Entering 

Re: [QE-users] Random seed for MD calculations

[徐翰]

Dear Professor Giannozzi

Thank you very much for the patch file. It works!! Now, the random numbers
generated in all rounds of calculations (including from_scratch and all the
restart calculations) are not always repetitive.

After more tests, however, I noticed that some "random numbers" are
generated much much more frequently than the others. Not surprisingly, this
occurs in both 'andersen' and 'svr' thermostats. What surprises me is that
this occurs on different machines at different institutes as well. I
therefore believe this should be reported, as further described below.

Again, all the test runs were performed using a 8-atom cell of fcc Si, with
dt = 20, tempw = 600. Each run lasts for 12000 steps, divided into a
from_scratch (steps 1–1000) and 11 restart calculations (steps 1001–2000,
2001–3000, ..., 11001–12000). Both 'andersen' and 'svr' are tested (nraise
= 100). A few representative results obtained on two different machines are
shown below (JobID on these two machines are very different). In order to
print out the random numbers, additional WRITE statements are added to the
patched dynamics_module.f90 file.


(1) svr (random numbers R1, R2, ..., R_Nf are generated in each step)

In run_4323, same sequence of random numbers were generated in rounds 2, 4,
6, 8, 10, and 12. In run_6632321 (on a different machine), same sequence of
random numbers were generated in rounds 3, 5, 6, 7, 8, and 12. Furthermore,
all these 12 rounds on two different machines got the same sequence of
random numbers...

===
run_4323
===
In pwscf.md.1.out, ...
 Entering Dynamics:iteration = 1
 rr**2 = R1**2 =  0.5153
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 26.0457
 Entering Dynamics:iteration = 2
 rr**2 = R1**2 =  0.0411
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 14.9090
 Entering Dynamics:iteration = 3
 rr**2 = R1**2 =  0.2940
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 20.9707
 Entering Dynamics:iteration = 4
 rr**2 = R1**2 =  1.0311
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 23.6863
 Entering Dynamics:iteration = 5
 rr**2 = R1**2 =  0.0971
 sum_of_gaussians2 = R2**2+...+R_Nf**2 =  8.4677
 Entering Dynamics:iteration = 6
 rr**2 = R1**2 =  0.1053
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 16.8005
 Entering Dynamics:iteration = 7
 rr**2 = R1**2 =  1.3552
In pwscf.md.2.out, ...
 Entering Dynamics:iteration =  1001
 rr**2 = R1**2 =  1.4715
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 24.4598
 Entering Dynamics:iteration =  1002
 rr**2 = R1**2 =  1.1178
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 19.6385
 Entering Dynamics:iteration =  1003
 rr**2 = R1**2 =  1.3458
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 15.8245
 Entering Dynamics:iteration =  1004
 rr**2 = R1**2 =  0.6379
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 18.7133
 Entering Dynamics:iteration =  1005
 rr**2 = R1**2 =  0.0146
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 11.9063
 Entering Dynamics:iteration =  1006
 rr**2 = R1**2 =  2.2273
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 23.3447
 Entering Dynamics:iteration =  1007
 rr**2 = R1**2 =  0.0441
In pwscf.md.3.out, ...
 Entering Dynamics:iteration =  2001
 rr**2 = R1**2 =  0.0243
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 15.0970
 Entering Dynamics:iteration =  2002
 rr**2 = R1**2 =  0.8330
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 16.2563
 Entering Dynamics:iteration =  2003
 rr**2 = R1**2 =  1.2013
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 23.2091
 Entering Dynamics:iteration =  2004
 rr**2 = R1**2 =  2.0586
 sum_of_gaussians2 = R2**2+...+R_Nf**2 =  8.7939
 Entering Dynamics:iteration =  2005
 rr**2 = R1**2 =  0.1827
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 15.2563
 Entering Dynamics:iteration =  2006
 rr**2 = R1**2 =  0.0228
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 12.6388
 Entering Dynamics:iteration =  2007
 rr**2 = R1**2 =  2.0902
In pwscf.md.4.out, ...
 Entering Dynamics:iteration =  3001
 rr**2 = R1**2 =  1.4715
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 24.4598
 Entering Dynamics:iteration =  3002
 rr**2 = R1**2 =  1.1178
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 19.6385
 Entering Dynamics:iteration =  3003
 rr**2 = R1**2 =  1.3458
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 15.8245
 Entering Dynamics:iteration =  3004
 rr**2 = R1**2 =  0.6379
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 18.7133
 Entering Dynamics:iteration =  3005
 rr**2 = R1**2 =  0.0146
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 11.9063
 Entering Dynamics:iteration =  3006
 rr**2 = R1**2 =  2.2273
 sum_of_gaussians2 = 

Re: [QE-users] Random seed for MD calculations

[Paolo Giannozzi]

Thank you for reporting this. Could you please try the attached patch 
and report if it fixes your problem?


Paolo

On 8/10/23 12:36, Han Hsu (徐翰) wrote:


You don't often get email from han...@ncu.edu.tw. Learn why this is 
important 



Dear QE developers and users,

I've been performing Born-Oppenheimer molecular dynamics (BOMD) 
calculations for hcp Fe at high-PT conditions using Quantum Espresso 
(pw.x, calculation= 'md'). For the thermostats, I am testing stochastic 
velocity rescaling (svr) and Andersen, as they are the only ones (in QE) 
that generate canonical ensembles. For these two thermostats, random 
numbers are generated in each time step. I find, however, that in 
restart calculations (restart_mode= 'restart'), the random numbers are 
not generated properly. In each restart calculation, the same sequence 
of random numbers are generated, as shown in detail in the next few 
paragraphs. Consequently, the calculation results (e.g pressure, 
temperature) exhibit a periodic behavior, with a period same as the 
length of each restart calculation. My question is: Can we avoid this 
problem by properly setting the random seed in restart calculations? If 
so, how to do that? (For restart_mode= 'from_scratch', such a problem 
does not occur, as also shown in the next few paragraphs).


Note: To print out the random numbers generated in each time step, a few 
more lines (write statements) have to be added to the file 
PW/src/dynamics_module.f90. Also, to obtain lots of data within a short 
time, I use fcc Si (8-atom cells) to run extensive tests. For these 
tests, each run consists of 12000 time steps divided into 12 rounds of 
calculations: A from_scratch round (steps 1–1000) followed by 11 restart 
rounds (steps 1001–2000, 2001–3000, ..., 11001–12000).


For the 'svr' thermostat, random numbers (R1, R2,..., R_Nf) from 
gaussian distribution are generated in teach time step. In the 
from_scratch calculation (steps 1–1000), the random numbers generated in 
steps 1–4 are are as below:


      Entering Dynamics:    iteration =     1
      rr**2 = R1**2 =  2.2273
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 23.3447
      Entering Dynamics:    iteration =     2
      rr**2 = R1**2 =  0.0441
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 17.0921
      Entering Dynamics:    iteration =     3
      rr**2 = R1**2 =  0.6032
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 16.4800
      Entering Dynamics:    iteration =     4
      rr**2 = R1**2 =  1.2989
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 13.2204

Next, in the restart calculation for steps 1001–2000, the random numbers 
generated in steps 1001–1004 are are as below:


      Entering Dynamics:    iteration =  1001
      rr**2 = R1**2 =  0.0008
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 18.0810
      Entering Dynamics:    iteration =  1002
      rr**2 = R1**2 =  0.8306
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 19.9568
      Entering Dynamics:    iteration =  1003
      rr**2 = R1**2 =  1.0171
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 14.5561
      Entering Dynamics:    iteration =  1004
      rr**2 = R1**2 =  0.0605
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 28.3579

Next, in the restart calculation for steps 2001–3000, the random numbers 
generated in steps 2001–2004 are are as below:


      Entering Dynamics:    iteration =  2001
      rr**2 = R1**2 =  0.0008
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 18.0810
      Entering Dynamics:    iteration =  2002
      rr**2 = R1**2 =  0.8306
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 19.9568
      Entering Dynamics:    iteration =  2003
      rr**2 = R1**2 =  1.0171
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 14.5561
      Entering Dynamics:    iteration =  2004
      rr**2 = R1**2 =  0.0605
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 28.3579

Next, in the restart calculation for steps 3001–4000, the random numbers 
generated in steps 3001–3004 are are as below:


      Entering Dynamics:    iteration =  3001
      rr**2 = R1**2 =  0.0008
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 18.0810
      Entering Dynamics:    iteration =  3002
      rr**2 = R1**2 =  0.8306
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 19.9568
      Entering Dynamics:    iteration =  3003
      rr**2 = R1**2 =  1.0171
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 14.5561
      Entering Dynamics:    iteration =  3004
      rr**2 = R1**2 =  0.0605
      sum_of_gaussians2 = R2**2+...+R_Nf**2 = 28.3579

Clearly, the same sequence of random numbers are generated in all 
restart rounds (steps 1001–2000, 2001–3000, 3001–4000, ..., etc.), while 
the from_scratch round (steps 1–1000) is an exception. I am thus 
wondering how to properly generate/seed the random numbers in restart 
calculations?


A few more remarks:
(1) I have done this test on different facilities at different 
institutes. The same pattern persists.
(2) For Andersen thermostat, the same 

[QE-users] Random seed for MD calculations

[徐翰]

Dear QE developers and users,

I've been performing Born-Oppenheimer molecular dynamics (BOMD)
calculations for hcp Fe at high-PT conditions using Quantum Espresso (pw.x,
calculation= 'md'). For the thermostats, I am testing stochastic velocity
rescaling (svr) and Andersen, as they are the only ones (in QE) that
generate canonical ensembles. For these two thermostats, random numbers are
generated in each time step. I find, however, that in restart calculations
(restart_mode= 'restart'), the random numbers are not generated properly.
In each restart calculation, the same sequence of random numbers are
generated, as shown in detail in the next few paragraphs. Consequently, the
calculation results (e.g pressure, temperature) exhibit a periodic
behavior, with a period same as the length of each restart calculation. My
question is: Can we avoid this problem by properly setting the random seed
in restart calculations? If so, how to do that? (For restart_mode=
'from_scratch', such a problem does not occur, as also shown in the next
few paragraphs).

Note: To print out the random numbers generated in each time step, a few
more lines (write statements) have to be added to the file
PW/src/dynamics_module.f90. Also, to obtain lots of data within a short
time, I use fcc Si (8-atom cells) to run extensive tests. For these tests,
each run consists of 12000 time steps divided into 12 rounds of
calculations: A from_scratch round (steps 1–1000) followed by 11 restart
rounds (steps 1001–2000, 2001–3000, ..., 11001–12000).

For the 'svr' thermostat, random numbers (R1, R2,..., R_Nf) from gaussian
distribution are generated in teach time step. In the from_scratch
calculation (steps 1–1000), the random numbers generated in steps 1–4 are
are as below:

 Entering Dynamics:iteration = 1
 rr**2 = R1**2 =  2.2273
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 23.3447
 Entering Dynamics:iteration = 2
 rr**2 = R1**2 =  0.0441
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 17.0921
 Entering Dynamics:iteration = 3
 rr**2 = R1**2 =  0.6032
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 16.4800
 Entering Dynamics:iteration = 4
 rr**2 = R1**2 =  1.2989
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 13.2204

Next, in the restart calculation for steps 1001–2000, the random numbers
generated in steps 1001–1004 are are as below:

 Entering Dynamics:iteration =  1001
 rr**2 = R1**2 =  0.0008
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 18.0810
 Entering Dynamics:iteration =  1002
 rr**2 = R1**2 =  0.8306
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 19.9568
 Entering Dynamics:iteration =  1003
 rr**2 = R1**2 =  1.0171
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 14.5561
 Entering Dynamics:iteration =  1004
 rr**2 = R1**2 =  0.0605
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 28.3579

Next, in the restart calculation for steps 2001–3000, the random numbers
generated in steps 2001–2004 are are as below:

 Entering Dynamics:iteration =  2001
 rr**2 = R1**2 =  0.0008
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 18.0810
 Entering Dynamics:iteration =  2002
 rr**2 = R1**2 =  0.8306
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 19.9568
 Entering Dynamics:iteration =  2003
 rr**2 = R1**2 =  1.0171
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 14.5561
 Entering Dynamics:iteration =  2004
 rr**2 = R1**2 =  0.0605
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 28.3579

Next, in the restart calculation for steps 3001–4000, the random numbers
generated in steps 3001–3004 are are as below:

 Entering Dynamics:iteration =  3001
 rr**2 = R1**2 =  0.0008
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 18.0810
 Entering Dynamics:iteration =  3002
 rr**2 = R1**2 =  0.8306
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 19.9568
 Entering Dynamics:iteration =  3003
 rr**2 = R1**2 =  1.0171
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 14.5561
 Entering Dynamics:iteration =  3004
 rr**2 = R1**2 =  0.0605
 sum_of_gaussians2 = R2**2+...+R_Nf**2 = 28.3579

Clearly, the same sequence of random numbers are generated in all restart
rounds (steps 1001–2000, 2001–3000, 3001–4000, ..., etc.), while the
from_scratch round (steps 1–1000) is an exception. I am thus wondering how
to properly generate/seed the random numbers in restart calculations?

A few more remarks:
(1) I have done this test on different facilities at different institutes.
The same pattern persists.
(2) For Andersen thermostat, the same pattern persists: In the restart
calculations, at the same respective time steps (e.g. steps 1001, 2001,
3001, ...), the same atoms are "randomly" chosen to interact with the heat
bath.

Finally, my understanding is that in Quantum Espresso, generation and seed
of random numbers are controlled by Modules/random_numbers.f90 and a
subroutine set_random_seed(). To properly set the random