AMBER Archive (2007)

Subject: Re: AMBER: solvateCap minin and mdin for qm-mm dftb

From: M. L. Dodson (mldodson_at_houston.rr.com)
Date: Fri Aug 24 2007 - 15:59:04 CDT

Francesco Pietra wrote:
> I forgot to examine the cap of solvent in xleap. Now, how are inpcrd
rest mdcrd
> files to be handled with VMD? This question because I am now trying to
> visualize the solute and cap with VMD. No problem for inpcrd and rst,
though
> the solute is completely out of the spherical assembly of TIP3PBOX. On
> resizing, the solute is not resized.
>
> With mdcrd there is a further problem, I believe a wrong choice of
filetype
> because the solute molecule is highly distorted with all options I tried.
>
> To make the cap I used
>
> solvateCap EQE TIP3PBOX {1.2 59.9 130.3} 20.0
>
> where the cartesian coordinates for a central carbon atom in the
solute were
> taken from the prepin file from antechamber.
>
> Does any shift of coordinates occur after the stage of prepin?
>
> If all that is not clear, I have to change option to define the
"position",
> perhaps one that Ross suggested.
>
> Thanks
>
> francesco
>
>

Let us use an experiment to answer the question: Take the
coordinate and prmtop files generated by LEaP and make a pdb file
with ambpdb. Look at the pdb file in a text editor and determine
whether the atom you originally chose for the center is "close" to
{1.2 59.9 130.3}.

Perhaps it would have been better if you had taken the center of
the box coordinates from the pdb (or xyz or ...) file that you
used to create the starting coordinates with LEaP.

Bud

