[AMBER] Problem with distance restraints in membrane system with pmemd.cuda

From: David Bickel <david.bickel.uni-duesseldorf.de>
Date: Wed, 6 Mar 2019 10:18:40 +0100

Dear Amber-team,

I encountered a problem running simulations of a membrane system with
distance restraints on GPUs. My system consists of a natural product
molecule (750 Da) in the center of a membrane. The goal is to perform a
steered MD, pulling the ligand from the membrane into the solvent.

The system was set up using packmol-memgen from AmberTools.

For equilibration, I first performed a minimization, a short NVT MD for
heating, and a NPT MD for the density equilibration on CPUs using pmemd.MPI.
For the MD steps, I used distance restraints (fxyz=0,0,1) between the
COM of the molecule (:1) and the phosphor atoms of the lipids (@P31), in
order to keep the molecule in the center of the membrane.

Now I wanted to switch to pmemd.cuda to equilibrate the system for
further 50 ns under NTP conditions, using the same distance restraints
as above. However, when running the system on GPUs the RESTRAINT energy
immediately overflows (within 1 step), while the corresponding
DUMPAVE-file still shows reasonable z-values. Eventually such a system
would blow up and cancel.

Starting the same calculation with pmemd.MPI works fine and returns
output as expected.

I tried equilibrating the system for further 10 ns on CPUs, before
switching to GPUs, but the issue remains the same.
I tried using different parameters for the constant pressure-condition
as well as disabling iwrap. No success either.
I also tried running the simulation on different clusters with different
Amber-builds (Amber18 & 17) under different UNIX-evironments. In any
case the result was the same, pmemd.cuda does not work while pmemd.MPI does.
Assuming, that the error is in my restraint definition, I also tried
defining the distance restraint towards the phosphor atoms of one
leaflet only. Still no, success.
By chance, I realized, that modifying the ligand atom selection in my
restraint file to the first 8 atoms (IGR1(1-8)) gets the system running
on GPUs. The initial reaction coordinate is ~4 A then and thus, the
restraint energy starts at ~80 kcal/mol. However, the system is running
just fine and the ligand gradually reorients to better fulfill the
restraints.

Does anyone have an idea, where the issue might be located? I would be
very grateful for any suggestions.

Below, I'll append all input/output files mentioned above.

Best regards,
David Bickel

#
# -- MD INPUT FILE
#
   > NPT MD for 9,000 ps for equilibration
   >  &cntrl
   >   ntx = 5,
   >   irest = 1, iwrap = 0,
   >
   >   ntpr = 1,
   >   ntwx = 5,
   >
   >   nstlim = 4500000,
   >   t = 1000.0,
   >   dt = 0.002,
   >
   >   ntt = 3,
   >   gamma_ln = 2.0,
   >   temp0 = 300.0, tempi = 300.0,
   >
   >   ntb = 2,
   >   barostat = 1, pres0 = 1.0, comp = 44.6, taup = 1.0,
   >   ntp = 3, csurften = 3, gamma_ten = 0.0, ninterface = 2,
   >
   >   ntc = 2, ntf = 2,
   >   tol = 0.000001,
   >
   >   nmropt = 1,
   >  /
   >  &wt type='DUMPFREQ',
   >   istep1=1,
   >  /
   >  &wt type='END',
   >  /
   >  DISANG=pxa_00_5kcal.rst
   >  DUMPAVE=06_ntp_equi.rstout
   >  LISTIN=POUT
   >  LISTOUT=POUT


#
# --- RESTRAINT FILE "pxa_00_5kcal.rst"
#
   >  &rst iat=-1,-1,0
   >    r1=-50.000000, r2=0.000000, r3=0.000000, r4=50.000000,
rk2=5.000000, rk3=5.000000,
   >    iresid=0,
   >    fxyz=0,0,1,
   >    outxyz=1,
   >
