@@ -458,44 +458,42 @@ Minimalist NVT input file
458
458
459
459
.. container:: justify
460
460
461
- Let us first perform a small (20 picoseconds)
461
+ Let us first perform a short (20 picoseconds)
462
462
equilibration in the NVT ensemble. In the NVT ensemble, the number of
463
463
atom (N) and volume (V) are maintained fixed, and the
464
464
temperature (T) is adjusted using a thermostat.
465
465
466
+ .. container:: justify
467
+
466
468
Let use write a new input script
467
- called nvt.mdp, and save it in the ' inputs/' folder.
469
+ called * nvt.mdp* , and save it in the * inputs/* folder.
468
470
Copy the following lines into it:
469
471
470
472
.. code-block:: bw
471
- :caption: *to be copied in inputs/nvt.mdp*
472
473
473
474
integrator = md
474
475
nsteps = 20000
475
476
dt = 0.001
476
477
477
478
.. container:: justify
478
479
479
- Here the molecular dynamics (md) integrator is used
480
+ Here, the molecular dynamics (md) integrator is used
480
481
(leapfrog algorithm), and a number of 20000 steps with
481
- timestep (dt) 0.001 ps is requested, so the total
482
- requested duration is 20 ps. Let us print the
483
- trajectory in a xtc file every 1 ps by adding:
482
+ a timestep *dt* equal of :math:`0.001 ~ \text{ps}` will be performed.
483
+ Let us print the trajectory in a compressed *xtc* file every 1 ps by adding:
484
484
485
485
.. code-block:: bw
486
- :caption: *to be copied in inputs/nvt.mdp*
487
486
488
487
nstxout-compressed = 1000
489
488
490
489
.. container:: justify
491
490
492
- Let us control the temperature over the course of the
491
+ Let us also control the temperature over the course of the
493
492
simulation using the v-rescale thermostat, which is
494
493
the Berendsen thermostat with an additional stochastic
495
- term (it is known to give proper canonical ensemble):
494
+ term :cite:`bussi2007canonical`. This thermostat is known to give proper canonical ensemble.
496
495
497
496
.. code-block:: bw
498
- :caption: *to be copied in inputs/nvt.mdp*
499
497
500
498
tcoupl = v-rescale
501
499
ref-t = 360
@@ -504,18 +502,17 @@ Minimalist NVT input file
504
502
505
503
.. container:: justify
506
504
507
- Here we also specified that the thermostat is
508
- applied to the entire system (we could choose to
509
- apply it only to a certain group of atom, which we
510
- will do later), and that the damping constant for
511
- the thermostat is 0.5 ps.
505
+ Here, we also specified that the thermostat is
506
+ applied to the entire system, and that the damping constant for
507
+ the thermostat is equal to 0.5 ps.
508
+
509
+ .. container:: justify
512
510
513
511
Note that the relatively high temperature of 360 K
514
512
has been chosen here in order to reduce the
515
513
viscosity of the solution and converge toward
516
- the desired result (i.e. the diffusion coefficient) faster.
517
- We now have a minimalist input file for performing
518
- the NVT step. Run it by typing in the terminal:
514
+ the desired result faster. We now have a minimalist input file for performing
515
+ the first NVT simulation. Run it by typing in the terminal:
519
516
520
517
.. code-block:: bw
521
518
@@ -524,9 +521,11 @@ Minimalist NVT input file
524
521
525
522
.. container:: justify
526
523
527
- Here ' -c min.gro' ensures that the previously
524
+ Here * -c min.gro* ensures that the previously
528
525
minimized configuration is used as a starting point.
529
526
527
+ .. container:: justify
528
+
530
529
After the completion of the simulation, we can
531
530
ensure that the system temperature indeed reached
532
531
the value of 360 K by using the energy command of
@@ -538,14 +537,16 @@ Minimalist NVT input file
538
537
539
538
.. container:: justify
540
539
541
- and choose ' temperature' .
540
+ and choose * temperature* .
542
541
543
- From the generated file, we can see that temperature
542
+ .. container:: justify
543
+
544
+ From the generated *Tnvt.xvg* file, we can see that temperature
544
545
started from 0, which was expected since the atoms
545
546
have no velocity during a minimization step, and
546
547
reaches a temperature slightly larger than the
547
548
requested 360 K after a duration of a few
548
- picoseconds:
549
+ picoseconds.
549
550
550
551
.. figure:: figures/bulksolution/temperature-light.png
551
552
:alt: Gromacs tutorial : temperature versus time.
@@ -555,27 +556,24 @@ Minimalist NVT input file
555
556
:alt: Gromacs tutorial : temperature versus time.
556
557
:class: only-dark
557
558
559
+ .. container:: figurelegend
560
+
558
561
Evolution of the temperature as a function of the time
559
562
during the NVT equilibration. Dashed line is the
560
563
requested temperature of 360 K.
561
564
562
- .. container:: justify
563
-
564
- A better control of the temperature is achieved in the next section.
565
-
566
565
Improving the NVT
567
566
=================
568
567
569
568
.. container:: justify
570
569
571
- So far, a very few commands have been placed in the
572
- NVT input file, meaning that most of the instruction
570
+ So far, very few commands have been placed in the
571
+ *.mdp* input file, meaning that most of the instruction
573
572
have been taken by GROMACS from the default
574
573
parameters. You can find what parameters were used
575
- during the last nvt run by opening the new nvt.mdp
576
- file that has been created (i.e. not the one in the
577
- ' inputs/' folder, but the one in the main folder).
578
- Exploring this new nvt.mdp file shows us that, for
574
+ during the last nvt run by opening the new *nvt.mdp*
575
+ file that has been created in the main folder.
576
+ Exploring this new *nvt.mdp* file shows us that, for
579
577
instance, plain cut-off Coulomb interactions have
580
578
been used:
581
579
@@ -586,9 +584,12 @@ Improving the NVT
586
584
587
585
.. container:: justify
588
586
589
- For this system, long range Coulomb interaction is a
587
+ For this system, long range Coulomb interaction would be a
590
588
better choice. We could also improve the thermostating of the
591
589
system by applying a separate thermostat to water molecules and ions.
590
+
591
+ .. container:: justify
592
+
592
593
Therefore, let us improve the NVT step by specifying
593
594
more options in the input file.
594
595
First, in the nvt.mdp file, let us impose the use of the
@@ -597,7 +598,6 @@ Improving the NVT
597
598
of 4, and cut-off of 4:
598
599
599
600
.. code-block:: bw
600
- :caption: *to be copied in inputs/nvt.mdp*
601
601
602
602
coulombtype = pme
603
603
fourierspacing = 0.1
@@ -611,7 +611,6 @@ Improving the NVT
611
611
Let us also specify the van der Waals interaction:
612
612
613
613
.. code-block:: bw
614
- :caption: *to be copied in inputs/nvt.mdp*
615
614
616
615
vdw-type = Cut-off
617
616
rvdw = 1.0
@@ -622,7 +621,6 @@ Improving the NVT
622
621
bonds of the water molecules:
623
622
624
623
.. code-block:: bw
625
- :caption: *to be copied in inputs/nvt.mdp*
626
624
627
625
constraint-algorithm = lincs
628
626
constraints = hbonds
@@ -631,11 +629,10 @@ Improving the NVT
631
629
.. container:: justify
632
630
633
631
Let us also use separate temperature baths for
634
- water and ions (here corresponding to the gromacs group called non-water) respectively
635
- by replacing:
632
+ water and ions (here corresponding to the GROMACS group called * non-water*)
633
+ respectively by replacing:
636
634
637
635
.. code-block:: bw
638
- :caption: *to be removed from inputs/nvt.mdp*
639
636
640
637
tcoupl = v-rescale
641
638
ref-t = 360
@@ -647,7 +644,6 @@ Improving the NVT
647
644
by:
648
645
649
646
.. code-block:: bw
650
- :caption: *to be copied in inputs/nvt.mdp*
651
647
652
648
tcoupl = v-rescale
653
649
tc-grps = Water non-Water
@@ -659,7 +655,6 @@ Improving the NVT
659
655
Let us specify neighbor searching parameters:
660
656
661
657
.. code-block:: bw
662
- :caption: *to be copied in inputs/nvt.mdp*
663
658
664
659
cutoff-scheme = Verlet
665
660
nstlist = 10
@@ -671,7 +666,6 @@ Improving the NVT
671
666
total velocities give the desired temperature instead of 0):
672
667
673
668
.. code-block:: bw
674
- :caption: *to be copied in inputs/nvt.mdp*
675
669
676
670
gen-vel = yes
677
671
gen-temp = 360
@@ -682,7 +676,6 @@ Improving the NVT
682
676
velocity of the whose system:
683
677
684
678
.. code-block:: bw
685
- :caption: *to be copied in inputs/nvt.mdp*
686
679
687
680
comm_mode = linear
688
681
comm_grps = system
@@ -703,6 +696,8 @@ Improving the NVT
703
696
:alt: Gromacs tutorial : temperature versus time.
704
697
:class: only-dark
705
698
699
+ .. contained:: figurelegend
700
+
706
701
Evolution of the temperature as a function of the time
707
702
during the NVT equilibration.
708
703
@@ -713,14 +708,13 @@ Adjust the density using NPT
713
708
714
709
Now that the system is properly equilibrated in the
715
710
NVT ensemble, let us perform an equilibration in the
716
- NPT ensemble, where the pressure is imposed and the
711
+ NPT ensemble where the pressure is imposed and the
717
712
volume of the box is free to relax. NPT relaxation ensures that the
718
713
density of the fluid converges toward its equilibrium value.
719
- Create a new input script, call it ' npt.mdp' , and
714
+ Create a new input script, call it * npt.mdp* , and
720
715
copy the following lines in it:
721
716
722
717
.. code-block:: bw
723
- :caption: *to be copied in inputs/npt.mdp*
724
718
725
719
integrator = md
726
720
nsteps = 50000
@@ -762,14 +756,14 @@ Adjust the density using NPT
762
756
763
757
.. container:: justify
764
758
765
- The main difference with the previous NVT script, is
759
+ The main difference with the previous NVT script is
766
760
the addition of the isotropic C-rescale pressure
767
761
coupling with a target pressure of 1 bar. Some other
768
- differences are the addition of the ' nstlog' and
769
- ' nstenergy' commands to control the frequency at
762
+ differences are the addition of the * nstlog* and
763
+ * nstenergy* commands to control the frequency at
770
764
which information are printed in the log file and in
771
- the energy file (edr), and the removing the
772
- ' gen-vel' commands. Run it using:
765
+ the energy file (* edr* ), and the removing the
766
+ * gen-vel* commands. Run it using:
773
767
774
768
.. code-block:: bash
775
769
@@ -779,7 +773,7 @@ Adjust the density using NPT
779
773
.. container:: justify
780
774
781
775
Let us have a look a both temperature, pressure and
782
- volume during the NPT step using the ' gmx energy'
776
+ volume during the NPT step using the * gmx energy*
783
777
command 3 times:
784
778
785
779
.. code-block:: bash
@@ -790,7 +784,7 @@ Adjust the density using NPT
790
784
791
785
.. container:: justify
792
786
793
- Choose respectively ' temperature' , ' pressure' and ' volume' .
787
+ Choose respectively * temperature*, * pressure* and * volume* .
794
788
This is what I see:
795
789
796
790
.. figure:: figures/bulksolution/NPT-light.png
@@ -801,6 +795,8 @@ Adjust the density using NPT
801
795
:alt: Gromacs tutorial : NPT equilibration
802
796
:class: only-dark
803
797
798
+ .. container:: figurelegend
799
+
804
800
From top to bottom: evolution of the temperature,
805
801
pressure, and volume of the simulation box as a
806
802
function of the time during the NPT equilibration.
@@ -817,6 +813,8 @@ Adjust the density using NPT
817
813
typical in molecular dynamics, particularly with
818
814
liquid water that is almost uncompressible.
819
815
816
+ .. container:: justify
817
+
820
818
Exact results may differ depending on the actual *.gro* file generated.
821
819
822
820
Measurement diffusion coefficient
@@ -828,15 +826,16 @@ Measurement diffusion coefficient
828
826
perform a longer simulation and extract quantities of
829
827
interest.
830
828
829
+ .. container:: justify
830
+
831
831
Here, as an illustration, the diffusion
832
832
coefficients of all 3 species (water and the two ions) will be
833
833
measured. First, let us perform a longer run in the
834
834
NVT ensemble. Create a new input file, call it
835
- ' pro.mdp' ( ' pro' is short for ' production' ), and copy the
835
+ * pro.mdp* (* pro* is short for * production* ), and copy the
836
836
following lines into it:
837
837
838
838
.. code-block:: bw
839
- :caption: *to be copied in inputs/pro.mdp*
840
839
841
840
integrator = md
842
841
nsteps = 200000
@@ -890,18 +889,24 @@ Measurement diffusion coefficient
890
889
891
890
.. container:: justify
892
891
893
- and select the SO4 ions by typing ' SO4' , and then press ' ctrl D' .
892
+ and select the SO4 ions by typing * SO4* , and then press * ctrl D* .
894
893
894
+ .. container:: justify
895
+
895
896
Fitting the slope of the
896
- MSD gives a value of 1.3e -5 cm\ :sup:`2`/s for the
897
+ MSD gives a value of :math:`1.3 \mathrm{e} -5 ~ \text{cm}^2/\text{s}` for the
897
898
diffusion coefficient.
899
+
900
+ .. container:: justify
898
901
899
902
Repeat the same for Na and water.
900
903
901
- For Na, the value is 1.5e-5
902
- cm\ :sup:`2`/s, and for water 5.2e-5 cm\ :sup:`2`/s
903
- (not too far from the experimental value of ~ 7e-5
904
- cm\ :sup:`2`/s at temperature T=360 K).
904
+ .. container:: justify
905
+
906
+ For Na, the value is :math:`1.5 \mathrm{e}-5 ~ \text{cm}^2/\text{s}`,
907
+ and for water :math:`5.2 \mathrm{e}-5 ~ \text{cm}^2/\text{s}`
908
+ (not too far from the experimental value of :math:`7.5 \mathrm{e}-5 ~ \text{cm}^2/\text{s}`
909
+ at temperature T=360 K).
905
910
906
911
.. admonition:: About MSD in molecular simulations
907
912
:class: info
@@ -935,5 +940,3 @@ Going further
935
940
Take advantage of the generated production run to extract more
936
941
equilibrium quantities. For instance, Gromacs allows you to
937
942
extract Radial Distribution Functions (RDF) using the *gmx rdf* commands.
938
-
939
- .. include:: ../contact/contactme.rst
0 commit comments