Re: [gmx-users] CG Lincs errors

[Nash, Anthony]

Hi Peter, 

Thanks for the reply, I know we spoke in length on this mater only just
recently. Many thanks for that.

I’ve taken the time step of collagen in vacu down to 0.0001 and I’ve
dropped the temp down to 280. I hope, running over 16 cores for two days
that this should relieve any tension in the backdone before in gradually
increase to 20-40fs. Again, all in vacu with no solvent.

I’ve actually thought about writing a script to modify the equilibrium
bond angles in the CG.itp file for the backbone, using the atomistic
structure as a template. After all, true collagen does not replicate
‘ideal’ collagen along the stretch of the protein e.g., the MMP1 binding
site is not very tightly coiled. Perhaps, there lies my problem.

I haven’t thought too much of the actual full fibril structure, I want to
capture the D-band gap in  type I collagen. When I give it some more
thought I will probably look into semi-isotropic pressure coupling first.

Thanks
Anthony 

Dr Anthony Nash
Department of Chemistry
University College London

Skeletal Tissue Dynamics Group
Committee member of London Matrix Group @LondonMatrixGrp





On 16/12/2016 09:02, "gromacs.org_gmx-users-boun...@maillist.sys.kth.se on
behalf of Peter Kroon"  wrote:

>As a note to Alex (and the rest of the list), the coarse-grained Martini
>forcefield is usually run with timesteps between 20-40 fs. 15fs is
>already rather low. I do agree that longer equilibration at low timestep
>(5 or 10) might help.
>
>Alternatively, Do you think a semiisotropic pressure coupling might be
>applicable in this case, since it's an infinite collagen polymer?
>
>
>Peter (PhD in the Martini group)
>
>
>On 16-12-16 00:21, Nash, Anthony wrote:
>> Alex and Mark, thanks for the information. I’ll drop dt down,
>> significantly, drop the temperature, and run it for a long time.
>>
>> Thanks for the ideas.
>> Anthony 
>>
>>
>> On 15/12/2016 21:54, "gromacs.org_gmx-users-boun...@maillist.sys.kth.se
>>on
>> behalf of Alex" > behalf of nedoma...@gmail.com> wrote:
>>
>>> Mark is right, no two ways about it. For initial equilibration and
>>> assessing preexisting structural strains try vacuum, _much_ smaller
>>> timesteps and possibly low temperatures in vacuum, only then transfer
>>>to
>>> solvent, etc. Algorithmically, LINCS requires convergence and you
>>>already
>>> are using a pretty high LINCS order... From what I see, dt = 15 fs at
>>>310K
>>> looks like a cowboy mode simulation in this case.
>>>
>>> Alex
>>>
>>> On Thu, Dec 15, 2016 at 2:32 PM, Mark Abraham
>>>
>>> wrote:
>>>
 Hi,

 If a simulation isn't stable with a small time step (as I think you
are
 saying) then moving to a larger time step is guaranteed to make that
 worse.
 Try an even smaller time step, for a long time, and see what happens.
Or
 take a subset of your protein and see what happens. Or simulate in
vacuo
 for a while. Your topology could be unsuited to your starting
structure,
 e.g. some part is under a lot of tension that gets released at some
 point
 and no finite time step can in practice deal with the velocity of the
 recoil...

 Mark

 On Thu, 15 Dec 2016 23:06 Nash, Anthony  wrote:

> Hi all,
>
> I¹m hoping for some help. I¹m very sorry, this is a bit of a long
>one.
>
> I¹ve been struggling for almost a month trying to run a CG
 representation
> of our all-atom model of a collagen protein (3 polypeptide chains in
>a
> protein). Our original AMBER all-atom model had been successful
 modelling
> using MD. I went on to use the latest version of Martinize.py with
>the
> latest version of the MARTINI forcefield fields.
>
> After a little tweaking (the way AMBER names histidine residues), I
> successful converted the molecule (approx 3100 amino acids) into a CG
> representation. I successfully energy min the protein in vacuum to a
> threshold of 500, and in solvent to a threshold of 750 using steepest
> descent. Looking for a system at an energy min of a threshold around
 300
 I