>
>
>
> --- "M. L. Dodson" <mldodson_at_houston.rr.com> wrote:
>
>> Francesco Pietra wrote:
>> > Hi Bud and Ross:
>> > Thanks
>> >
>> > When launching the preliminary minimization with adjusted
>> >
>> > Initial min for qmmm with dftb
>> > &cntrl
>> > imin=1, maxcyc=500, ncyc=200,
>> > cut=20.0, ntb=0, ntc=2, ntf=2,
>> > ivcap=0, fcap=10, igb=10,
>> > ifqnt=1
>> > /
>> > &qmmm
>> > qmmask=':1',
>> > qmcharge=0,
>> > qmtheory=7,
>> > qmshake=0,
>> > qm_ewald=0
>> > /
>> >
>> >
>> > sander bombed. min_qmmm_dftb.out ended with:
>> >
>> > LOADING THE QUANTUM ATOMS AS GROUPS
>> > Mask :1; matches 98 atoms
>> > SANDER BOMB in subroutine read_qmmm_namelist
>> > qmtheory = DFTB but qmgb > 1. DFTB does not currently
>> > support full QM GB. Only qmgb values of 0 or 1 are supported.
>> >
>> >
>> > Once flag "qmgb=0" is added just after "qmtheory=7", the
minimization
>> is going
>> > on. From tail -f it is at NSTEP 50.
>> >
>> > I was unable to detect where GB is stated in the minin. I also
wonder
>> about the
>> > difference between "qmgb=0" and qmgb=1" which I did not fin in
either
>> Amber9
>> > manual or the pdf included in Marcus Elnster files.
>> >
>> > If there is nothing wrong about GB statements, I assume that qmgb=0
>> (or qmgb=1)
>> > flag should also be added to mdin, modified (for a trial) as follows
>> >
>> > 300K constant temp for qmmd with dftb
>> > &cntrl
>> > imin=0, ntb=0
>> > cut=20.0, ntc=2, ntf=2,
>> > tempi=300.0, temp0=300.0,
>> > ntt=3, gamma_ln=1.0,
>> > nstlim=1000, dt=0.002,
>> > ivcap=0, fcap=10, igb=10
>> > ntpr=10, ntwx=10,ifqnt=1
>> > /
>> > &qmmm
>> > qmmask=':1',
>> > qmcharge=0,
>> > qmtheory=7,
>> > qmgb=0,
>> > qmshake=0,
>> > qm_ewald=0
>> > /
>> >
>> >
>> >
>> >
>> > Thanks
>> >
>> > francesco
>> >
>>
>> I'm sorry, my bad. I forgot that QMMM in amber9 does not support
>> GB, either. You will not be able to use igb=10. (That is what
>> specifies a GB simulation.)
>>
>> Bud
>>
>> >
>> >
>> > --- "M. L. Dodson" <mldodson_at_houston.rr.com> wrote:
>> >
>> >> Francesco Pietra wrote:
>> >>> I am trying to build prmtop and inpcrd for qm-mm with theory=7
>> (dftb). I am
>> >>> stuck at solvateCap, argument #3 (position).
>> >>>
>> >>> Using *.prepin and *frcmod from previous successful runs with
theory=2,
>> >>>
>> >>> startx
>> >>>
>> >>> $AMBERHOME/exe/xleap -s -f $AMBERHOME/dat/leap/cmod/leaprc.ff99
>> >>>
>> >>> source leaprc.gaff
>> >>>
>> >>> loadamberprep EQE.prepin
>> >>>
>> >>> loadamberparams EQE.frcmod
>> >>>
>> >>> check EQE (unit is OK)
>> >>>
>> >>> solvateCap EQE TIP3PBOX EQE 20.0
>> >>>
>> >>> returns
>> >>>
>> >>> Argument #3 is type Unit must be of type: [molecule residue
atom list]
>> >>>
>> >>> Using as Argument #3: "1" (residue number in prepin), "C17" (a
central
>> >> carbon
>> >>> atom number in prepin, "c" the C17 atom type, or, finally (and
blindly)
>> >>> EQE.2.CA the error is: Argument #3 is type String.
>> >>>
>> >>> I am short of imagination (even consulting the manual and the
>> mailing list)
>> >> of
>> >>> what should be given in my case as Argument #3.
>> >> Without getting into the structured design of LEaP, perhaps the
>> >> easiest way to for you to proceed is to pick the atom you want to
>> >> make the center of your spherical cap. For purposes of
>> >> illustration, let us assume the Cartesian coordinates of that
>> >> atom are:
>> >>
>> >> x = 1.0, y = 2.0, z = 3.0
>> >>
>> >> Then a solvate cap command using these coordinates might be:
>> >>
>> >> solvateCap EQE TIP3PBOX {1.0 2.0 3.0} 20.0
>> >>
>> >> The desc command is your friend when trying to learn the syntax of
>> >> LEaP. You just have to practice.
>> >>
>> >>> ________________--
>> >>> That said, on summing up info from the manual and suggestions from
>> Bud and
>> >>> Gustavo, I have a number of doubts about min.in and md.in, in
>> particular
>> >> if
>> >>> qm_pme should be omitted at all. The two files I have set up read:
>> >>>
>> >>> Initial min for qmmm with dftb
>> >>> &cntrl
>> >>> imin=1, maxcyc=500, ncyc=200,
>> >>> cut=20.0, ntb=1, ntc=2, ntf=2,
>> >>> ivcap=0, fcap=10
>> >>> ifqnt=1
>> >>> /
>> >>> &qmmm
>> >>> qmmask=':1',
>> >>> qmcharge=0,
>> >>> qmtheory=7,
>> >>> qmshake=0,
>> >>> qm_ewald=0
>> >>> /
>> >>>
>> >>> 300K constant temp for qmmd with dftb
>> >>> &cntrl
>> >>> imin=0, ntb=1
>> >>> cut=20.0, ntc=2, ntf=2,
>> >>> tempi=300.0, temp0=300.0,
>> >>> ntt=3, gamma_ln=1.0,
>> >>> nstlim=20000, dt=0.002,
>> >>> ivcap=0, fcap=10,
>> >>> ntpr=100, ntwx=100,ifqnt=1
>> >>> /
>> >>> &qmmm
>> >>> qmmask=':1',
>> >>> qmcharge=0,
>> >>> qmtheory=7,
>> >>> qmshake=0,
>> >>> qm_ewald=0
>> >>> /
>> >>>
>> >>>
>> >>> I want to try if, at the current implementation of qm-mmm/dftb
>> >>> in Amber9, speed is acceptable for a 98-atoms molecule. If
>> >>> acceptable, to compare ESCF and EPtot with AM1.
>> >>>
>> >> Since amber9 sander cannot do pme with QMMM and qmtheory=7, you
>> >> should not turn on qm_pme. Also, for solvent cap simulations,
>> >> ntb=0, not ntb=1 as you have specified. And you should think
>> >> about the suggestion on page 99 of the manual to set igb=10. I
>> >> also think you may have to reduce the timestep. I use 0.5fs for
>> >> QMMM steered MD, but you may be able to get by with 1.0fs, 2.0fs
>> >> is generally workable, but is said to be borderline for classical
>> >> MD. See the email list archives.
>> >>
>> >> Of course the definition of acceptable speed is personal, but I
>> >> have found 90+ atom QM regions to be workable with qmtheory=7.
>> >>
>> >>>
>> >>> Thanks for checking
>> >>>
>> >>> francesco pietra
>> >> HTH,
>> >>
>> >> Bud Dodson
>> >> --
>> >> M. L. Dodson
>> >> Email: mldodson-at-houston-dot-rr-dot-com
>> >> Phone: eight_three_two-five_63-386_one 

-- 
M. L. Dodson
Email:	mldodson-at-houston-dot-rr-dot-com
Phone:	eight_three_two-56_three-386_one
-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber_at_scripps.edu
To unsubscribe, send "unsubscribe amber" to majordomo_at_scripps.edu