Hi Gustaf, thanks for the files.
So I haven't run anything yet but I think I've spotted an issue right
off the bat. In your original leap script you have the following:
```
source leaprc.water.opc
source leaprc.protein.ff19SB
WAT = OPC
```
This is problematic because it allows the 19SB frcmod to overwrite
parameters set by the OPC frcmod. The manual gives an example of
loading different parameter files at the beginning of section 3
(Molecular mechanics force fields) that gives a typical load order,
but I don't think it emphasizes how much the order can matter - that's
something that needs to be fixed.
Anyway, if I reorder the parameter files so the protein FF comes first
and the solvent model last, like so:
```
source leaprc.protein.ff19SB
source leaprc.water.opc
go = loadpdb ../ff19SB_opc/waterPeptBox.pdb
setbox go vdw
saveamberparm go test.waterPeptBox-opc.prmtop test.waterPeptBox-opc.inpcrd
quit
```
And then compare the two topologies using the 'comparetop' command from cpptraj:
```
parm waterPeptBox-opc.prmtop
parm test.waterPeptBox-opc.prmtop
comparetop parmindex 0 parmindex 1 out parm.diff
```
I can see that in my "test" topology (denoted by '>') the correct OPC
water parameters are assigned compared to your original topology, e.g.
```
# LJ params
< OW Rmin= 1.7683 Eps= 0.152
> OW Rmin= 1.77717 Eps= 0.212801
```
I sent the full diff and topology/coordintae files to you off-list.
Because the parameters assigned are different, the output coordinates
end up being a little different as well. I still need to run tests on
the "new" topology, but I suspect this may be the problem right here.
Let me know if this makes sense and/or if you have more questions.
-Dan
On Fri, Oct 23, 2020 at 7:28 AM Daniel Roe <daniel.r.roe.gmail.com> wrote:
>
> If you’d like, send me off-list your topology and the problematic coordinates and I can check it out and try to reproduce the bad behavior, and hopefully smoke out any bugs.
>
> -Dan
>
> On Fri, Oct 23, 2020 at 7:13 AM Gustaf Olsson <gustaf.olsson.lnu.se> wrote:
>>
>> Thank you very much for the tips, I will look at them carefully as soon as I have time. I am sure I will have use for these at some point.
>>
>> Just to avoid confusion. This is a very small system with an initially linear, capped hexapeptide. I am doing energy minimisation, temperature equilibration and pressure equilibration with no restraints. I am aware of the basic idea for proteins/membranes, alternating restraining the movement of some structures while letting others move and so forth. My reservation is that this only system seems to have everything constrained though no restraints have been used whatsoever. The other systems using ff14SB with SPCE/TIP3P and the plain OPC water system are equilibrated using identical input so unless there is some input option not compatible with ff19SB, then the protocol should not be the cause of the observed behaviour.
>>
>> If the system did not equilibrate properly, that would be “fine”. Though I am getting nothing in the “process_m*.perl” scripts output regarding energy minimisation or pressure equilibration. Visually, there is no “diffusion”/spreading of waters. Everything seems to be “locked” in place, no water is displaying any expected degree of rotation or movement. The system is basically incompressible.
>>
>> I could try ff14SB-OPC and ff19SB-tip3p/spce to see if anything changes though was hoping that someone would notice some obvious mistake.
>>
>> Thank you again for suggestions and tips, it is greatly appreciated.
>>
>> I’ll keep trying
>> Best regards
>> // Gustaf
>>
>>
>> > On 23 Oct 2020, at 11:36, Daniel Roe <daniel.r.roe.gmail.com> wrote:
>> >
>> > Hi Gustaf,
>> >
>> > I’m not sure what equilibration protocol you’re using, but the one I use
>> > generally does the trick, detailed here:
>> > https://aip.scitation.org/doi/abs/10.1063/5.0013849
>> >
>> > I also have a script that implements it for amber (it’s a bit rough around
>> > the edges so use with care):
>> > https://github.com/drroe/AmberMdPrep
>> >
>> > Hopefully it is helpful.
>> >
>> > -Dan
>> >
>> > On Fri, Oct 23, 2020 at 5:26 AM Gustaf Olsson <gustaf.olsson.lnu.se> wrote:
>> >
>> >> I am experiencing some other problem, it seems that the OPC water
>> >> simulation works as intended though a box with a capped hexapeptide and a
>> >> Cl- ion does not.
>> >>
>> >> I have run 6 short simulations 3 with only water (OPC, SPCE and TIP3P) as
>> >> well as 3 with the same water models as well as a hexapeptide and a
>> >> chloride ion (FF14SP-SPCE, FF14SB-TIP3P and FF19SB-OPC).
>> >>
>> >> All systems with only water equilibrates nicely inside the timeframe. The
>> >> pressure equilibration of OPC water is finished faster then the SPCE and
>> >> TIP3P models and all simulation produce usable plots with the
>> >> “process_m*.perl” scripts. Everything looks fine.
>> >>
>> >> The FF19SB-OPC system does NOT behave as intended. I can keep hammering
>> >> this system with higher pressure and shorter taup though looking at the
>> >> simulation, it seems like I am trying to equilibrate a block of ice with
>> >> air trapped inside: https://gifyu.com/image/865p (Do note that I have NOT
>> >> autoimaged with CPPTRAJ, the restricted motion of the water prevents
>> >> changing the shape).
>> >>
>> >> Running only OPC water works just fine, so it seems that it is the FF19SB
>> >> part which does not work as intended. For some reason, the entire system
>> >> appears to be suffering from some motion restriction though I have
>> >> implemented none. I have tried to see if the manual contained any
>> >> information regarding particular options for input files or similar but did
>> >> not find anything.
>> >>
>> >> I used the same input for all simulations:
>> >>
>> >> Energy minimization, 5000 steps SD then 5000 CG
>> >> &cntrl
>> >> imin=1,
>> >> maxcyc=10000,
>> >> ncyc=5000,
>> >> ntb=1,
>> >> cut=9.0,
>> >> /
>> >>
>> >> 125 ps to constant T
>> >> &cntrl
>> >> ig =-1,
>> >> imin = 0,
>> >> iwrap = 0,
>> >> irest = 0,
>> >> ntx = 1,
>> >> ntb = 1,
>> >> cut = 9.0,
>> >> ntc = 2,
>> >> ntf = 2,
>> >> tempi = 0.0,
>> >> temp0 = 293.15,
>> >> ntt = 3,
>> >> gamma_ln = 1.0,
>> >> nstlim = 62500, dt = 0.002,
>> >> ntpr = 1000, ntwx = 1000, ntwr = 2000
>> >> /
>> >>
>> >> Constant T, pressure eq: 500 ps MD
>> >> &cntrl
>> >> ig =-1,
>> >> imin = 0,
>> >> iwrap = 0,
>> >> irest = 1,
>> >> ntx = 5,
>> >> ntb = 2,
>> >> pres0 = 1,
>> >> ntp = 1,
>> >> barostat = 1,
>> >> taup = 2.0,
>> >> cut = 9.0,
>> >> ntc = 2,
>> >> ntf = 2,
>> >> tempi = 293.15,
>> >> temp0 = 293.15,
>> >> ntt = 3,
>> >> gamma_ln = 1.0,
>> >> nstlim = 250000, dt = 0.002,
>> >> ntpr = 1000, ntwx = 1000, ntwr = 2000
>> >> /
>> >>
>> >> All the intended output files are produced and seemingly (except for the
>> >> restricted motion) looks normal and reads just fine in cpptraj and vmd.
>> >>
>> >> Any ideas regarding what I am doing wrong?
>> >>
>> >> Best regards
>> >> // Gustaf
>> >>
>> >>
>> >>> On 22 Oct 2020, at 15:06, Carlos Simmerling <carlos.simmerling.gmail.com>
>> >> wrote:
>> >>>
>> >>> Dave - no I haven't tried that during equilibration, since tighter
>> >> coupling
>> >>> fixes it I just made it part of our equilibration protocol long ago and
>> >>> moved on. I'll try that since I have a sample setup with a bubble and see
>> >>> if it works better.
>> >>> carlos
>> >>>
>> >>> On Thu, Oct 22, 2020 at 8:48 AM David A Case <david.case.rutgers.edu>
>> >> wrote:
>> >>>
>> >>>> On Wed, Oct 21, 2020, Carlos Simmerling wrote:
>> >>>>
>> >>>>> what is your pressure coupling constant? I often get vacuum bubbles
>> >> unless
>> >>>>> I include a step in my equilibration with tight coupling, taup around
>> >> 0.1
>> >>>>> or smaller.
>> >>>>
>> >>>> Carlos: have you ever tried barostat=2? It would be interesting to know
>> >>>> if that makes a difference. Also, what is your starting density? If
>> >>>> you use a closeness parameter of 0.7--0.8 (in solvateBox or solvateOct),
>> >>>> you'll get a higher starting density, which might(?) reduce the tendency
>> >>>> to get bubbles.
>> >>>>
>> >>>> ....dac
>> >>>>
>> >>>>
>> >>>> _______________________________________________
>> >>>> AMBER mailing list
>> >>>> AMBER.ambermd.org
>> >>>> http://lists.ambermd.org/mailman/listinfo/amber
>> >>>>
>> >>> _______________________________________________
>> >>> AMBER mailing list
>> >>> AMBER.ambermd.org
>> >>> http://lists.ambermd.org/mailman/listinfo/amber
>> >>
>> >> _______________________________________________
>> >> AMBER mailing list
>> >> AMBER.ambermd.org
>> >> http://lists.ambermd.org/mailman/listinfo/amber
>> >>
>> > _______________________________________________
>> > AMBER mailing list
>> > AMBER.ambermd.org
>> > http://lists.ambermd.org/mailman/listinfo/amber
>>
>> _______________________________________________
>> AMBER mailing list
>> AMBER.ambermd.org
>> http://lists.ambermd.org/mailman/listinfo/amber
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Fri Oct 23 2020 - 08:00:02 PDT