> begin to see LINCS warnings. Observing the initial structure, there
>is
> nothing obviously wrong with the bond network (both protein and
 polarized
> CG water).
>
> I take the system that energy mins at 750 (protein-water mix, with no
> fault reported), and went straight to NPT, 20fs step. Blew up. After
>a
 bit
> of chatting with the MARTINI community, I¹ve started with an NVT
 ensemble,
> beginning at 5s then through 10fs, 15fs, and 20fs. I only run for
>1000
> steps before switching. Keeping any of the simulations running for
 longer
> throws lincs warnings followed by a segmentation fault from the
 warning:

Re: [gmx-users] CG Lincs errors

[Peter Kroon]

As a note to Alex (and the rest of the list), the coarse-grained Martini
forcefield is usually run with timesteps between 20-40 fs. 15fs is
already rather low. I do agree that longer equilibration at low timestep
(5 or 10) might help.

Alternatively, Do you think a semiisotropic pressure coupling might be
applicable in this case, since it's an infinite collagen polymer? 


Peter (PhD in the Martini group)


On 16-12-16 00:21, Nash, Anthony wrote:
> Alex and Mark, thanks for the information. I’ll drop dt down,
> significantly, drop the temperature, and run it for a long time.
>
> Thanks for the ideas.
> Anthony 
>
>
> On 15/12/2016 21:54, "gromacs.org_gmx-users-boun...@maillist.sys.kth.se on
> behalf of Alex"  behalf of nedoma...@gmail.com> wrote:
>
>> Mark is right, no two ways about it. For initial equilibration and
>> assessing preexisting structural strains try vacuum, _much_ smaller
>> timesteps and possibly low temperatures in vacuum, only then transfer to
>> solvent, etc. Algorithmically, LINCS requires convergence and you already
>> are using a pretty high LINCS order... From what I see, dt = 15 fs at 310K
>> looks like a cowboy mode simulation in this case.
>>
>> Alex
>>
>> On Thu, Dec 15, 2016 at 2:32 PM, Mark Abraham 
>> wrote:
>>
>>> Hi,
>>>
>>> If a simulation isn't stable with a small time step (as I think you are
>>> saying) then moving to a larger time step is guaranteed to make that
>>> worse.
>>> Try an even smaller time step, for a long time, and see what happens. Or
>>> take a subset of your protein and see what happens. Or simulate in vacuo
>>> for a while. Your topology could be unsuited to your starting structure,
>>> e.g. some part is under a lot of tension that gets released at some
>>> point
>>> and no finite time step can in practice deal with the velocity of the
>>> recoil...
>>>
>>> Mark
>>>
>>> On Thu, 15 Dec 2016 23:06 Nash, Anthony  wrote:
>>>
 Hi all,

 I¹m hoping for some help. I¹m very sorry, this is a bit of a long one.

 I¹ve been struggling for almost a month trying to run a CG
>>> representation
 of our all-atom model of a collagen protein (3 polypeptide chains in a
 protein). Our original AMBER all-atom model had been successful
>>> modelling
 using MD. I went on to use the latest version of Martinize.py with the
 latest version of the MARTINI forcefield fields.

 After a little tweaking (the way AMBER names histidine residues), I
 successful converted the molecule (approx 3100 amino acids) into a CG
 representation. I successfully energy min the protein in vacuum to a
 threshold of 500, and in solvent to a threshold of 750 using steepest
 descent. Looking for a system at an energy min of a threshold around
>>> 300
>>> I
 begin to see LINCS warnings. Observing the initial structure, there is
 nothing obviously wrong with the bond network (both protein and
>>> polarized
 CG water).

 I take the system that energy mins at 750 (protein-water mix, with no
 fault reported), and went straight to NPT, 20fs step. Blew up. After a
>>> bit
 of chatting with the MARTINI community, I¹ve started with an NVT
>>> ensemble,
 beginning at 5s then through 10fs, 15fs, and 20fs. I only run for 1000
 steps before switching. Keeping any of the simulations running for
>>> longer
 throws lincs warnings followed by a segmentation fault from the
