I need to constrain two water droplets, each having 5000 molecules with different residue names (SOL1 and SOL2) and atom names. I have recently tried pull code with the following MDP options:
The water droplet on the right side is evaporating fast and the distance between the COMs are not being constrained as the initial configuration. The final configuration is shown at right side after 1380 ps.
Kindly help me constrain the COMs of both the droplets at their initial configurations. I plan to find the free energy at such constratined distances between their COMs.
Please see the image below which is after 6.5 ns. Here, the smaller droplet is droplet which was towards the right side. It has evaporated and all the water molecules have been scattered inside the box. I have checked the pullx.xvg file, the distances are almost constant at the start but increases towards the end. The main problem is the evaporation. Kindly suggest regarding this.
If a droplet evaporates, its COM is no longer well defined as the distance of evaporated molecules jumps up and down with periodic boundary conditions.
So, how can I prevent this undesired evaporation? Only the right side droplet is evaporating with huge fluctuations in positions of water molecules, while the left side droplet is moving about the right droplet with comparatively much less fluctuations and no evaporation. I am using OPC water model for the simulations. Kindly suggest any changes to be made in the MDP file.
I remember that we also experienced problems with evaporating droplets when system sizes got very large. If a box had a long extension in, say, the x-direction, then water droplets near the origin would behave as expected, while droplets with high x positions would heat up, start to rotate and evaporate. Also, if you take an exact same droplet that runs fine in a small box, and you put it into a far corner of a much (factor 100) larger box, it would start to evaporate. This all happens in single precision, but not in double precision. I think single precision might not be enough to calculate differences of positions (as e.g. needed to fulfil the constraints within each water molecule) accurately enough for molecules with large x-, y-, or z coordinates when using floating point. So you could try to check whether running in double precision solves the issue. Are you running with or without thermostat and barostat?
Hi, if you are running single precision and setting tc-grps=system, you could also look at the temperature of each group (SOL1 and SOL2) separately (e.g. by using gmx traj). Another ‘signature’ of the single-precision issue that Carsten highlighted is the unphisical transfer of thermal energy from one group to the other, i.e. you should see a droplet heating up and the other cooling down.
Yes, the droplets were created individually by NVT for 50ns and later merged together into the same GRO file. The radial density of the droplets was around 1 g/cc.
you have to use the GROMACS version that is compiled in double precision. If it is provided on the cluster you use, it is probably named gmx_d instead of just gmx. Or did you compile it yourself? Then you will have to make an extra double precision build.
Sure, but I had expereinced quite extreme heating/cooling of water droplets when running single-precision in large domain. So I wouldn’t be surpised if the temperature of one of the two droplets is much greated than 300K (and consequently the other may be even frozen).
We have installed the GROMACS version 2021.3 in single precision. Would you kindly share how to make an extra double precision build or share any available source in internet. It will be really helpful for me.
I have simulated two systems based on your suggestions:
With the mdp option “tc-grps=System” and “ref_t = 300”
With the mdp option “tc-grps=SOL1 SOL2” & “ref_t = 300 300”, i.e., separate for both droplets.
Kindly observe the results of the temperature of both droplets with time. For case 1, as the temperature of one droplet is rising, the temperature of the other droping. Whereas, the temperature of both remains almost same (about 150K) thoguhout the simulation for case 2. However, the sum of the temperatures of both the droplets is about 300K.
I have encountered recently that as I use separate tc-grps for SOL1 and SOL2, the droplet with SOL1 residues breaks into 2 droplets after 4.3 ns of the simulation. The splitting takes place within 1.2 ns. I have used the following mdp options:
;TEMPERATURE coupling
Tcoupl = nose-hoover
tc-grps = SOL1 SOL2 ; Groups to couple separately
nsttcouple = -1
tau_t = 0.5 0.5 ; Time constant (ps) and
ref_t = 300.0 300.0 ; reference temperature (K)
; PULL CODE
pull = yes
pull_ngroups = 2
pull_ncoords = 1
pull_group1_name = SOL1
pull_group2_name = SOL2
pull_coord1_type = umbrella ; harmonic biasing force
pull_coord1_geometry = distance ; simple distance increase
pull_coord1_groups = 1 2
pull_coord1_dim = Y Y Y
;pull_coord1_rate = 0 ; nm/ps (-ve means in the -ve Z direction)
pull_coord1_k = 10 ; kJ mol^-1 nm^-2
pull_coord1_start = yes ; define initial COM distance > 0