From: Götz, Alexander (alexander.goetz_at_mytum.de)
Date: Mon Mar 27 2017 - 11:56:09 CDT
Dear NAMD users,
I am currently facing some troubles with my simulation and couldn’t find any more solution on the mailing list. So I hope you have a solution, idea of a good hint.
I want to simulate a 35 residue peptide in TFE/TIP3 (equal numbers) in a NPT ensemble. The box has an initial size of 104*60*60 Angstrom and the perfectly alpha-helical peptide is aligned along the X-Axis. We used this box in the past with peptides of up to 31 residues without any problems for several years. However, the new peptides will certainly unfold a little bite more than the old ones (approx. 6-8 residues) and therefore the 60 Angstrom of the y and z-axis will very likely become too short.
Originally the peptide is heavily tumbling in the box. To keep the box dimensions, I wanted to get rid of the rotational and translational movements of the peptide by using colvars. In the beginning I introduced a rotational and a translational restrains. As a result, I ended up with a heavily fluctuating/breathing box, what finnaly lead to the well-known partitioning error. Also with margin 2.0 this happens but what is much more problematic is that my box sometimes gets so small on one site, that there are only a few Angstroms between the peptide and the boundaries. This also stops me from just increasing the marign to values above 2.0.
Googling the mailing list, I found some tips for using slower rates of pressure coupling during equilibration. So I have set up the following equilibration protocol now (I only write down what I changed from the original CHARMM GUI protocols for NVT equilibration and NPT run in water solvents):
1. NVT, 295K, reassignFreq 500, reassingTemp 500, timestep 1.0, langevinDamping 5.0, run 500000
2. NPT 295K, FlexibleCell yes, wrapNearest off, margin 2.0, reassignFreq 500, reassingTemp 295, timestep 1.0, langevinDamping 5.0, langevinPistonPeriod 1000, lagevinPistonDecay 500, run 500000
3. NPT 295K, FlexibleCell yes, wrapNearest off, margin 2.0, reassignFreq 500, reassingTemp 295, timestep 1.0, langevinDamping 5.0, langevinPistonPeriod 750, lagevinPistonDecay 250, run 500000
4. NPT 295K, FlexibleCell yes, wrapNearest off, margin 2.0, reassignFreq 500, reassingTemp 295, timestep 1.0, langevinDamping 5.0, langevinPistonPeriod 200, lagevinPistonDecay 100, run 500000
5. NPT 295K, FlexibleCell yes, wrapNearest off, margin 2.0, reassignFreq 500, reassingTemp 295, timestep 2.0, langevinDamping 5.0, langevinPistonPeriod 100, lagevinPistonDecay 50, run 500000
The Production run works with the same settings like the last step of the equilibration, without reassigning the temperature every 500 steps. I commonly end up with a box fluctuation error after 5ns. I tried several parameters and also made a lot of changes in my equilibration protocol (It has been much worse in the beginning, so the protocol now is the best I could get until now) and also switched of the translational restrain because translation would not be a problem. For the equilibration run, I use the standard CHARMM GUI colvars on the backbone and side-chains together with my rotational ones. The rotational restrain is only on the backbone CA atoms of a core region of the peptide (approx. 24 residues in the center). For production I only stay with the rotational constrain (script shown in the end of this message). Changing the integration timestep from 2fs to 1fs didn’t made any difference. From what I can see in VMD (didn’t checked the absolute values), the box immediately starts breathing/oscillating after switching from NVT to NPT and the amplitude of the oscillation increases with time and piston decay and period. Both is not unexpected to me; however, I would expect the oscillation to slow down with ongoing time but I seems to increase.
I don’t have any more ideas what causes this problem. The worst case solution would be to increase the box size to an equally spaced one which would not require any constrain. However, this will significantly blow up (110k atoms vs 28k) and slow down my simulation which currently runs for approx. 110ns/day without fine tuning. I hope you have a tip or solution for me.
Best regards
Alex
#Peptide Translation and Rotation###
Colvarstrajfrequency 100
Colvarsrestartfrequency 100
### ROTATIONAL RESTRAINING ###
colvar {
name rot
orientation {
atoms {
atomsFile peptide.pdb
atomsCol B
atomsColValue 1.0
}
refPositionsFile peptide.pdb
refPositionsCol B
refPositionsColValue 1.0
}
}
harmonic {
colvars rot # Acting on colvar "rot"
centers (1.0, 0.0, 0.0, 0.0) # Unit quaternion (no rotation)
forceConstant 1.0
}
This archive was generated by hypermail 2.1.6 : Mon Dec 31 2018 - 23:20:11 CST