IGR1(1)=1,IGR1(2)=2,IGR1(3)=3,IGR1(4)=4,IGR1(5)=5,IGR1(6)=6,IGR1(7)=7,IGR1(8)=8,IGR1(9)=9,IGR1(10)=10,IGR1(11)=11,IGR1(12)=12,IGR1(13)=13,IGR1(14)=14,IGR1(15)=15,IGR1(16)=16,IGR1(17)=17,IGR1(18)=18,IGR1(19)=19,IGR1(20)=20,IGR1(21)=21,IGR1(22)=22,IGR1(23)=23,IGR1(24)=24,IGR1(25)=25,IGR1(26)=26,IGR1(27)=27,IGR1(28)=28,IGR1(29)=29,IGR1(30)=30,IGR1(31)=31,IGR1(32)=32,IGR1(33)=33,IGR1(34)=34,IGR1(35)=35,IGR1(36)=36,IGR1(37)=37,IGR1(38)=38,IGR1(39)=39,IGR1(40)=40,IGR1(41)=41,IGR1(42)=42,IGR1(43)=43,IGR1(44)=44,IGR1(45)=45,IGR1(46)=46,IGR1(47)=47,IGR1(48)=48,IGR1(49)=49,IGR1(50)=50,IGR1(51)=51,IGR1(52)=52,IGR1(53)=53,IGR1(54)=54,IGR1(55)=55,IGR1(56)=56,IGR1(57)=57,IGR1(58)=58,IGR1(59)=59,IGR1(60)=60,IGR1(61)=61,IGR1(62)=62,IGR1(63)=63,IGR1(64)=64,IGR1(65)=65,IGR1(66)=66,IGR1(67)=67,IGR1(68)=68,IGR1(69)=69,IGR1(70)=70,IGR1(71)=71,IGR1(72)=72,IGR1(73)=73,IGR1(74)=74,IGR1(75)=75,IGR1(76)=76,IGR1(77)=77,IGR1(78)=78,IGR1(79)=79,IGR1(80)=80,IGR1(81)=81,IGR1(82)=82,IGR1(83)=83,IGR1(84)=84,IGR1(85)=85,IGR1(86)=86,IGR1(87)=87,IGR1(88)=88,IGR1(89)=89,IGR1(90)=90,IGR1(91)=91,IGR1(92)=92,
   >
IGR2(1)=151,IGR2(2)=285,IGR2(3)=419,IGR2(4)=553,IGR2(5)=687,IGR2(6)=821,IGR2(7)=955,IGR2(8)=1089,IGR2(9)=1223,IGR2(10)=1357,IGR2(11)=1491,IGR2(12)=1625,IGR2(13)=1759,IGR2(14)=1893,IGR2(15)=2027,IGR2(16)=2161,IGR2(17)=2295,IGR2(18)=2429,IGR2(19)=2563,IGR2(20)=2697,IGR2(21)=2831,IGR2(22)=2965,IGR2(23)=3099,IGR2(24)=3233,IGR2(25)=3367,IGR2(26)=3501,IGR2(27)=3635,IGR2(28)=3769,IGR2(29)=3903,IGR2(30)=4037,IGR2(31)=4171,IGR2(32)=4305,IGR2(33)=4439,IGR2(34)=4573,IGR2(35)=4707,IGR2(36)=4841,IGR2(37)=4975,IGR2(38)=5109,IGR2(39)=5243,IGR2(40)=5377,IGR2(41)=5807,IGR2(42)=5932,IGR2(43)=6057,IGR2(44)=6182,IGR2(45)=6307,IGR2(46)=6432,IGR2(47)=6557,IGR2(48)=6682,IGR2(49)=6807,IGR2(50)=6932,IGR2(51)=7057,IGR2(52)=7182,IGR2(53)=7307,IGR2(54)=7432,IGR2(55)=7557,IGR2(56)=7682,IGR2(57)=7807,IGR2(58)=7932,IGR2(59)=8057,IGR2(60)=8182,IGR2(61)=8307,IGR2(62)=8441,IGR2(63)=8575,IGR2(64)=8709,IGR2(65)=8843,IGR2(66)=8977,IGR2(67)=9111,IGR2(68)=9245,IGR2(69)=9379,IGR2(70)=9513,IGR2(71)=9647,IGR2(72)=9781,IGR2(73)=9915,IGR2(74)=10049,IGR2(75)=10183,IGR2(76)=10317,IGR2(77)=10451,IGR2(78)=10585,IGR2(79)=10719,IGR2(80)=10853,IGR2(81)=10987,IGR2(82)=11121,IGR2(83)=11255,IGR2(84)=11389,IGR2(85)=11523,IGR2(86)=11657,IGR2(87)=11791,IGR2(88)=11925,IGR2(89)=12059,IGR2(90)=12193,IGR2(91)=12327,IGR2(92)=12461,IGR2(93)=12595,IGR2(94)=12729,IGR2(95)=12863,IGR2(96)=12997,IGR2(97)=13131,IGR2(98)=13265,IGR2(99)=13399,IGR2(100)=13533,IGR2(101)=13963,IGR2(102)=14088,IGR2(103)=14213,IGR2(104)=14338,IGR2(105)=14463,IGR2(106)=14588,IGR2(107)=14713,IGR2(108)=14838,IGR2(109)=14963,IGR2(110)=15088,IGR2(111)=15213,IGR2(112)=15338,IGR2(113)=15463,IGR2(114)=15588,IGR2(115)=15713,IGR2(116)=15838,IGR2(117)=15963,IGR2(118)=16088,IGR2(119)=16213,IGR2(120)=16338,
   >    nstep1=0, nstep2=0,
   >  &end