>>> warning:
 "3 particles communicated to PME rank 7 are more than 2/3 times the
 cut-off out of the domain decomposition cell of their charge group in
 dimension x."

 Observing the trajectories of any of the extended simulations shows
>>> the
 protein snapping like a rope, and always at the same place. I have
>>> watched
 every trajectory at this point, using numerous energy min start
>>> points,
>>> to
 try and understand why it is blowing up. I can¹t see any obvious
>>> reason.
>>> I
 was told to consider how the temperature is changing. Below is an
>>> example
 of the temperature and pressure from an NPT of 20fs step continued
>>> from
 the very short 20fs step NVT simulation (hoping that perhaps CG
>>> without
 pressure just doesn¹t behave happily; I was wrong).


 TEMP:
 Š
 6.63  311.000336
 6.645000  311.371643
 6.66  311.724213
 6.675000  313.878693
 6.69  681558.937500


 PRESSURE:
 Š
 6.633.559879
 6.6450003.901433
 6.663.589078
 6.6750004.158611
 6.69  81762.437500

 The final LINCS warning from this same run:

 Step 300, time 4.5 (ps)  LINCS WARNING
 relative constraint deviation after LINCS:
 rms 0.35, max 0.003386 (between atoms 2125 and 2126)
 bonds that rotated more than 45 degrees:
  atom 1 atom 2  angle  

Re: [gmx-users] CG Lincs errors

[Nash, Anthony]

Alex and Mark, thanks for the information. I’ll drop dt down,
significantly, drop the temperature, and run it for a long time.

Thanks for the ideas.
Anthony 


On 15/12/2016 21:54, "gromacs.org_gmx-users-boun...@maillist.sys.kth.se on
behalf of Alex"  wrote:

>Mark is right, no two ways about it. For initial equilibration and
>assessing preexisting structural strains try vacuum, _much_ smaller
>timesteps and possibly low temperatures in vacuum, only then transfer to
>solvent, etc. Algorithmically, LINCS requires convergence and you already
>are using a pretty high LINCS order... From what I see, dt = 15 fs at 310K
>looks like a cowboy mode simulation in this case.
>
>Alex
>
>On Thu, Dec 15, 2016 at 2:32 PM, Mark Abraham 
>wrote:
>
>> Hi,
>>
>> If a simulation isn't stable with a small time step (as I think you are
>> saying) then moving to a larger time step is guaranteed to make that
>>worse.
>> Try an even smaller time step, for a long time, and see what happens. Or
>> take a subset of your protein and see what happens. Or simulate in vacuo
>> for a while. Your topology could be unsuited to your starting structure,
>> e.g. some part is under a lot of tension that gets released at some
>>point
>> and no finite time step can in practice deal with the velocity of the
>> recoil...
>>
>> Mark
>>
>> On Thu, 15 Dec 2016 23:06 Nash, Anthony  wrote:
>>
>> > Hi all,
>> >
>> > I¹m hoping for some help. I¹m very sorry, this is a bit of a long one.
>> >
>> > I¹ve been struggling for almost a month trying to run a CG
>>representation
>> > of our all-atom model of a collagen protein (3 polypeptide chains in a
>> > protein). Our original AMBER all-atom model had been successful
>>modelling
>> > using MD. I went on to use the latest version of Martinize.py with the
>> > latest version of the MARTINI forcefield fields.
>> >
>> > After a little tweaking (the way AMBER names histidine residues), I
>> > successful converted the molecule (approx 3100 amino acids) into a CG
>> > representation. I successfully energy min the protein in vacuum to a
>> > threshold of 500, and in solvent to a threshold of 750 using steepest
>> > descent. Looking for a system at an energy min of a threshold around
>>300
>> I
>> > begin to see LINCS warnings. Observing the initial structure, there is
>> > nothing obviously wrong with the bond network (both protein and
>>polarized
>> > CG water).
>> >
>> > I take the system that energy mins at 750 (protein-water mix, with no
>> > fault reported), and went straight to NPT, 20fs step. Blew up. After a
>> bit
>> > of chatting with the MARTINI community, I¹ve started with an NVT
>> ensemble,
>> > beginning at 5s then through 10fs, 15fs, and 20fs. I only run for 1000
>> > steps before switching. Keeping any of the simulations running for
>>longer
>> > throws lincs warnings followed by a segmentation fault from the
>>warning:
>> >
>> > "3 particles communicated to PME rank 7 are more than 2/3 times the
>> > cut-off out of the domain decomposition cell of their charge group in
>> > dimension x."
>> >
>> > Observing the trajectories of any of the extended simulations shows
>>the
>> > protein snapping like a rope, and always at the same place. I have
>> watched
>> > every trajectory at this point, using numerous energy min start
>>points,
>> to
>> > try and understand why it is blowing up. I can¹t see any obvious
>>reason.
>> I
>> > was told to consider how the temperature is changing. Below is an
>>example
>> > of the temperature and pressure from an NPT of 20fs step continued
>>from
>> > the very short 20fs step NVT simulation (hoping that perhaps CG
>>without
>> > pressure just doesn¹t behave happily; I was wrong).
>> >
>> >
>> > TEMP:
>> > Š
>> > 6.63  311.000336
>> > 6.645000  311.371643
>> > 6.66  311.724213
>> > 6.675000  313.878693
>> > 6.69  681558.937500
>> >
>> >
>> > PRESSURE:
>> > Š
>> > 6.633.559879
>> > 6.6450003.901433
>> > 6.663.589078
>> > 6.6750004.158611
>> > 6.69  81762.437500
>> >
>> > The final LINCS warning from this same run:
>> >
>> > Step 300, time 4.5 (ps)  LINCS WARNING
>> > relative constraint deviation after LINCS:
>> > rms 0.35, max 0.003386 (between atoms 2125 and 2126)
>> > bonds that rotated more than 45 degrees:
>> >  atom 1 atom 2  angle  previous, current, constraint length
>> >2125   2126   68.30.2781   0.2691  0.2700
>> >2125   2127   45.90.2789   0.2701  0.2700
>> >
>> >
>> > At this stage the structure ruptures as described above.
>> >
>> >
>> > My NVT settings (with NPT included to save space) are:
>> >
>> > -
>> > title= Martini
>> >
>> > integrator   = md
>> > dt   = 0.015
>> > nsteps   = 1000
>> > nstcomm  = 100
>> > ;comm-grps  

