Let assume that I would like to create amorphopus polymer sample composed of C32H66 (dotriacontane) chains. I am trying to do it in the following way:

1. Create a single molecule (e.g. in Avogadro), save the yxz coordinates.

2. Use e.g. fortran code and create from the xyz file the pdb file. Preview of some middle lines of the file is below (I see that formating is lost):

HETATM 31 C31 DTCN 1 31.632 16.485 -14.222 1.00 0.00 C

HETATM 32 C32 DTCN 1 33.154 16.573 -14.167 1.00 0.00 C

HETATM 33 H011 DTCN 1 0.328 -0.862 0.798 1.00 0.00 H

HETATM 34 H012 DTCN 1 0.344 -0.416 -0.951 1.00 0.00 H

3. Copy to my working catalog and modify the oplsaa.ff force field folder (as I have chosen this potential), especially aminoacids.rtp and aminoacids.hdb files (the names probably depends on the Gromacs distributions, I am using version 4.5.3):

aminoacids.rtp

; C32H66 Dotriacontane

[ DTCN ]

[ atoms ]

C01 opls_135 -0.180 1

C02 opls_136 -0.120 2

C03 opls_136 -0.120 3

C04 opls_136 -0.120 4

C05 opls_136 -0.120 5

C06 opls_136 -0.120 6

C07 opls_136 -0.120 7

C08 opls_136 -0.120 8

C09 opls_136 -0.120 9

C10 opls_136 -0.120 10

C11 opls_136 -0.120 11

C12 opls_136 -0.120 12

C13 opls_136 -0.120 13

C14 opls_136 -0.120 14

C15 opls_136 -0.120 15

C16 opls_136 -0.120 16

C17 opls_136 -0.120 17

C18 opls_136 -0.120 18

C19 opls_136 -0.120 19

C20 opls_136 -0.120 20

C21 opls_136 -0.120 21

C22 opls_136 -0.120 22

C23 opls_136 -0.120 23

C24 opls_136 -0.120 24

C25 opls_136 -0.120 25

C26 opls_136 -0.120 26

C27 opls_136 -0.120 27

C28 opls_136 -0.120 28

C29 opls_136 -0.120 29

C30 opls_136 -0.120 30

C31 opls_136 -0.120 31

C32 opls_135 -0.180 32

H011 opls_156 0.060 1

H012 opls_156 0.060 1

H013 opls_156 0.060 1

H021 opls_156 0.060 2

H022 opls_156 0.060 2

H031 opls_156 0.060 3

H032 opls_156 0.060 3

H041 opls_156 0.060 4

H042 opls_156 0.060 4

H051 opls_156 0.060 5

H052 opls_156 0.060 5

H061 opls_156 0.060 6

H062 opls_156 0.060 6

H071 opls_156 0.060 7

H072 opls_156 0.060 7

H081 opls_156 0.060 8

H082 opls_156 0.060 8

H091 opls_156 0.060 9

H092 opls_156 0.060 9

H101 opls_156 0.060 10

H102 opls_156 0.060 10

H111 opls_156 0.060 11

H112 opls_156 0.060 11

H121 opls_156 0.060 12

H122 opls_156 0.060 12

H131 opls_156 0.060 13

H132 opls_156 0.060 13

H141 opls_156 0.060 14

H142 opls_156 0.060 14

H151 opls_156 0.060 15

H152 opls_156 0.060 15

H161 opls_156 0.060 16

H162 opls_156 0.060 16

H171 opls_156 0.060 17

H172 opls_156 0.060 17

H181 opls_156 0.060 18

H182 opls_156 0.060 18

H191 opls_156 0.060 19

H192 opls_156 0.060 19

H201 opls_156 0.060 20

H202 opls_156 0.060 20

H211 opls_156 0.060 21

H212 opls_156 0.060 21

H221 opls_156 0.060 22

H222 opls_156 0.060 22

H231 opls_156 0.060 23

H232 opls_156 0.060 23

H241 opls_156 0.060 24

H242 opls_156 0.060 24

