The gfortran compiler (and several others) implemented RANDOM_SEED(), with no arguments, as a way of seeding the random number generator non-repeatably. This seemed the best approach for our Fortran example codes, when it became available. However, this was not part of the Fortran standard, and so is implementation-dependent. In Fortran 2018 the RANDOM_INIT subroutine was introduced, with a REPEATABLE=.FALSE. argument, and this has been in gfortran since v9.1. Plenty of time has elapsed, so we might replace RANDOM_SEED() by this soon.

I'm using ubuntu 16.04 with fftw, blas, lapack installed. However, I'm failing to compile the program even after adding the correct path to the libraries in SConstruct file.

Please let me know what I might be doing wrong. I did the following changes to the SConstruct file to provide proper path to the fftw and lapack libraries.

FFTW_F90FLAGS='-fdefault-real-8 -fall-intrinsics -std=f2008 -Wall -I/usr/include'
LAPACK_LINKFLAGS='-L/usr/lib/ -llapack'
FFTW_LINKFLAGS='-L/usr/lib/x86-64-linux-gnu/ -lfftw3'

errors.txt

I'm attaching the error.txt file that has the last few lines that were printed on screen including the compilation error, when it failed to build.

This would not be too difficult, and would probably be useful, although it would be challenging to run the code long enough to get results comparable with the Fortran version.

In the row 123, column 31,of the code file md_nvt_lj.f90, the number of dim is 2. But I think it should be equal to 1. I am not sure.
The detailed content is as the folling picture describes:

Hi Allen,
It is more conceptual question. In MD we add lrc beyond cutoff radius for configurational energy why dont we consider or add any thing for forces coming from lj like interactions beyond cutoff radius? I can understand taking derivative of lj potential energy would increase one extra power in denominator for each term which will decrease the forces between particles with increasing separation. Please explain or give me some reference for this.

Is there a bug in md_nve_hs.f90 on lines 158-165? I understand that this strategy stores a collision partner in the partner array using the lower of the two collision indices. However, the collision between i and j can open up collisions of partner(k) with any other hard sphere, correct? This should require calling update with 'gt' and 'lt' in the loop making it unnecessary to call update ( i, box, lt ) and update ( j, box, lt ) outside the loop.

! Update collision lists DO k = 1, n IF ( k==i .OR. k==j .OR. partner(k) == i .OR. partner(k) == j ) THEN CALL update ( k, box, gt ) ! Search for partners >k END IF END DO CALL update ( i, box, lt ) ! Search for partners <i CALL update ( j, box, lt ) ! Search for partners <j

In the Fortran example codes, a handful of FORALL statements appear. This statement was declared obsolescent in the Fortran 2018 standard, so we should soon replace those instances.

The files link_list_module.f90, mc_lj_ll_module.f90, md_lj_ll_module.f90, and md_lj_llle_module.f90 were all written on the assumption of fixed box dimensions, for simplicity. In particular the number of cells sc is set, and the head array dimensioned, just once at the start of the program. Nonetheless, the SConstruct file includes a build_mc_npt_lj_ll/mc_npt_lj program using linked lists, and the GUIDE includes test results obtained for this code. This program will give incorrect results if the cell size becomes smaller than the cutoff, and will also become needlessly inefficient if the box size grows significantly. If it is not possible to modify the MC codes in a clean and simple way, we should remove that build, and those results. On the other hand, it is fairly easy to modify the MD codes along those lines, so we should do that, and include a build of md_npt_lj with linked lists.

The same problems do not apply to the Python example based on md_lj_ll_module.py, but we should include a test which terminates execution if sc<3 at any point.

It should not be difficult to add a Python version of the Gibbs ensemble code, and may be helpful.

The file grint_data.zip is larger than all the others combined, and only of interest to people who want to test the grint example. Consider moving it to a separate data repository.

In a few of the Python examples, the NumPy type np.int was used instead of np.int_ (or a more specific type such as np.int64). This is most likely a typo. np.int was actually an alias for the builtin int, and was deprecated in NumPy 1.20. Thus, the error has been inconsequential until recently, but the deprecations expired in NumPy 1.24, so it needs fixing promptly. The affected files are grint.py, initialize.py, mc_chain_nvt_sw.py, and mc_chain_wl_sw.py.

This example would probably benefit from adding a Python version. At the same time we could do some tidying up: calculating atomic positions outside the potential routines might be an improvement.

Dear Dr. Allen

I get errors when I compile your source code. I tried to comment such lines in SConstruct to localize problematic source files.

Here is the list of errors:
1.