Re: [gmx-users] CG Lincs errors

[Alex]

Mark is right, no two ways about it. For initial equilibration and
assessing preexisting structural strains try vacuum, _much_ smaller
timesteps and possibly low temperatures in vacuum, only then transfer to
solvent, etc. Algorithmically, LINCS requires convergence and you already
are using a pretty high LINCS order... From what I see, dt = 15 fs at 310K
looks like a cowboy mode simulation in this case.

Alex

On Thu, Dec 15, 2016 at 2:32 PM, Mark Abraham 
wrote:

> Hi,
>
> If a simulation isn't stable with a small time step (as I think you are
> saying) then moving to a larger time step is guaranteed to make that worse.
> Try an even smaller time step, for a long time, and see what happens. Or
> take a subset of your protein and see what happens. Or simulate in vacuo
> for a while. Your topology could be unsuited to your starting structure,
> e.g. some part is under a lot of tension that gets released at some point
> and no finite time step can in practice deal with the velocity of the
> recoil...
>
> Mark
>
> On Thu, 15 Dec 2016 23:06 Nash, Anthony  wrote:
>
> > Hi all,
> >
> > I¹m hoping for some help. I¹m very sorry, this is a bit of a long one.
> >
> > I¹ve been struggling for almost a month trying to run a CG representation
> > of our all-atom model of a collagen protein (3 polypeptide chains in a
> > protein). Our original AMBER all-atom model had been successful modelling
> > using MD. I went on to use the latest version of Martinize.py with the
> > latest version of the MARTINI forcefield fields.
> >
> > After a little tweaking (the way AMBER names histidine residues), I
> > successful converted the molecule (approx 3100 amino acids) into a CG
> > representation. I successfully energy min the protein in vacuum to a
> > threshold of 500, and in solvent to a threshold of 750 using steepest
> > descent. Looking for a system at an energy min of a threshold around 300
> I
> > begin to see LINCS warnings. Observing the initial structure, there is
> > nothing obviously wrong with the bond network (both protein and polarized
> > CG water).
> >
> > I take the system that energy mins at 750 (protein-water mix, with no
> > fault reported), and went straight to NPT, 20fs step. Blew up. After a
> bit
> > of chatting with the MARTINI community, I¹ve started with an NVT
> ensemble,
> > beginning at 5s then through 10fs, 15fs, and 20fs. I only run for 1000
> > steps before switching. Keeping any of the simulations running for longer
> > throws lincs warnings followed by a segmentation fault from the warning:
> >
> > "3 particles communicated to PME rank 7 are more than 2/3 times the
> > cut-off out of the domain decomposition cell of their charge group in
> > dimension x."
> >
> > Observing the trajectories of any of the extended simulations shows the
> > protein snapping like a rope, and always at the same place. I have
> watched
> > every trajectory at this point, using numerous energy min start points,
> to
> > try and understand why it is blowing up. I can¹t see any obvious reason.
> I
> > was told to consider how the temperature is changing. Below is an example
> > of the temperature and pressure from an NPT of 20fs step continued from
> > the very short 20fs step NVT simulation (hoping that perhaps CG without
> > pressure just doesn¹t behave happily; I was wrong).
> >
> >
> > TEMP:
> > Š
> > 6.63  311.000336
> > 6.645000  311.371643
> > 6.66  311.724213
> > 6.675000  313.878693
> > 6.69  681558.937500
> >
> >
> > PRESSURE:
> > Š
> > 6.633.559879
> > 6.6450003.901433
> > 6.663.589078
> > 6.6750004.158611
> > 6.69  81762.437500
> >
> > The final LINCS warning from this same run:
> >
> > Step 300, time 4.5 (ps)  LINCS WARNING
> > relative constraint deviation after LINCS:
> > rms 0.35, max 0.003386 (between atoms 2125 and 2126)
> > bonds that rotated more than 45 degrees:
> >  atom 1 atom 2  angle  previous, current, constraint length
> >2125   2126   68.30.2781   0.2691  0.2700
> >2125   2127   45.90.2789   0.2701  0.2700
> >
> >
> > At this stage the structure ruptures as described above.
> >
> >
> > My NVT settings (with NPT included to save space) are:
> >
> > -
> > title= Martini
> >
> > integrator   = md
> > dt   = 0.015
> > nsteps   = 1000
> > nstcomm  = 100
> > ;comm-grps   =
> >
> > nstxout  = 1000
> > nstvout  = 1000
> > nstfout  = 0
> > nstlog   = 1
> > nstenergy= 1
> > nstxout-compressed   = 0
> > compressed-x-precision   = 0
> > ;compressed-x-grps=
> > energygrps   = collagen solvent
> >
> > cutoff-scheme= Verlet
> > nstlist  = 20
> > ns_type  = grid
> >
> > pbc  = xyz
> > 