#
# --- STARTING CONDITIONS FROM MDOUT
#
   > R1 = -99.000 R2 =   0.000 R3 =   0.000 R4 =  99.000 RK2 =   5.000
RK3 =    5.000
   > Rcurr:    0.037  Rcurr-(R2+R3)/2:    0.037
MIN(Rcurr-R2,Rcurr-R3):    0.037


#
# --- RESULTS FROM MDOUT
#
   >
--------------------------------------------------------------------------------
   >    4.  RESULTS
   >
--------------------------------------------------------------------------------
   >
   >
   >
   >  NSTEP =        1   TIME(PS) =     925.002  TEMP(K) = 298.48  PRESS
=    24.9
   >  Etot   = **************  EKtot   =     25421.1250 EPtot      =
**************
   >  BOND   =      2076.2292  ANGLE   =      8505.3570 DIHED     
=      5083.6430
   >  1-4 NB =      1858.0038  1-4 EEL =    -22464.4248 VDWAALS   
=      4273.3586
   >  EELEC  =    -98833.1557  EHBOND  =         0.0000 RESTRAINT  =
**************
   >  EAMBER (non-restraint)  =    -99500.9889
   >  EKCMT  =      7057.5051  VIRIAL  =      6851.8744 VOLUME     =   
382196.4740
   > SURFTEN    =        30.6828
   > Density    =         1.0174
   >
------------------------------------------------------------------------------
   >
   >  NMR restraints: Bond =    0.015   Angle =     0.000 Torsion =    
0.000
   >
===============================================================================
   >
   >  NSTEP =        2   TIME(PS) =     925.004  TEMP(K) = 297.72  PRESS
=    -4.3
   >  Etot   = **************  EKtot   =     25356.4648 EPtot      =
**************
   >  BOND   =      2095.6856  ANGLE   =      8523.1258 DIHED     
=      5092.7207
   >  1-4 NB =      1861.4326  1-4 EEL =    -22455.4085 VDWAALS   
=      4265.0974
   >  EELEC  =    -98822.5718  EHBOND  =         0.0000 RESTRAINT  =
**************
   >  EAMBER (non-restraint)  =    -99439.9183
   >  EKCMT  =      7037.1397  VIRIAL  =      7072.4544 VOLUME     =   
382197.2895
   > SURFTEN    =        -1.2366
   > Density    =         1.0174
   >
------------------------------------------------------------------------------
   >
   >  NMR restraints: Bond =    0.016   Angle =     0.000 Torsion =    
0.000
   >
===============================================================================
   >
[...]


#
# --- DUMPAVE 06_ntp_equi.rstout
#
   >        0     x:    0.679    y:     0.728    z: 0.056        0.997
   >        1     x:    0.679    y:     0.729    z: 0.056        0.998
   >        2     x:    0.680    y:     0.730    z: 0.056        0.999
   >        3     x:    0.680    y:     0.731    z: 0.057        1.000
   >        4     x:    0.681    y:     0.732    z: 0.057        1.001
   >        5     x:    0.682    y:     0.733    z: 0.058        1.002
   >        6     x:    0.683    y:     0.733    z: 0.059        1.003
   >        7     x:    0.683    y:     0.733    z: 0.060        1.004
   >        8     x:    0.684    y:     0.732    z: 0.062        1.004
   >        9     x:    0.684    y:     0.732    z: 0.063        1.004
[...]

-- 
David Bickel
PhD Student
Computational Pharmaceutical Chemistry and Molecular Bioinformatics
Heinrich-Heine-Universität Düsseldorf
Institut für Pharmazeutische und Medizinische Chemie
Universitätsstr. 1
D-40225 Düsseldorf
Germany
Office: 26.23.00.22
Phone:  +49 211 8112528
URL:    http://cpclab.uni-duesseldorf.de
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Wed Mar 06 2019 - 01:30:03 PST
Custom Search