$ scons
scons: Reading SConscript files ...
scons: done reading SConscript files.
scons: Building targets ...
scons: building associated VariantDir targets: build_md_nve_lj_omp build_smc_nvt_lj build_md_chain_nve_lj build_qmc_pi_sho build_grint build_test_pot_qq build_md_chain_mts_lj build_test_pot_gb build_md_nve_lj build_adjust build_test_pot_twist build_mc_npt_lj_ll build_mc_zvt_lj_ll build_bd_nvt_lj build_pair_distribution build_test_pot_at build_mc_chain_wl_sw build_test_pot_dq build_md_nve_hs build_mc_chain_nvt_sw build_md_nvt_lj build_mc_nvt_lj_re build_mc_nvt_poly_lj build_md_nvt_lj_llle build_md_nve_lj_ll build_test_pot_bend build_md_nvt_lj_le build_test_pot_dd build_cluster build_mc_zvt_lj build_error_calc build_mc_nvt_lj_ll build_initialize build_eos_hs build_diffusion build_mc_nvt_sc build_t_tensor build_mc_npt_hs build_mc_npt_lj build_mc_nvt_lj build_mc_nvt_hs build_eos_lj build_hit_and_miss build_qmc_walk_sho build_sample_mean build_mc_npt_sc build_md_nve_lj_vl build_md_lj_mts build_mc_chain_nvt_cbmc_lj build_md_npt_lj build_dpd build_mc_gibbs_lj build_corfun build_mesh build_fft3dwrap build_wl_hist build_qmc_pi_lj build_diffusion_test
gfortran -o build_diffusion_test/diffusion_test.o -c -fdefault-real-8 -fall-intrinsics -std=f2008 -Wall -Ibuild_diffusion_test -Jbuild_diffusion_test build_diffusion_test/diffusion_test.f90
build_diffusion_test/diffusion_test.f90:134:73:

     CALL random_normals ( 0.0, SQRT(temperature), zeta ) ! Random momenta
                                                                         1
Error: There is no specific subroutine for the generic ‘random_normals’ at (1)
build_diffusion_test/diffusion_test.f90:84:51:

   CALL random_normals ( 0.0, SQRT(temperature), v )
                                                   1
Error: There is no specific subroutine for the generic ‘random_normals’ at (1)
scons: *** [build_diffusion_test/diffusion_test.o] Error 1
scons: building terminated because of errors. 

$ scons
scons: Reading SConscript files ...
scons: done reading SConscript files.
scons: Building targets ...
scons: building associated VariantDir targets: build_md_nve_lj_omp build_smc_nvt_lj build_md_chain_nve_lj build_qmc_pi_sho build_grint build_test_pot_qq build_md_chain_mts_lj build_test_pot_gb build_md_nve_lj build_adjust build_test_pot_twist build_mc_npt_lj_ll build_mc_zvt_lj_ll build_bd_nvt_lj build_pair_distribution build_test_pot_at build_mc_chain_wl_sw build_test_pot_dq build_md_nve_hs build_mc_chain_nvt_sw build_md_nvt_lj build_mc_nvt_lj_re build_mc_nvt_poly_lj build_md_nvt_lj_llle build_md_nve_lj_ll build_test_pot_bend build_md_nvt_lj_le build_test_pot_dd build_cluster build_mc_zvt_lj build_error_calc build_mc_nvt_lj_ll build_initialize build_eos_hs build_diffusion build_mc_nvt_sc build_t_tensor build_mc_npt_hs build_mc_npt_lj build_mc_nvt_lj build_mc_nvt_hs build_eos_lj build_hit_and_miss build_qmc_walk_sho build_sample_mean build_mc_npt_sc build_md_nve_lj_vl build_md_lj_mts build_mc_chain_nvt_cbmc_lj build_md_npt_lj build_dpd build_mc_gibbs_lj build_corfun build_mesh build_fft3dwrap build_wl_hist build_qmc_pi_lj
gfortran -o build_qmc_pi_lj/qmc_pi_lj.o -c -fdefault-real-8 -fall-intrinsics -std=f2008 -Wall -Ibuild_qmc_pi_lj -Jbuild_qmc_pi_lj build_qmc_pi_lj/qmc_pi_lj.f90
build_qmc_pi_lj/qmc_pi_lj.f90:151:49:

               rik(:) = random_translate_vector ( dr_max/box, r(:,i,k) ) ! Trial move to new position (in box=1 units) 
                                                 1