H251 opls_156 0.060 25

H252 opls_156 0.060 25

H261 opls_156 0.060 26

H262 opls_156 0.060 26

H271 opls_156 0.060 27

H272 opls_156 0.060 27

H281 opls_156 0.060 28

H282 opls_156 0.060 28

H291 opls_156 0.060 29

H292 opls_156 0.060 29

H301 opls_156 0.060 30

H302 opls_156 0.060 30

H311 opls_156 0.060 31

H312 opls_156 0.060 31

H321 opls_156 0.060 32

H322 opls_156 0.060 32

H323 opls_156 0.060 32

[ bonds ]

C01 C02

C02 C03

C03 C04

C04 C05

C05 C06

C06 C07

C07 C08

C08 C09

C09 C10

C10 C11

C11 C12

C12 C13

C13 C14

C14 C15

C15 C16

C16 C17

C17 C18

C18 C19

C19 C20

C20 C21

C21 C22

C22 C23

C23 C24

C24 C25

C25 C26

C26 C27

C27 C28

C28 C29

C29 C30

C30 C31

C31 C32

H011 C01

H012 C01

H013 C01

H021 C02

H022 C02

H031 C03

H032 C03

H041 C04

H042 C04

H051 C05

H052 C05

H061 C06

H062 C06

H071 C07

H072 C07

H081 C08

H082 C08

H091 C09

H092 C09

H101 C10

H102 C10

H111 C11

H112 C11

H121 C12

H122 C12

H131 C13

H132 C13

H141 C14

H142 C14

H151 C15

H152 C15

H161 C16

H162 C16

H171 C17

H172 C17

H181 C18

H182 C18

H191 C19

H192 C19

H201 C20

H202 C20

H211 C21

H212 C21

H221 C22

H222 C22

H231 C23

H232 C23

H241 C24

H242 C24

H251 C25

H252 C25

H261 C26

H262 C26

H271 C27

H272 C27

H281 C28

H282 C28

H291 C29

H292 C29

H301 C30

H302 C30

H311 C31

H312 C31

H321 C32

H322 C32

H323 C32

aminoacids.hdb

DTCN 32

3 4 H01 C01 C02 C03

2 4 H02 C02 C03 C04

2 4 H03 C03 C04 C05

2 4 H04 C04 C05 C06

2 4 H05 C05 C06 C07

2 4 H06 C06 C07 C08

2 4 H07 C07 C08 C09

2 4 H08 C08 C09 C10

2 4 H09 C09 C10 C11

2 4 H10 C10 C11 C12

2 4 H11 C11 C12 C13

2 4 H12 C12 C13 C14

2 4 H13 C13 C14 C15

2 4 H14 C14 C15 C16

2 4 H15 C15 C16 C17

2 4 H16 C16 C17 C18

2 4 H17 C17 C18 C19

2 4 H18 C18 C19 C20

2 4 H19 C19 C20 C21

2 4 H20 C20 C21 C22

2 4 H21 C21 C22 C23

2 4 H22 C22 C23 C24

2 4 H23 C23 C24 C25

2 4 H24 C24 C25 C26

2 4 H25 C25 C26 C27

2 4 H26 C26 C27 C28

2 4 H27 C27 C28 C29

2 4 H28 C28 C29 C30

2 4 H29 C29 C30 C31

2 4 H30 C30 C31 C32

2 4 H31 C31 C30 C29

3 4 H32 C32 C31 C30

4. Use pdb2gmx command to create .gro file:

pdb2gmx -f X.pdb

and select forcefield from the working catalog and no water

5. Use genbox command to create larger random sample:

genbox -ci conf.gro -nmol 30 -box 4.8 4.8 12 -o box.gro

6. Create (e.g. in Fortran) pdb file for this larger system:

HETATM 1 C01 DTCN 1 16.290 5.740 69.430 1.00 0.00 C

HETATM 2 C02 DTCN 1 17.270 6.150 68.320 1.00 0.00 C

HETATM 3 C03 DTCN 1 16.780 5.670 66.950 1.00 0.00 C

