[gmx-users] stacked bilayers and pbc errors

[Senthil Kandasamy]

I created a larger bilayer by stacking an existing small bilayer.

genconf -f small.gro -o large.gro -nbox 2 2 1.

Now, I wanted to energy minimize the larger bilayer and run md for 
atleast 100 ps to equilibrate before I do anuything else with it.

on running md, I get the following errors .( I have posted only the 
relevant errors in chronological order in md.log, chopping out 
irrelevant sections)


_________

first sign of error:


12748  12749   46.9    1.1193   0.3750      0.1530
Constraint error in algorithm Lincs at step -2
Going to use C-settle (15384 waters)
wo = 0.888096, wh =0.0559521, wohh = 18.0154, rc = 0.081665, ra = 0.00645837
rb = 0.051255, rc2 = 0.16333, rone = 1, dHH = 0.16333, dOH = 0.1

.
.
.
deleted a whole bunch of output
.
.
.


   12749  12750   60.0    0.1721   0.2150      0.1530
Constraint error in algorithm Lincs at step -1
Started mdrun on node 0 Sat Dec 14 17:09:55 2002
Initial temperature: 1.00068e+09 K
Grid: 18 x 16 x 19 cells

Testing x86 processor CPUID...

Testing x86 SSE capabilities...
CPU and OS support SSE.
Using Gromacs SSE single precision assembly innerloops.

            Step           Time         Lambda      Annealing
               0        0.00000        0.00000        1.00000

There are 58952 atoms in your xtc output selection
    Rel. Constraint Deviation:  Max    between atoms     RMS
        Before LINCS         1333075.000000  12263  12264   28219.177734
         After LINCS         816884.687500  12263  12264   17073.349609


.
.
.
deleted a whole lot  more of the output
.
.
.
   12770  12771   83.5    0.1530  80.9084      0.1530
   12771  12772   92.7    0.1530 316.1893      0.1530
   12772  12773   88.9    0.1530 670.0599      0.1530
   12773  12774   91.2    0.1530 435.2308      0.1530
   12774  12775   89.7    0.1530 639.5835      0.1530
   12775  12776   91.0    0.1530  65.1449      0.1530
   12776  12777   86.3    0.1530  13.2053      0.1530
   12777  12778   95.3    0.1530   1.7792      0.1530
   12778  12779   73.8    0.1530   0.5516      0.1530
   12779  12780   31.3    0.1530   0.1794      0.1530
Constraint error in algorithm Lincs at step 0

t = 0.000 ps: Water molecule starting at atom 52749 can not be settled.
Check for bad contacts and/or reduce the timestep.Wrote pdb files with 
previous and current coordinates
    Energies (kJ/mol)
           Angle    Proper Dih. Ryckaert-Bell.  Improper Dih.          LJ-14
     4.34342e+05    1.48698e+04    3.58255e+04    4.46091e+04    4.56380e+05
      Coulomb-14        LJ (SR)   Coulomb (SR)      Potential    Kinetic En.
     3.01652e+04    2.38870e+12   -9.84130e+05    2.38870e+12            nan
    Total Energy    Temperature Pressure (bar)
             nan            nan            nan
_________________________________________ these are from md.log

the onscreen output has:

warning : Only triclinic boxes with the first verctor parallel to the 
x-axis and box (3x3):
Box [ 0] = {nan, nan, nan}
...
can not fix pbc.
_______________________________________________________


Reducing the timestep does not do any good.

My question is , how can simple stacking cause bad contacts?
if the inital configuration is well equilibrated, shouldn't the larger 
stacked box also be well equilibrated?
how does this affect the pbc?
should i have used some sort of a pbc option when I stacked the bilayers 
using genconf?

my mdp file looks like this:

> ;tieleman's dppc64 modified by senthilk to run dppc256 
> ;
> ;	Input file
> ;
> title               =  dppc256 bilayer
> cpp                 =  /lib/cpp
> constraints         =  all-bonds
> integrator          =  md
> dt                  =  0.001	; ps !
> nsteps              =  25000	; total 25 ps.
> nstxout             =  2500
> nstvout             =  2500
> nstfout             =  2500
> nstlog              =  2500
> nstenergy           =  500
> nstxtcout           =  500
> nstlist             =  10
> ns_type             =  grid
> rlist               =  1.0
> rcoulomb            =  1.0
> rvdw                =  1.0
> pbc			=  xyz
> ; Berendsen temperature coupling is on in two groups
> Tcoupl              =  berendsen
> tc-grps		    =  DPC	SOL 
> tau_t               =  0.1	0.1	
> ref_t               =  323	323 
> ; Energy monitoring
> energygrps          = DPC    SOL 
> ; Isotropic pressure coupling is now on
> Pcoupl              = berendsen
> Pcoupltype          = isotropic
> tau_p               =  0.5
> compressibility     =  4.5e-5
> ref_p               =  1.0
> ; Generate velocites is on at 300 K.
> gen_vel             =  yes
> gen_temp            =  323.0
> gen_seed            =  173529