Date: Thu, 25 Apr 2002 15:25:40 +0200
From: David Smith
Subject: Re: Number of WATBOX

yuann wrote:
>
> Hello all,
> When adding 'cap' of water, I'd like to
> increase my number of water, instead of WATBOX216.
> Does there any command of 'solvateCap' work, or how
> to create a new WATBOX? Thanks!
>
> sychen

I have an answer and a related question.

ANSWER:

Although WATBOX216 contains 216 waters LEaP is smart enough to put as
many waters as you need. This is acheived simply by specifying the
centre and the radius of a sphere of waters you would like to add.

eg: For example I have a protein that I gave the unit name OXM in LEaP.
The central atom of the substrate is at coordinates { 70.7 , 39.3 ,
117.9 }. I wanted to solvate this system with a cap of 41 Angstrom
radius. In LeAP (after loading my pdb file and paramters etc.) I typed:

solvatecap OXM WATBOX216 { 70.7 , 39.3 , 117.9 } 41.0

This added 4279 waters to my system. Whereas setting the closeness value
(see manual) to 0.5 instead of the default 1.0:

solvatecap OXM WATBOX216 { 70.7 , 39.3 , 117.9 } 41.0 0.5

gave 5597 waters.

(PS. You don't have to use cartesian cordinates, you can use the residue
or atom name if you want to)

QUESTION:

I was not overly satisfied with using the cap option with my particular
protein as it was taking a long time for water to diffuse properly into
my system. I found I had to continually add water. to the outside to
make up for that which had diffused into my barrel. I noticed that if I
did constant pressure simulations with periodic boundaries and
equilibrated until I got a good density, everything was a bit nicer. My
brilliant idea was to then extract my cap (non-periodic sphere) from
this system that I knew had a good density and run from there.

So I used the CUTRES command in CARNAL (nice code) to make a group of
appropriate radius and subsequently to generate a pdb file. This works
nicely and I can look at my roughly spherical ball in rasmol etc. Now,
of course, I want to make cordinate and topology files of this new
system and go merrily on my way.

However, I seem to have got stuck with LeAP.

The new file has 7144 waters using the same 41 angstrom cutoff as above
(quite different to that produced using solvatecap and the default
closeness). However, the numbering system is still that of the periodic
box which had some 25000 waters. When I try to read this new pdb file
into LEaP I get a whole bunch of:

-- residue 3071: duplicate H1 atoms (total 6)
-- residue 3071: duplicate H2 atoms (total 6)
-- residue 3071: duplicate O atoms (total 6)
-- residue 3327: duplicate H1 atoms (total 4)
-- residue 3327: duplicate H2 atoms (total 4)
-- residue 3327: duplicate O atoms (total 4)
-- residue 3583: duplicate H1 atoms (total 6)
-- residue 3583: duplicate H2 atoms (total 6)
-- residue 3583: duplicate O atoms (total 6)

ATOM NAMES IN EACH RESIDUE MUST BE UNIQUE:
(same-name atoms are reduced to a single atom)

and I only get 2950 waters if I write a pdb file from here. I thought
the problem might be with TER cards so I added one for every water but
it didn't help. Its hard to tell exactly what's wrong as LEaP seems to
(sensibly) renumber the residues. It seems to me, however, that when
reading a residue number with more than 4 digits, the first one is lost.
That is, residue 5306, residue 15306, and residue 25306 are all read as
residue 5306. My c is terrible so I thought I might ask for somebody
more qualified to help before I destroyed the very nice code.

So finally the question(s) are:

Is this the correct explanation ?
Is there a workaround in LEaP ?
Is there a workaround in CARNAL ?

Thanks in advance for all of your ingenious solutions.

Dave.