Re: [gmx-users] CG Lincs errors

[Mark Abraham]

Hi,

If a simulation isn't stable with a small time step (as I think you are
saying) then moving to a larger time step is guaranteed to make that worse.
Try an even smaller time step, for a long time, and see what happens. Or
take a subset of your protein and see what happens. Or simulate in vacuo
for a while. Your topology could be unsuited to your starting structure,
e.g. some part is under a lot of tension that gets released at some point
and no finite time step can in practice deal with the velocity of the
recoil...

Mark

On Thu, 15 Dec 2016 23:06 Nash, Anthony  wrote:

> Hi all,
>
> I¹m hoping for some help. I¹m very sorry, this is a bit of a long one.
>
> I¹ve been struggling for almost a month trying to run a CG representation
> of our all-atom model of a collagen protein (3 polypeptide chains in a
> protein). Our original AMBER all-atom model had been successful modelling
> using MD. I went on to use the latest version of Martinize.py with the
> latest version of the MARTINI forcefield fields.
>
> After a little tweaking (the way AMBER names histidine residues), I
> successful converted the molecule (approx 3100 amino acids) into a CG
> representation. I successfully energy min the protein in vacuum to a
> threshold of 500, and in solvent to a threshold of 750 using steepest
> descent. Looking for a system at an energy min of a threshold around 300 I
> begin to see LINCS warnings. Observing the initial structure, there is
> nothing obviously wrong with the bond network (both protein and polarized
> CG water).
>
> I take the system that energy mins at 750 (protein-water mix, with no
> fault reported), and went straight to NPT, 20fs step. Blew up. After a bit
> of chatting with the MARTINI community, I¹ve started with an NVT ensemble,
> beginning at 5s then through 10fs, 15fs, and 20fs. I only run for 1000
> steps before switching. Keeping any of the simulations running for longer
> throws lincs warnings followed by a segmentation fault from the warning:
>
> "3 particles communicated to PME rank 7 are more than 2/3 times the
> cut-off out of the domain decomposition cell of their charge group in
> dimension x."
>
> Observing the trajectories of any of the extended simulations shows the
> protein snapping like a rope, and always at the same place. I have watched
> every trajectory at this point, using numerous energy min start points, to
> try and understand why it is blowing up. I can¹t see any obvious reason. I
> was told to consider how the temperature is changing. Below is an example
> of the temperature and pressure from an NPT of 20fs step continued from
> the very short 20fs step NVT simulation (hoping that perhaps CG without
> pressure just doesn¹t behave happily; I was wrong).
>
>
> TEMP:
> Š
> 6.63  311.000336
> 6.645000  311.371643
> 6.66  311.724213
> 6.675000  313.878693
> 6.69  681558.937500
>
>
> PRESSURE:
> Š
> 6.633.559879
> 6.6450003.901433
> 6.663.589078
> 6.6750004.158611
> 6.69  81762.437500
>
> The final LINCS warning from this same run:
>
> Step 300, time 4.5 (ps)  LINCS WARNING
> relative constraint deviation after LINCS:
> rms 0.35, max 0.003386 (between atoms 2125 and 2126)
> bonds that rotated more than 45 degrees:
>  atom 1 atom 2  angle  previous, current, constraint length
>2125   2126   68.30.2781   0.2691  0.2700
>2125   2127   45.90.2789   0.2701  0.2700
>
>
> At this stage the structure ruptures as described above.
>
>
> My NVT settings (with NPT included to save space) are:
>
> -
> title= Martini
>
> integrator   = md
> dt   = 0.015
> nsteps   = 1000
> nstcomm  = 100
> ;comm-grps   =
>
> nstxout  = 1000
> nstvout  = 1000
> nstfout  = 0
> nstlog   = 1
> nstenergy= 1
> nstxout-compressed   = 0
> compressed-x-precision   = 0
> ;compressed-x-grps=
> energygrps   = collagen solvent
>
> cutoff-scheme= Verlet
> nstlist  = 20
> ns_type  = grid
>
> pbc  = xyz
> verlet-buffer-tolerance  = 0.005
>
> coulombtype  = PME ;reaction-field
> rcoulomb = 1.1
> fourierspacing   = 0.16 ;0.2  ;0.12
>
> epsilon_r= 2.5 ;15  ; 2.5 (with polarizable water)
> epsilon_rf   = 0
> vdw_type = cutoff
> vdw-modifier = Potential-shift-verlet
> rvdw = 1.1
>
> tcoupl   = v-rescale ;berendsen ;v-rescale
> tc-grps  = collagen solvent
> tau_t= 0.5 0.5 ;1.0 1.0
> ref_t= 310 310
>
> Pcoupl   = berendsen   ;parrinello-rahman
> Pcoupltype   = isotropic
> tau_p= 12.0  ; 

