[gmx-users] walls...still struggling (not)

[Sherwin J. Singer]

I can report 1) a problem solved, and 2) unexpected behavior
concerning genbox input options.

1) I've written recently about constructing flexible silanol groups on
an otherwise frozen silicate wall.  Gromacs md would move the silanol
oxygens and hydrogens, but would claim the kinetic energy was zero,
gave "nan" for the number of degrees of freedom, and would not
thermostat or energy minimize the structure.

My problem turned out to be the center of mass constraints.  I had
given the two groups for center of mass removal as "comm_grps = WALL1
WALL2", where the WALL1 group is frozen, and WALL2 is the silanols.
When I simply specify "comm_grps = SYSTEM", all is well.  I'm not sure
why removing 3 center of mass degree of freedom among 448 silanol
atoms messed me up, but it did.

2) I then went to solvate my wall.
"genbox -cp walls_run1.gro -p walls.top" added 0 waters, while
"genbox -cp walls_run1.gro -cs spc216.gro -p walls.top" added 9180
waters.  Isn't spc216.gro the default for the -cs option?

Regards,
  Sherwin Singer


On Wed, 2003-11-12 at 03:31, David van der Spoel wrote:
> On Wed, 2003-11-12 at 00:30, Sherwin J. Singer wrote:
> > David:
> > 
> > Thanks for your reply!  I just got back into town, so here are the
> > answers to your questions.
> > 
> > * I have harmonic bonds between O and H of the silanol groups, and a
> > harmonic Si-O-H angle function.  There is no dihedral around the Si-O
> > bond.  I agree that the OH could undergo azimuthal rotation around the
> > Si-O axis, although I would expect non-bonded interactions to make it
> > less than perfectly free.
> > 
> > * I stayed away from dihedrals, since I thought that would require a
> > periodic wall or dummy atoms to do that right.  I was trying to get
> > something simple working, and then build from there.  That is why at the
> > moment the silanol hydrogens are the only atoms that are not frozen.
> > 
> OK, check the rotation anyway once you get it going.
> 
> > * Yes, there are charges on all the wall atoms, both frozen and
> > un-frozen.  SI atoms have a charge of +.3, OS atoms between silicons
> > have a charge of -.15.  OA and HO atoms of the silanol groups have
> > charges of -.475 and .4, respectively.
> > 
> > I am curious why are you asking about the potentials and the charges. 
> > Even if my model was poorly chosen, wouldn't I still be able to minimize
> > or cool down the initial configuration?  Would the model cause, for all
> > groups, the output "Number of degrees of freedom in T-Coupling group ...
> > is nan"?
> > 
> The number of DOF for a group is computed from 3*natoms - Nconstr. Now
> if the atoms are frozen, the nDOF is zero, it may be that something is
> divided by zero somewhere. Check what it says in your tpr file. You
> might actually consider to let the oxygens fluctuate as well, or, for
> debugging purposes, to make the OH bond a flexible one.
> 
> 
> > I would think that my input of group memberships was the reason that the
> > silanol H's never minimize or cool down and for the "nan" mentioned
> > above, but the frozen/unfrozen specification seems to work and visually
> > ngmx seems to recognize my choice of groups.
> > 
> > Again, help would be very much appreciated.
> > 
> > Regards,
> >   Sherwin
> > 
> > 
> > 
> > On Sat, 2003-11-08 at 06:15, David wrote:
> > > On Sat, 2003-11-08 at 00:53, Sherwin J. Singer wrote:
> > > > I am attempting to simulate a wall consisting of a few layers of
> > > > silicate running across system in the xy-plane.  Eventually, the wall
> > > > will be in contact with water, but for now I'm trying to just set up
> > > > the wall.  The desired results has some similarities to Wensink,
> > > > Hoffmann, Apol and Berendsen, Langmuir 16:7392 (2000).
> > > > 
> > > > The strategy I took is to make the wall a giant molecule containing
> > > > two residues.  The first residue (WALL1) contains inner frozen inner
> > > > atoms, and the second residue (WALL2) contains atoms that can move.
> > > > To simplify as much as possible now, the moving WALL2 atoms only
> > > > include H-atoms of silanol groups (...Si-OH) at the wall surface.  The
> > > > WALL1 and WALL2 atoms are identified in an index file, and the total
> > > > number of atoms is 4704.
> > > > 
> > > > Some things are working right.  Under the "md" integrator the WALL1
> > > > atoms are really frozen, and only the WALL2 atoms move.
> > > > 
> > > > However, I cannot relax the moving "WALL2" atoms.  I cannot get
> > > > conjugate gradient or steepest descents working, and I cannot even
> > > > cool down the initial structure by running molecular dynamics with
> > > > strong coupling to a thermostat and low temperature.  Those silanol
> > > > H's keep flailing violently.
> > > How are the H connected? Did you define dihedrals around the Si-O bond?
> > > Otherwise they have in principle free rotation. Is there water in your
> > > system? Do you have charges on the Si-O atoms?
> > > 
> > > 
> > > > 
> > > > * With steepest descents, I find 
> > > > 
> > > > Steepest Descents:
> > > >    Tolerance         =  1.00000e+02
> > > >    Number of steps   =         1000
> > > > 
> > > > Stepsize too small (6.10352e-07 nm)Converged to machine precision,
> > > > but not to the requested precision (100)
> > > > Steepest Descents did not converge in 15 steps
> > > >   Potential Energy  =  1.52274e+08
> > > > Maximum force:  2.37367e+02
> > > > 
> > > > * With md at 20K, dt=.0005, tcoupl=berendsen, tau_t=.01, ref_t=20, I
> > > > generate a trajectory but no apparent cooling even after thousands of
> > > > steps.  The following output from grompp makes me think I've done
> > > > something wrong.  Please note the "nan" for the number of degrees of
> > > > freedom in the T-Coupling groups.
> > > > 
> > > > ==========================================================================
> > > > checking input for internal consistency...
> > > > calling /lib/cpp...
> > > > processing topology...
> > > > Generated 1125 of the 1326 non-bonded parameter combinations
> > > > Excluding 3 bonded neighbours for WALL 1
> > > > processing coordinates...
> > > > double-checking input for internal consistency...
> > > > renumbering atomtypes...
> > > > converting bonded parameters...
> > > > #      BONDS:   1344
> > > > #     ANGLES:   1788
> > > > initialising group options...
> > > > processing index file...
> > > > Making dummy/rest group for Acceleration containing 4704 elements
> > > > Making dummy/rest group for Freeze containing 448 elements
> > > > Making dummy/rest group for Energy Mon. containing 4704 elements
> > > > Number of degrees of freedom in T-Coupling group WALL1 is nan
> > > > Number of degrees of freedom in T-Coupling group WALL2 is nan
> > > > Making dummy/rest group for User1 containing 4704 elements
> > > > Making dummy/rest group for User2 containing 4704 elements
> > > > Making dummy/rest group for XTC containing 4704 elements
> > > > Making dummy/rest group for Or. Res. Fit containing 4704 elements
> > > > T-Coupling       has 2 element(s): WALL1 WALL2
> > > > Energy Mon.      has 1 element(s): rest
> > > > Acceleration     has 1 element(s): rest
> > > > Freeze           has 2 element(s): WALL1 rest
> > > > User1            has 1 element(s): rest
> > > > User2            has 1 element(s): rest
> > > > VCM              has 2 element(s): WALL1 WALL2
> > > > XTC              has 1 element(s): rest
> > > > Or. Res. Fit     has 1 element(s): rest
> > > > ==========================================================================
> > > > 
> > > > Help would be greatly appreciated.
> > > > 
> > > > Thanks,
> > > >   Sherwin Singer
-- 
Sherwin J. Singer
Department of Chemistry
Ohio State University
Columbus OH, 43214
singer at chemistry.ohio-state.edu
614-292-8909
614-292-1685 (fax)