HETATM 4 C04 DTCN 1 17.770 6.080 65.850 1.00 0.00 C

HETATM 5 C05 DTCN 1 17.280 5.590 64.480 1.00 0.00 C

HETATM 6 C06 DTCN 1 18.250 6.000 63.370 1.00 0.00 C

HETATM 7 C07 DTCN 1 17.770 5.520 61.990 1.00 0.00 C

HETATM 8 C08 DTCN 1 18.760 5.930 60.890 1.00 0.00 C

HETATM 9 C09 DTCN 1 18.270 5.450 59.520 1.00 0.00 C

HETATM 10 C10 DTCN 1 19.250 5.860 58.420 1.00 0.00 C

HETATM 11 C11 DTCN 1 18.760 5.380 57.050 1.00 0.00 C

HETATM 12 C12 DTCN 1 19.740 5.790 55.940 1.00 0.00 C

HETATM 13 C13 DTCN 1 19.260 5.310 54.560 1.00 0.00 C

HETATM 14 C14 DTCN 1 20.230 5.720 53.460 1.00 0.00 C

HETATM 15 C15 DTCN 1 19.750 5.240 52.090 1.00 0.00 C

HETATM 16 C16 DTCN 1 20.720 5.650 50.980 1.00 0.00 C

HETATM 17 C17 DTCN 1 20.240 5.170 49.600 1.00 0.00 C

HETATM 18 C18 DTCN 1 21.220 5.590 48.500 1.00 0.00 C

HETATM 19 C19 DTCN 1 20.740 5.100 47.130 1.00 0.00 C

HETATM 20 C20 DTCN 1 21.710 5.520 46.020 1.00 0.00 C

HETATM 21 C21 DTCN 1 21.230 5.030 44.660 1.00 0.00 C

HETATM 22 C22 DTCN 1 22.200 5.450 43.550 1.00 0.00 C

HETATM 23 C23 DTCN 1 21.720 4.960 42.170 1.00 0.00 C

HETATM 24 C24 DTCN 1 22.700 5.380 41.070 1.00 0.00 C

HETATM 25 C25 DTCN 1 22.210 4.890 39.700 1.00 0.00 C

HETATM 26 C26 DTCN 1 23.190 5.310 38.590 1.00 0.00 C

HETATM 27 C27 DTCN 1 22.710 4.810 37.220 1.00 0.00 C

HETATM 28 C28 DTCN 1 23.690 5.230 36.110 1.00 0.00 C

HETATM 29 C29 DTCN 1 23.210 4.730 34.740 1.00 0.00 C

HETATM 30 C30 DTCN 1 24.180 5.160 33.630 1.00 0.00 C

HETATM 31 C31 DTCN 1 23.710 4.650 32.270 1.00 0.00 C

HETATM 32 C32 DTCN 1 24.670 5.070 31.170 1.00 0.00 C

HETATM 33 H011 DTCN 1 16.670 6.100 70.400 1.00 0.00 H

HETATM 34 H012 DTCN 1 16.200 4.640 69.480 1.00 0.00 H

HETATM 35 H013 DTCN 1 15.300 6.180 69.240 1.00 0.00 H

HETATM 36 H021 DTCN 1 18.270 5.720 68.550 1.00 0.00 H

HETATM 37 H022 DTCN 1 17.360 7.250 68.320 1.00 0.00 H

HETATM 38 H031 DTCN 1 15.780 6.110 66.740 1.00 0.00 H

HETATM 39 H032 DTCN 1 16.690 4.560 66.970 1.00 0.00 H

HETATM 40 H041 DTCN 1 18.770 5.640 66.060 1.00 0.00 H

HETATM 41 H042 DTCN 1 17.860 7.180 65.830 1.00 0.00 H

HETATM 42 H051 DTCN 1 16.280 6.030 64.260 1.00 0.00 H

HETATM 43 H052 DTCN 1 17.190 4.490 64.500 1.00 0.00 H

HETATM 44 H061 DTCN 1 19.260 5.570 63.580 1.00 0.00 H