Error: Type mismatch in argument ‘dr_max’ at (1); passed REAL(8) to REAL(4)
build_qmc_pi_lj/qmc_pi_lj.f90:164:35:

                  IF ( metropolis ( delta ) ) THEN ! Accept Metropolis test
                                   1
Error: Type mismatch in argument ‘delta’ at (1); passed REAL(8) to REAL(4)
scons: *** [build_qmc_pi_lj/qmc_pi_lj.o] Error 1
scons: building terminated because of errors.

$ scons
scons: Reading SConscript files ...
scons: done reading SConscript files.
scons: Building targets ...
scons: building associated VariantDir targets: build_md_nve_lj_omp build_smc_nvt_lj build_md_chain_nve_lj build_qmc_pi_sho build_grint build_test_pot_qq build_md_chain_mts_lj build_test_pot_gb build_md_nve_lj build_adjust build_test_pot_twist build_mc_npt_lj_ll build_mc_zvt_lj_ll build_bd_nvt_lj build_pair_distribution build_test_pot_at build_mc_chain_wl_sw build_test_pot_dq build_md_nve_hs build_mc_chain_nvt_sw build_md_nvt_lj build_mc_nvt_lj_re build_mc_nvt_poly_lj build_md_nvt_lj_llle build_md_nve_lj_ll build_test_pot_bend build_md_nvt_lj_le build_test_pot_dd build_cluster build_mc_zvt_lj build_error_calc build_mc_nvt_lj_ll build_initialize build_eos_hs build_diffusion build_mc_nvt_sc build_t_tensor build_mc_npt_hs build_mc_npt_lj build_mc_nvt_lj build_mc_nvt_hs build_eos_lj build_hit_and_miss build_qmc_walk_sho build_sample_mean build_mc_npt_sc build_md_nve_lj_vl build_md_lj_mts build_mc_chain_nvt_cbmc_lj build_md_npt_lj build_dpd build_mc_gibbs_lj build_corfun build_mesh build_fft3dwrap build_wl_hist
gfortran -o build_wl_hist/wl_hist.o -c -fdefault-real-8 -fall-intrinsics -std=f2008 -Wall -Ibuild_wl_hist -Jbuild_wl_hist build_wl_hist/wl_hist.f90
gfortran -o build_wl_hist/wl_hist build_wl_hist/wl_hist.o
gfortran -o build_fft3dwrap/fft3dwrap.o -c -fdefault-real-8 -fall-intrinsics -std=f2008 -Wall -I/opt/local/include -Ibuild_fft3dwrap -Jbuild_fft3dwrap build_fft3dwrap/fft3dwrap.f90
f951: Warning: Nonexistent include directory ‘/opt/local/include’ [-Wmissing-include-dirs]
build_fft3dwrap/fft3dwrap.f90:9: Error: Can't open included file 'fftw3.f03'
scons: *** [build_fft3dwrap/fft3dwrap.o] Error 1
scons: building terminated because of errors.

The last one is quite trivial, however, some people use math libraries from repositories. I think, it's possible to include such condition in the file.

P.S. This is fascinating that your new book is coming.

I suggest adding topics such as molecular-dynamics, monte-carlo in the About section.

There's nothing in the GUIDE about this example. At the same time, we might consider adding a Python version.

Hi professor ,
in my mentioned link , the line 137 ,the np.real is used to get the real part of the complex ,but in your book ,P498 ,there is "where Im(z) is the imaginary part of a complex number".
is there some wrong with that?

This should be possible, and may be helpful. However, since the production runs need to be quite long, the resulting program will be of limited use.

This looks relatively straightforward and may be informative and useful to students, especially those unfamiliar with Fortran, since replica exchange is so widespread.

Not originally included in the production release, since the user would have to have a suitable MPI/Python package installed, and the interface is less standard than for Fortran and C. Examples are pypar and mpi4py. We would have to pick one, and could not help much with enquiries from users who could not get their own MPI/Python installation to work.

There appears to be a bug in the link-list MC code, possibly in mc_lj_ll_module.f90. This manifests as a drop in the acceptance rate of volume-change moves in mc_npt_lj_ll, in some cases to zero, giving unreliable results. While this is investigated, all the link-list codes must be considered suspect.

Currently the mc_zvt_lj.f90 program, with or without neighbour lists, handles variation of N in a very crude way. Arrays are dimensioned with size 2N at the point of reading in the initial starting configuration, and are unchanged thereafter. If the value of N approaches the upper limit of the array, the program stops with a diagnostic message. In the GUIDE, this is mentioned; the reason given is "for simplicity". However, it is not very satisfactory, and also not very complicated to introduce an array reallocation routine at the appropriate place. We should do this.