[AMBER] Thermodynamic Integration problems

From: Prokopios Andrikopoulos via AMBER <amber.ambermd.org>
Date: Mon, 24 Jul 2023 16:27:17 +0000

Dear AMBER Users,
I've been trying to do a TI for a mutation (res position 362) on a quite large trimer (1125 res per chain with a ligand at each chain at positions 3376-3378 - the ligand associated to the mutation is at 3378). I'm using ff14sb with TIP3P via the Amber 22 pmemd.cuda and the parmed from Amber 18 (v22 parmed does not work on my cluster installation atm). I've adapted the relaxation regime from tutorial 3.1, by changing ntmin=2, ntf=1 and ntp=1 to make it acceptable for TI. I do a minimization > slow heating > and a series of 1ns constant pressure runs (with one more minimization in-between) releasing slowly the restrains as in tutorial 3.1. The relaxation regime works without any problems without the TI keywords (icfe, clambda etc.) on the normal prmtop using either a box or truncated octahedron for solvation. When I prepare the dual topology parm+rst7 is where the problems start. My parmed command is: tiMerge :1-3377 :3378-6754 :362 :3739. Parmed moves the mutated res in the merged files right after the 1st chain, so in the relaxation inputs I have to point to 1126 instead, with the associated ligand at 3379.

While the normal relaxation w/out TI reports a normal Density of ~1.00, the TI calculations density is reported in the 0.12-0.25 range in the constant pressure runs. A difference in energy is also evident in the initial minimization between the normal system and the TI complex with ë = 0.0: i.e. after 1000 min steps, this is EAMBER = -2180402.3012 and -2171818.4624 for the normal and TI calculation, respectively.
Below are inputs for the complex with ë = 0.0. The input for the constant pressure run is the 4th in the regime with reduced backbone restrains. This fails for all ë values at 65-%75% of the total steps with the ERROR: Calculation halted. Periodic box dimensions have changed too much from their initial values. I would appreciate any help with this problem, and I can provide more details for the calculations.

A secondary problem I have encountered is when I use a truncated octahedron for the solvation of the system with solvateOct. I use this command in leap for the complex and the protein: setBox protein/complex vdw { x = y = z }, and then I manually edit the last line of protein/complex.rst7 and merged_protein/complex.rst7 with the angles a = b = c = 109.47. When I run the calculation, it halts with Error: ifbox=2 in prmtop but angles are not correct. Incidentally, when I set ifbox=3 it runs without this error. Is this an acceptable solution?

minimisation
 &cntrl
   imin = 1, ntmin = 2, maxcyc = 10000,
   ntpr = 50, ntwe = 0, ntwr = 500,
   ntb = 1, ntp = 0, ntf = 1, ntc = 2,
   ntx = 1, ntr = 1, ioutfm = 1, ntxo = 2,
   cut = 10.0,
   restraintmask = ':1-3379',
   restraint_wt = 100.,


   icfe = 1, ifsc = 1, clambda = 0.0, scalpha = 0.5, scbeta = 12.0,
   logdvdl = 0,
   timask1 = ':362', timask2 = ':1126',
   scmask1 = ':362', scmask2 = ':1126',
 /

 &ewald
 /


heating to 298K
 &cntrl
   imin = 0, nstlim = 1000000, irest = 0, ntx = 1, dt = 0.001,
   ntt = 3, temp0 = 298.0, tempi = 100.0, tautp = 1.0, ig = -1,
   ntc = 2, ntf = 1,
   ntb = 1, ntp = 0, tol = 0.00001,
   ioutfm = 1, iwrap = 0, ntxo = 2, nscm = 0, cut = 8.0,
   ntwe = 0, ntwx = 10000, ntpr = 1000, ntwr = 1000,

   nmropt = 1, gamma_ln = 1.,
   ntr = 1, restraint_wt = 100.,
   restraintmask=':1-3379',

   icfe = 1, ifsc = 1, clambda = 0.0, scalpha = 0.5, scbeta = 12.0,
   logdvdl = 0,
   timask1 = ':362', timask2 = ':1126',
   scmask1 = ':362', scmask2 = ':1126',
 /

 &ewald
 /

 &wt
   type='TEMP0',
   istep1 = 0, istep2 = 1000000,
   value1 = 100., value2 = 298.
 /

 &wt type = 'END'
 /

4th constant pressure run of 1ns
&cntrl
   imin = 0, nstlim = 1000000, dt = 0.001,
   irest = 1, ntx = 5, ig = -1,
   temp0 = 298.0,
   ntc = 2, ntf = 1,
   ntwx = 10000, ntwe = 0, ntwr = 1000, ntpr = 1000,
   cut = 8.0, iwrap = 0,
   ntt = 3, gamma_ln = 1.0, ntb = 2, ntp = 1, barostat = 2,
   nscm = 0,
   ntr = 1, restraintmask = "@CA,N,C", restraint_wt = 1.
   ioutfm = 1, ntxo = 2,

   icfe = 1, ifsc = 1, clambda = 0.0, scalpha = 0.5, scbeta = 12.0,
   logdvdl = 0,
   timask1 = ':362', timask2 = ':1126',
   scmask1 = ':362', scmask2 = ':1126',
 /

 &ewald
 /


Thanks for any suggestions
P. A.

_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Mon Jul 24 2023 - 10:00:03 PDT
Custom Search