HETATM 45 H062 DTCN 1 18.350 7.110 63.360 1.00 0.00 H

HETATM 46 H071 DTCN 1 16.770 5.960 61.790 1.00 0.00 H

HETATM 47 H072 DTCN 1 17.680 4.410 62.010 1.00 0.00 H

HETATM 48 H081 DTCN 1 19.750 5.490 61.110 1.00 0.00 H

HETATM 49 H082 DTCN 1 18.850 7.040 60.880 1.00 0.00 H

HETATM 50 H091 DTCN 1 17.270 5.890 59.310 1.00 0.00 H

HETATM 51 H092 DTCN 1 18.170 4.340 59.530 1.00 0.00 H

HETATM 52 H101 DTCN 1 20.250 5.420 58.620 1.00 0.00 H

HETATM 53 H102 DTCN 1 19.340 6.970 58.400 1.00 0.00 H

HETATM 54 H111 DTCN 1 17.760 5.820 56.830 1.00 0.00 H

HETATM 55 H112 DTCN 1 18.670 4.270 57.050 1.00 0.00 H

HETATM 56 H121 DTCN 1 20.740 5.350 56.150 1.00 0.00 H

HETATM 57 H122 DTCN 1 19.830 6.900 55.920 1.00 0.00 H

HETATM 58 H131 DTCN 1 18.250 5.750 54.350 1.00 0.00 H

HETATM 59 H132 DTCN 1 19.160 4.200 54.580 1.00 0.00 H

HETATM 60 H141 DTCN 1 21.240 5.280 53.660 1.00 0.00 H

HETATM 61 H142 DTCN 1 20.320 6.830 53.450 1.00 0.00 H

HETATM 62 H151 DTCN 1 18.740 5.680 51.880 1.00 0.00 H

HETATM 63 H152 DTCN 1 19.660 4.130 52.090 1.00 0.00 H

HETATM 64 H161 DTCN 1 21.730 5.220 51.190 1.00 0.00 H

HETATM 65 H162 DTCN 1 20.810 6.760 50.960 1.00 0.00 H

HETATM 66 H171 DTCN 1 19.240 5.610 49.390 1.00 0.00 H

HETATM 67 H172 DTCN 1 20.150 4.060 49.620 1.00 0.00 H

HETATM 68 H181 DTCN 1 22.210 5.150 48.710 1.00 0.00 H

HETATM 69 H182 DTCN 1 21.300 6.690 48.490 1.00 0.00 H

HETATM 70 H191 DTCN 1 19.740 5.540 46.920 1.00 0.00 H

HETATM 71 H192 DTCN 1 20.650 3.990 47.140 1.00 0.00 H

HETATM 72 H201 DTCN 1 22.710 5.090 46.230 1.00 0.00 H

HETATM 73 H202 DTCN 1 21.790 6.630 46.010 1.00 0.00 H

HETATM 74 H211 DTCN 1 20.230 5.470 44.450 1.00 0.00 H

HETATM 75 H212 DTCN 1 21.140 3.920 44.670 1.00 0.00 H

HETATM 76 H221 DTCN 1 23.210 5.020 43.750 1.00 0.00 H

HETATM 77 H222 DTCN 1 22.290 6.560 43.530 1.00 0.00 H

HETATM 78 H231 DTCN 1 20.720 5.390 41.960 1.00 0.00 H

HETATM 79 H232 DTCN 1 21.640 3.850 42.200 1.00 0.00 H

HETATM 80 H241 DTCN 1 23.700 4.950 41.280 1.00 0.00 H

HETATM 81 H242 DTCN 1 22.780 6.490 41.050 1.00 0.00 H

HETATM 82 H251 DTCN 1 21.210 5.320 39.480 1.00 0.00 H

HETATM 83 H252 DTCN 1 22.140 3.780 39.710 1.00 0.00 H

HETATM 84 H261 DTCN 1 24.190 4.880 38.800 1.00 0.00 H

HETATM 85 H262 DTCN 1 23.270 6.410 38.570 1.00 0.00 H