[gmx-users] CG Lincs errors

[Nash, Anthony]

Hi all,

I¹m hoping for some help. I¹m very sorry, this is a bit of a long one.

I¹ve been struggling for almost a month trying to run a CG representation
of our all-atom model of a collagen protein (3 polypeptide chains in a
protein). Our original AMBER all-atom model had been successful modelling
using MD. I went on to use the latest version of Martinize.py with the
latest version of the MARTINI forcefield fields.

After a little tweaking (the way AMBER names histidine residues), I
successful converted the molecule (approx 3100 amino acids) into a CG
representation. I successfully energy min the protein in vacuum to a
threshold of 500, and in solvent to a threshold of 750 using steepest
descent. Looking for a system at an energy min of a threshold around 300 I
begin to see LINCS warnings. Observing the initial structure, there is
nothing obviously wrong with the bond network (both protein and polarized
CG water).

I take the system that energy mins at 750 (protein-water mix, with no
fault reported), and went straight to NPT, 20fs step. Blew up. After a bit
of chatting with the MARTINI community, I¹ve started with an NVT ensemble,
beginning at 5s then through 10fs, 15fs, and 20fs. I only run for 1000
steps before switching. Keeping any of the simulations running for longer
throws lincs warnings followed by a segmentation fault from the warning:

"3 particles communicated to PME rank 7 are more than 2/3 times the
cut-off out of the domain decomposition cell of their charge group in
dimension x."

Observing the trajectories of any of the extended simulations shows the
protein snapping like a rope, and always at the same place. I have watched
every trajectory at this point, using numerous energy min start points, to
try and understand why it is blowing up. I can¹t see any obvious reason. I
was told to consider how the temperature is changing. Below is an example
of the temperature and pressure from an NPT of 20fs step continued from
the very short 20fs step NVT simulation (hoping that perhaps CG without
pressure just doesn¹t behave happily; I was wrong).


TEMP:
Š
6.63  311.000336
6.645000  311.371643
6.66  311.724213
6.675000  313.878693
6.69  681558.937500


PRESSURE:
Š
6.633.559879
6.6450003.901433
6.663.589078
6.6750004.158611
6.69  81762.437500

The final LINCS warning from this same run:

Step 300, time 4.5 (ps)  LINCS WARNING
relative constraint deviation after LINCS:
rms 0.35, max 0.003386 (between atoms 2125 and 2126)
bonds that rotated more than 45 degrees:
 atom 1 atom 2  angle  previous, current, constraint length
   2125   2126   68.30.2781   0.2691  0.2700
   2125   2127   45.90.2789   0.2701  0.2700


At this stage the structure ruptures as described above.


My NVT settings (with NPT included to save space) are:

-
title= Martini

integrator   = md
dt   = 0.015
nsteps   = 1000
nstcomm  = 100
;comm-grps   =

nstxout  = 1000
nstvout  = 1000
nstfout  = 0
nstlog   = 1
nstenergy= 1
nstxout-compressed   = 0
compressed-x-precision   = 0
;compressed-x-grps=
energygrps   = collagen solvent

cutoff-scheme= Verlet
nstlist  = 20
ns_type  = grid

pbc  = xyz
verlet-buffer-tolerance  = 0.005

coulombtype  = PME ;reaction-field
rcoulomb = 1.1
fourierspacing   = 0.16 ;0.2  ;0.12

epsilon_r= 2.5 ;15  ; 2.5 (with polarizable water)
epsilon_rf   = 0
vdw_type = cutoff
vdw-modifier = Potential-shift-verlet
rvdw = 1.1

tcoupl   = v-rescale ;berendsen ;v-rescale
tc-grps  = collagen solvent
tau_t= 0.5 0.5 ;1.0 1.0
ref_t= 310 310

Pcoupl   = berendsen   ;parrinello-rahman
Pcoupltype   = isotropic
tau_p= 12.0  ; parrinello-rahman is more stable with
larger tau-p, DdJ, 20130422
compressibility  = 10e-4
ref_p= 1.0
refcoord_scaling = com

gen_vel  = no
gen_temp = 310
gen_seed = 473529

continuation = yes
constraints  = none
constraint_algorithm = lincs
lincs-warnangle = 45
lincs-order=8
lincs-iter=4


‹‹

Every setting bar the lincs iter, order, warnangle were supplied with the
latest version of MARTINI. During many NVT runs I have adjusted the tau-t
to try and keep the thermostat from oscillating its way into infinity.

I¹m curious, will an out of control thermostat break a structure, or will
a structure breaking (for what ever reason this structure is breaking)
cause the thermostat to go out of control?

My only