HETATM 86 H271 DTCN 1 21.710 5.230 37.010 1.00 0.00 H

HETATM 87 H272 DTCN 1 22.630 3.700 37.240 1.00 0.00 H

HETATM 88 H281 DTCN 1 24.700 4.810 36.330 1.00 0.00 H

HETATM 89 H282 DTCN 1 23.760 6.340 36.090 1.00 0.00 H

HETATM 90 H291 DTCN 1 22.210 5.150 34.530 1.00 0.00 H

HETATM 91 H292 DTCN 1 23.140 3.630 34.760 1.00 0.00 H

HETATM 92 H301 DTCN 1 25.190 4.740 33.850 1.00 0.00 H

HETATM 93 H302 DTCN 1 24.260 6.260 33.620 1.00 0.00 H

HETATM 94 H311 DTCN 1 22.700 5.080 32.040 1.00 0.00 H

HETATM 95 H312 DTCN 1 23.640 3.550 32.280 1.00 0.00 H

HETATM 96 H321 DTCN 1 24.750 6.180 31.110 1.00 0.00 H

HETATM 97 H322 DTCN 1 24.300 4.690 30.180 1.00 0.00 H

HETATM 98 H323 DTCN 1 25.680 4.640 31.340 1.00 0.00 H

HETATM 99 C01 DTCN 2 17.400 33.830 107.160 1.00 0.00 C

HETATM 100 C02 DTCN 2 18.030 32.720 106.310 1.00 0.00 C

...

7. Use pdb2gmx to obtain .gro and .topol files for this larger system

8. Prepare .mdp file for energy minimisation:

; LINES STARTING WITH ';' ARE COMMENTS

title = Minimization of Hen Egg White Lysozyme (1AKI.pdb) ADD YOUR NAME AND DATE HERE ; Title of run

; The following lines tell the program the standard locations where to find certain files

cpp = /lib/cpp ; Preprocessor

include = -I../top ; Directories to include in the topology format

; Parameters describing what to do, when to stop and what to save

integrator = steep ; Algorithm options

; steep = steepest descent minimization

; MD = Leap Frog algorith for integrating Newtons equations of motion )

emtol = 10.0 ; Stop minimization when the energy changes by less than emtol kJ/mol.

nsteps = 1000 ; Maximum number of (minimization) steps to perform

nstenergy = 10 ; Write energies to disk every nstenergy steps

nstxtcout = 10 ; Write coordinates to disk every nstxtcout steps

xtc_grps = DTCN ; Which coordinate group(s) to write to disk

energygrps = DTCN ; Which energy group(s) to write to disk

; Parameters describing how to find the neighbors of each atom and how to calculate the interactions

nstlist = 5 ; Frequency to update the neighbor list and long range forces

ns_type = grid ; Method to determine neighbor list (simple, grid)

rlist = 1.4 ; Cut-off for making neighbor list (short range forces)

coulombtype = cut-off ; Treatment of long range electrostatic interactions

rcoulomb = 1.4 ; long range electrostatic cut-off

rvdw = 1.4 ; long range Van der Waals cut-off

constraints = none ; Bond types to replace by constraints

pbc = no ; Periodic Boundary Conditions (yes/no)

9. Modify .topol file:

1 opls_135 1 DTCN C01 1 -0.18 12.01

2 opls_136 1 DTCN C02 2 -0.12 12.01

3 opls_136 1 DTCN C03 3 -0.12 12.01

4 opls_136 1 DTCN C04 4 -0.12 12.01

5 opls_136 1 DTCN C05 5 -0.12 12.01

6 opls_136 1 DTCN C06 6 -0.12 12.01

7 opls_136 1 DTCN C07 7 -0.12 12.01

8 opls_136 1 DTCN C08 8 -0.12 12.01

9 opls_136 1 DTCN C09 9 -0.12 12.01

10 opls_136 1 DTCN C10 10 -0.12 12.01

11 opls_136 1 DTCN C11 11 -0.12 12.01

12 opls_136 1 DTCN C12 12 -0.12 12.01

13 opls_136 1 DTCN C13 13 -0.12 12.01

14 opls_136 1 DTCN C14 14 -0.12 12.01

15 opls_136 1 DTCN C15 15 -0.12 12.01

16 opls_136 1 DTCN C16 16 -0.12 12.01

17 opls_136 1 DTCN C17 17 -0.12 12.01

18 opls_136 1 DTCN C18 18 -0.12 12.01

19 opls_136 1 DTCN C19 19 -0.12 12.01

20 opls_136 1 DTCN C20 20 -0.12 12.01

21 opls_136 1 DTCN C21 21 -0.12 12.01

22 opls_136 1 DTCN C22 22 -0.12 12.01

23 opls_136 1 DTCN C23 23 -0.12 12.01

24 opls_136 1 DTCN C24 24 -0.12 12.01

25 opls_136 1 DTCN C25 25 -0.12 12.01

26 opls_136 1 DTCN C26 26 -0.12 12.01

27 opls_136 1 DTCN C27 27 -0.12 12.01

28 opls_136 1 DTCN C28 28 -0.12 12.01

29 opls_136 1 DTCN C29 29 -0.12 12.01

30 opls_136 1 DTCN C30 30 -0.12 12.01

31 opls_136 1 DTCN C31 31 -0.12 12.01

32 opls_135 1 DTCN C32 32 -0.18 12.01

33 opls_156 1 DTCN H011 1 0.06 1.01 (without modification I got here 33 instead of 1 in the 6th column)

34 opls_156 1 DTCN H012 1 0.06 1.01

35 opls_156 1 DTCN H013 1 0.06 1.01

36 opls_156 1 DTCN H021 2 0.06 1.01

37 opls_156 1 DTCN H022 2 0.06 1.01

10. Use grompp command to obtain tpr file:

grompp -f minim.mdp -c mol30.gro -p chg.top

11. Perform minimisation by:

mdrun

command

12. Copy confout.gro, chg.top and Jobscript file to another folder and prepare .mdp file for NVT relaxation:

; VARIOUS PREPROCESSING OPTIONS =

title = Yo

cpp = /lib/cpp

include = -I../top

define =

; RUN CONTROL PARAMETERS =

integrator = md

; start time and timestep in ps =

tinit = 0

dt = 0.002

nsteps = 25000

; number of steps for center of mass motion removal =

nstcomm = 1

comm-grps =

; LANGEVIN DYNAMICS OPTIONS =

; Temperature, friction coefficient (amu/ps) and random seed =

bd-fric = 0

ld-seed = 1993

; OUTPUT CONTROL OPTIONS =

; Output frequency for coords (x), velocities (v) and forces (f) =

nstxout = 1000

nstvout = 0

nstfout = 0

; Output frequency for energies to log file and energy file =

nstlog = 1000

nstenergy = 1000

; Output frequency and precision for xtc file =

; This selects the subset of atoms for the xtc file. You can =

; select multiple groups. By default all atoms will be written. =

; Selection of energy groups =

energygrps = DTCN

; NEIGHBORSEARCHING PARAMETERS =

; nblist update frequency =

nstlist = 10

; ns algorithm (simple or grid) =

ns-type = grid

; Periodic boundary conditions: xyz or none =

pbc = xy

; nblist cut-off =

rlist = 1.4

; OPTIONS FOR ELECTROSTATICS AND VDW =

; Method for doing electrostatics =

coulombtype = Cut-off

rcoulomb = 1.4

; Dielectric constant (DC) for cut-off or DC of reaction field =

epsilon-r = 1

; Method for doing Van der Waals =

vdw-type = Cut-off

; cut-off lengths =

rvdw-switch = 0

rvdw = 1.4

; Apply long range dispersion corrections for Energy and Pressure =

DispCorr = no

; Spacing for the PME/PPPM FFT grid =

fourierspacing = 0.12

; FFT grid size, when a value is 0 fourierspacing will be used =

fourier-nx = 0

fourier-ny = 0

fourier-nz = 0

; EWALD/PME/PPPM parameters =

pme-order = 4

ewald-rtol = 1e-05

epsilon-surface = 0

; OPTIONS FOR WEAK COUPLING ALGORITHMS =

; Temperature coupling =

tcoupl = V-rescale

; Groups to couple separately =

tc-grps = System

; Time constant (ps) and reference temperature (K) =

tau-t = 0.1

ref-t = 300.0

; Pressure coupling =

Pcoupl = no

; SIMULATED ANNEALING CONTROL =

annealing = no

; Time at which temperature should be zero (ps) =

; GENERATE VELOCITIES FOR STARTUP RUN =

gen-vel = yes

gen-temp = 300

gen-seed = -1

; OPTIONS FOR BONDS =

constraints = all-bonds

; Type of constraint algorithm =

constraint-algorithm = Lincs

; Do not constrain the start configuration =

unconstrained-start = no

; Relative tolerance of shake =

shake-tol = 0.0001

; Highest order in the expansion of the constraint coupling matrix =

lincs-order = 4

; Lincs will write a warning to the stderr if in one step a bond =

; rotates over more degrees than =

lincs-warnangle = 30

; Convert harmonic bonds to morse potentials =

morse = no

; NMR refinement stuff =

; Distance restraints type: No, Simple or Ensemble =

disre = No

; Force weighting of pairs in one distance restraint: Equal or Conservative =

disre-weighting = Equal

; Use sqrt of the time averaged times the instantaneous violation =

disre-mixed = no

disre-fc = 1000

disre-tau = 0

; Output frequency for pair distances to energy file =

nstdisreout = 100

; Free energy control stuff =

free-energy = no

init-lambda = 0

delta-lambda = 0

sc-alpha = 0

sc-sigma = 0.3

; Non-equilibrium MD stuff =

acc-grps =

accelerate =

freezegrps =

freezedim =

cos-acceleration = 0

energygrp-excl =

; Electric fields =

; Format is number of terms (int) and for all terms an amplitude (real) =

; and a phase angle (real) =

E-x =

E-xt =

E-y =

E-yt =

E-z =

E-zt =

; User defined thingies =

user1-grps =

user2-grps =

userint1 = 0

userint2 = 0

userint3 = 0

userint4 = 0

userreal1 = 0

userreal2 = 0

userreal3 = 0

userreal4 = 0

13. Copy confout.gro, chg.top and Jobscript file to another folder and prepare .mdp file for NPT relaxation:

; VARIOUS PREPROCESSING OPTIONS =

title = Yo

cpp = /lib/cpp

include = -I../top

define =

; RUN CONTROL PARAMETERS =

integrator = md

; start time and timestep in ps =

tinit = 0

dt = 0.002

nsteps = 50000

; number of steps for center of mass motion removal =

nstcomm = 1

comm-grps =

; LANGEVIN DYNAMICS OPTIONS =

; Temperature, friction coefficient (amu/ps) and random seed =

bd-fric = 0

ld-seed = 1993

; OUTPUT CONTROL OPTIONS =

; Output frequency for coords (x), velocities (v) and forces (f) =

nstxout = 1000

nstvout = 0

nstfout = 0

; Output frequency for energies to log file and energy file =

nstlog = 1000

nstenergy = 1000

; Output frequency and precision for xtc file =

; This selects the subset of atoms for the xtc file. You can =

; select multiple groups. By default all atoms will be written. =

; Selection of energy groups =

energygrps = DTCN

; NEIGHBORSEARCHING PARAMETERS =

; nblist update frequency =

nstlist = 10

; ns algorithm (simple or grid) =

ns-type = grid

; Periodic boundary conditions: xyz or none =

pbc = xy

; nblist cut-off =

rlist = 1.4

nwall = 2

wall-atomtype = C C

wall-type = 12-6

wall-r-linpot = -1.0

wall-ewald-zfac = 3

; OPTIONS FOR ELECTROSTATICS AND VDW =

; Method for doing electrostatics =

coulombtype = Cut-off

rcoulomb = 1.4

; Dielectric constant (DC) for cut-off or DC of reaction field =

epsilon-r = 1

; Method for doing Van der Waals =

vdw-type = Cut-off

; cut-off lengths =

rvdw-switch = 0

rvdw = 1.4

; Apply long range dispersion corrections for Energy and Pressure =

DispCorr = EnerPres

; Spacing for the PME/PPPM FFT grid =

fourierspacing = 0.16

; FFT grid size, when a value is 0 fourierspacing will be used =

fourier-nx = 0

fourier-ny = 0

fourier-nz = 0

; EWALD/PME/PPPM parameters =

pme-order = 4

ewald-rtol = 1e-05

epsilon-surface = 0

ewald-geometry = 3dc

; OPTIONS FOR WEAK COUPLING ALGORITHMS =

; Temperature coupling =

tcoupl = V-rescale

; Groups to couple separately =

tc-grps = System

; Time constant (ps) and reference temperature (K) =

tau-t = 0.1

ref-t = 300.0

; Pressure coupling =

Pcoupl = Parrinello-Rahman

Pcoupltype = Isotropic

; Time constant (ps), compressibility (1/bar) and reference P (bar) =

tau-p = 2.0

compressibility = 4.5e-5

ref-p = 1.0

; SIMULATED ANNEALING CONTROL =

annealing = no

; Time at which temperature should be zero (ps) =

; GENERATE VELOCITIES FOR STARTUP RUN =

gen-vel = no

gen-temp = 300

gen-seed = -1

; OPTIONS FOR BONDS =

constraints = all-bonds

; Type of constraint algorithm =

constraint-algorithm = Lincs

; Do not constrain the start configuration =

unconstrained-start = yes

; Relative tolerance of shake =

shake-tol = 0.0001

; Highest order in the expansion of the constraint coupling matrix =

lincs-order = 4

; Lincs will write a warning to the stderr if in one step a bond =

; rotates over more degrees than =

lincs-warnangle = 30

; Convert harmonic bonds to morse potentials =

morse = no

; NMR refinement stuff =

; Distance restraints type: No, Simple or Ensemble =

disre = No

; Force weighting of pairs in one distance restraint: Equal or Conservative =

disre-weighting = Equal

; Use sqrt of the time averaged times the instantaneous violation =

disre-mixed = no

disre-fc = 1000

disre-tau = 0

; Output frequency for pair distances to energy file =

nstdisreout = 100

; Free energy control stuff =

free-energy = no

init-lambda = 0

delta-lambda = 0

sc-alpha = 0

sc-sigma = 0.3

; Non-equilibrium MD stuff =

acc-grps =

accelerate =

freezegrps =

freezedim =

cos-acceleration = 0

energygrp-excl =

; Electric fields =

; Format is number of terms (int) and for all terms an amplitude (real) =

; and a phase angle (real) =

E-x =

E-xt =

E-y =

E-yt =

E-z =

E-zt =

; User defined thingies =

user1-grps =

user2-grps =

userint1 = 0

userint2 = 0

userint3 = 0

userint4 = 0

userreal1 = 0

userreal2 = 0

userreal3 = 0

userreal4 = 0

Remarks:

A. I am using only pbc=xy periodic boundary conditions as I would like next make a bombardment of my sample, so I need surface and not broken (in the case of pbc) molecules in Z-direction.

B. Till the NPT relaxation this procedure works well. With this npt.mdp file the NPT relaxation also works, but I am not sure if such content of this file is the optimal one. I encountered a lot of error when trying some other options (e.g. PME algoritm): problem with spliting calclations to too many nodes (or threads), some segmentation fault errors, matrix inversion errors - all of them I found in the websites about Gromacs but I managed to run only with this configuration. The only error I have here is the "The domain decomposition grid has shifted too much in the Z-direction around cell 1 0 7" error (also reported somewhere in the web) but after restarting from checkpoint Gromacs managed to finish the calculations.

So, are my .mdp files (for minimisation, NVT and NPT relaxation) "ok" or something needs to be changed?

Thank you for reading all of this :)