Dear Vojtech,
The QM region looks fine - this seems to be a case (for the specific Hamiltonian chosen) that needs some additional tricks for SCF convergence beyond what is currently implemented in AMBER. I am sorry I can't help further at the moment.
Thanks for providing examples of coordinate sets that work / don't work which we can use to try and improve the code.
All the best,
Andy
On Nov 15, 2012, at 9:01 AM, <vojtech.mlynsky.upol.cz> <vojtech.mlynsky.upol.cz> wrote:
> Dear all again,
>
>
>
> .. and sorry for prolonging ..
>
> I have finally run the sqm test for only QM-region system in geometry just
> before disruption during QM/MM-MD occurs .. and SCF convergence of Single-
> Point calculation fails (using default AMBER settings)! What is quite
> interesting for me that actual geometry of QM-region (attached at the end)
> is in fact reasonable and *very* similar to that starting one (I have posted
> and tested before..) ..
>
>
>
>
> With best regards,
>
> Vojtech.
>
>
> --
> ==============================
> Department of Physical Chemistry
> Palacky University Olomouc
> tr. 17. listopadu 12, 779 00 Olomouc
> Czech Republic
> ==============================
>
>
>
>
>
>
> QM-region pdb just before disruption ..
>
>
>
>
>
>
> REMARK
>
> ATOM 1 C QMR 1 1.749 1.728 -1.108
>
> ATOM 2 H QMR 1 1.972 1.196 -2.132
>
> ATOM 3 O QMR 1 3.047 1.621 -0.520
>
> ATOM 4 C QMR 1 3.176 0.344 0.224
>
> ATOM 5 H QMR 1 4.059 -0.221 -0.251
>
> ATOM 6 C QMR 1 0.863 0.817 -0.259
>
> ATOM 7 H QMR 1 0.577 1.377 0.689
>
> ATOM 8 C QMR 1 1.845 -0.324 0.035
>
> ATOM 9 H QMR 1 1.528 -0.842 0.998
>
> ATOM 10 O QMR 1 1.975 -1.219 -1.125
>
> ATOM 11 H QMR 1 1.152 -1.897 -1.173
>
> ATOM 12 O QMR 1 -0.367 0.438 -1.005
>
> ATOM 13 P QMR 1 -1.207 -0.805 -0.182
>
> ATOM 14 O QMR 1 -1.337 -0.315 1.244
>
> ATOM 15 O QMR 1 -0.496 -2.040 -0.682
>
> ATOM 16 O QMR 1 -2.685 -0.651 -0.775
>
> ATOM 17 C QMR 1 -3.743 -1.621 -0.402
>
> ATOM 18 H QMR 1 -3.869 -1.443 0.710
>
> ATOM 19 H QMR 1 -3.302 -2.687 -0.669
>
> ATOM 20 N QMR 1 1.552 0.211 7.611
>
> ATOM 21 C QMR 1 2.107 1.530 7.707
>
> ATOM 22 H QMR 1 2.704 1.906 8.624
>
> ATOM 23 N QMR 1 1.933 2.252 6.563
>
> ATOM 24 C QMR 1 1.083 1.433 5.777
>
> ATOM 25 C QMR 1 0.685 1.658 4.448
>
> ATOM 26 O QMR 1 0.999 2.687 3.707
>
> ATOM 27 N QMR 1 -0.059 0.638 3.820
>
> ATOM 28 H QMR 1 -0.341 0.690 2.831
>
> ATOM 29 C QMR 1 -0.209 -0.585 4.508
>
> ATOM 30 N QMR 1 -0.816 -1.556 3.827
>
> ATOM 31 H QMR 1 -1.020 -1.407 2.891
>
> ATOM 32 H QMR 1 -0.972 -2.547 4.132
>
> ATOM 33 N QMR 1 0.241 -0.891 5.782
>
> ATOM 34 C QMR 1 0.915 0.148 6.393
>
> ATOM 35 N QMR 1 -1.365 -0.970 -7.884
>
> ATOM 36 C QMR 1 -0.179 -1.678 -7.779
>
> ATOM 37 H QMR 1 0.407 -2.136 -8.624
>
> ATOM 38 N QMR 1 0.101 -1.969 -6.501
>
> ATOM 39 C QMR 1 -0.870 -1.303 -5.705
>
> ATOM 40 C QMR 1 -1.113 -1.439 -4.287
>
> ATOM 41 N QMR 1 -0.338 -2.081 -3.428
>
> ATOM 42 H QMR 1 -0.451 -2.052 -2.412
>
> ATOM 43 H QMR 1 0.649 -2.266 -3.723
>
> ATOM 44 N QMR 1 -2.324 -0.795 -3.854
>
> ATOM 45 H QMR 1 -2.652 -1.033 -2.901
>
> ATOM 46 C QMR 1 -3.206 -0.223 -4.796
>
> ATOM 47 H QMR 1 -4.059 0.230 -4.278
>
> ATOM 48 N QMR 1 -2.966 -0.107 -6.082
>
> ATOM 49 C QMR 1 -1.836 -0.727 -6.586
>
> ATOM 50 H QMR 1 1.465 2.778 -1.172
>
> ATOM 51 H QMR 1 3.417 0.507 1.274
>
> ATOM 52 H QMR 1 -4.676 -1.450 -0.938
>
> ATOM 53 H QMR 1 1.688 -0.587 8.341
>
> ATOM 54 H QMR 1 -1.900 -0.680 -8.788
>
> END
>
>
>
>
>
>
>
>
>
> ---------- Původní zpráva ----------
> Od: Andreas Goetz <agoetz.sdsc.edu>
> Datum: 14. 11. 2012
> Předmět: Re: [AMBER] Problem with running QM/MM - MD simulation using AM1/d-
> PhoT method in AMBER12
> "Dear all,
>
> It might be an SCF convergence problem in AMBER. Once the SCF fails to
> converge the forces will be wrong and it is not surprising that the system
> blows up.
>
> Vojtech and Mark, as an additional test, could you extract the QM region
> (including link atoms) and check whether you can converge a pure QM single
> point SCF with sqm? You can use 'writepdb = 1' in the &qmmm namelist, from
> which you can extract the coordinates and prepare an sqm input file. Thanks.
>
> All the best,
> Andy
>
> On Nov 14, 2012, at 12:46 PM, Marc van der Kamp wrote:
>
>> Hi all, esp. Gustavo and Brian,
>>
>> I have been helping Vojtech with his project. Let me point out first of
> all
>> that, as far as I'm aware, AM1/d-PhoT was developed in the group of Darrin
>> York especially for the same system that he is working on:
>> J Am Chem Soc. <http://www.ncbi.nlm.nih.gov/pubmed/18345664#
> (http://www.ncbi.nlm.nih.gov/pubmed/18345664#)> 2008 Apr
>> 9;130(14):4680-91. Epub 2008 Mar 18.
>> Quantum mechanical/molecular mechanical simulation study of the mechanism
>> of hairpin ribozyme catalysis.
>> Nam K<http://www.ncbi.nlm.nih.gov/pubmed?term=Nam%20K%5BAuthor%5D&cauthor=
> true&cauthor_uid=18345664
> (http://www.ncbi.nlm.nih.gov/pubmed?term=Nam%20K%5BAuthor%5D&cauthor=true&cauthor_uid=18345664)
>>
>> , Gao J<http://www.ncbi.nlm.nih.gov/pubmed?term=Gao%20J%5BAuthor%5D&
> cauthor=true&cauthor_uid=18345664
> (http://www.ncbi.nlm.nih.gov/pubmed?term=Gao%20J%5BAuthor%5D&cauthor=true&cauthor_uid=18345664)
>>
>> , York DM<http://www.ncbi.nlm.nih.gov/pubmed?term=York%20DM%5BAuthor%5D&
> cauthor=true&cauthor_uid=18345664
> (http://www.ncbi.nlm.nih.gov/pubmed?term=York%20DM%5BAuthor%5D&cauthor=true&cauthor_uid=18345664)
>>
>>
>> http://www.ncbi.nlm.nih.gov/pmc/articles/PMC2655239/
> (http://www.ncbi.nlm.nih.gov/pmc/articles/PMC2655239/)
>> (and reference therein)
>>
>> So the method *should* be perfectly applicable to this system, if
> anything.
>> Minimizations in AMBER do indeed give a reasonable potential energy
>> surface, but the reason to use a semi-empirical method is of course to do
>> some QM/MM MD.
>>
>> Also, we have been using (almost) the same QM region as the above paper,
> as
>> well as slightly larger ones, and the problem still occurs.
>> It seems that in the context Nam & York have used the AM1/d-PhoT
>> parameters, it works fine. This was in CHARMM (with a spherical system
> with
>> stochastic boundary conditions), with the CHARMM27 parameters, *and* with
>> GHO-boundary atoms.
>> This is the reason that we also tried with the CHARMM27 parameters (using
> a
>> CHAMBER prmtop). However, the instability problems still persist. It would
>> be surprising to me if the GHO-boundary atoms vs. link-atoms would make
> the
>> crucial difference... But perhaps they are?
>> The boundary conditions are a difference as well of course, but again not
>> very likely that this causes the extreme QM-region instability. (Vojtech
> is
>> currently running an extra test with qm_ewald=0 and a largish qmcut value
>> to see if that makes a difference...)
>>
>> Any ideas welcome!
>>
>> All the best,
>> Marc
>>
>>
>> On 14 November 2012 19:37, Brian Radak <radak004.umn.edu> wrote:
>>
>>> Hi all,
>>>
>>> Having used AM1/d-PhoT rather extensively, I have found it to be somewhat
>>> numerically sensitive to physically large (as in spread out) QM regions.
> I
>>> don't know if this is unique to the d-orbital implementation (which is
>>> different from the normal sp-orbital NDDO code) or a general
> characteristic
>>> of NDDO. I would give PM6 or MNDO/d a try, as those should also use
>>> d-orbitals for phosphorous (although I would have to double check that in
>>> the code to be sure).
>>>
>>> Regards,
>>> Brian
>>>
>>> P.S. The are virtually no assurances that AM1/d-PhoT is a good model for
>>> nucleobases, as these are entirely outside of the training set. Although
>>> the nitrogen parameters are untouched from AM1, the changing of the C, O,
>>> and H parameters certainly could have unbalanced some very important
>>> effects when it is present. That being said, I'm pretty sure that
>>> AM1/d-PhoT predicts superior proton affinities (wrt experiment) for small
>>> organic compounds compared to AM1, which might salvage some of the
> effects
>>> you are looking for.
>>>
>>> On Wed, Nov 14, 2012 at 1:41 PM, Gustavo Seabra <gustavo.seabra.gmail.com
>>>> wrote:
>>>
>>>> Hi,
>>>>
>>>> On Nov 14, 2012, at 2:15 PM, Vojtech wrote:
>>>>
>>>>> I have also tested another two SE methods (SCC-DFTB and AM1) and did
>>> not
>>>>> observe that problem.
>>>>
>>>> Ok, this *may* be an indication that the SE method you've chose does not
>>>> work very well for your system.
>>>>
>>>>> Furthermore, could you please explain me little bit more what do you
>>>> mean by
>>>>> 'parametrization of QM region'?
>>>>
>>>> Although its not very clear, sander *does* need parameters for the QM
>>>> region. First, because there are steps (as I mentioned in the other
>>> e-mail)
>>>> that are usually done with the system full-MM. Specifically, the initial
>>>> energy minimization and thermalization.
>>>>
>>>> Secondly, there are QM/MM interactions that use MM parameters: the vdW
>>>> interactions between QM-MM atoms utilizes the MM parameters for the QM
>>>> atoms and, if your QM-MM border crosses bonds, any bond, angle or
>>> dihedral
>>>> term involving at least one QM atom is calculated using the regular MM
>>>> equations and parameters. So, the QM region does need MM parameters.
>>>>
>>>> Finally, I had the impression that your QM region is not a standard
>>>> residue. So, for any of those parts, you would need to create the
>>>> appropriate parameters for it.
>>>>
>>>> Gustavo Seabra
>>>> Professor Adjunto
>>>> Departamento de Química Fundamental
>>>> Universidade Federal de Pernambuco
>>>> Recife - PE - Brasil
>>>> +55-81-2126-7450 x5023
>>>>> ---------- Forwarded message ----------
>>>>> From: Gustavo Seabra
>>>>> Date: 14 November 2012 12:40
>>>>> Subject: Re: [AMBER] Problem with running QM/MM - MD simulation using
>>>> AM1/d-
>>>>> PhoT method in AMBER12
>>>>> To: AMBER Mailing List
>>>>>
>>>>>
>>>>> Hi Vojtech,
>>>>>
>>>>> Just to eliminate possibilities: Have you tried a different QM
>>>> hamiltonian?
>>>>> Does the same thing happen?
>>>>>
>>>>> Anyways, how did you parametrize the QM region?
>>>>>
>>>>> What I usually do when running QM/MM calculations is to first energy-
>>>>> minimise the system, then slowly heat it to the desired temperature,
>>> all
>>>>> with full-MM. Only after the system is reasonably thermalised at the
>>>> final
>>>>> temperature I turn on the QM region.
>>>>>
>>>>> Cheers,
>>>>> Gustavo Seabra
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> On Wed Nov 14 2012 08:05:29 Vojtech said:
>>>>>> Dear all,
>>>>>>
>>>>>> I would like to mention here one strange issue that I have recently
>>>>>> encountered in AMBER12.
>>>>>>
>>>>>>
>>>>>>
>>>>>> I try to run QM/MM MD simulation using recently implemented SE -
>>> AM1/d-
>>>>> PhoT
>>>>>> method (abbreviation AM1/d) as the qm_theory for the QM part of
>>>> molecule.
>>>>>> The system is RNA-enzym (Ribozyme) and contains just NA-nucleobases,
>>> Na+
>>>>>> counter-ions and rectangular box of TIP3P molecules (around 30.000
>>> atoms
>>>>> in
>>>>>> total). The QM region contains just one phosphate with the ribose-ring
>>>> at
>>>>>> one site, methoxy group on the other and 2 nucleobases within H-bond
>>>>>> distance from both non-bridging oxygens of that phosphate (typically
>>> 49
>>>>>> atoms). The link atoms are thereby placed instead of C-C bonds (or C-N
>>>>> bonds
>>>>>> in case of nucleobases).
>>>>>>
>>>>>> Problem is that my QM region suddenly 'explodes' during initial
>>>> (typically
>>>>>> 50 ps-long) heating phase (constant volume, example of script at the
>>>> end)
>>>>> as
>>>>>> soon as the temperature of the system reaches ~250 K. Firstly hydrogen
>>>>> atoms
>>>>>> and them the other QM atoms just fly apart. I have also tried
>>> different
>>>>>> modifications of my input script, charmm parameters and charges
>>>>> (converting
>>>>>> starting pdb to charmm psf file, then backwards using chamber module),
>>>>> more
>>>>>> sequentional pre-heating with some equilibration steps and even larger
>>>> QM
>>>>>> regions. All modifications resulted in same error during pre-heating
>>> or
>>>>> just
>>>>>> at the beginning of subsequent equilibration. Notably, initial QM
>>>>>> minimalization fails (SCF convergency) for the system with largest QM
>>>>> region
>>>>>> (100 atoms). I have also tried to heat the system up to 200 K only
>>> and
>>>>>> thereby it 'survises' the pre-heating phase. But during the
>>>> equilibration
>>>>>> phase (constant pressure) same issue occurs (also suddenly, there is
>>> no
>>>>> hint
>>>>>> of something going wrong / any possible trend in both QM/Total
>>> Energies
>>>> ..
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> I am very confused and will much appreciate any reply or suggestion.
>>>>>>
>>>>>> Thanks in advance,
>>>>>> Vojtech
>>>>
>>>> _______________________________________________
>>>> AMBER mailing list
>>>> AMBER.ambermd.org
>>>> http://lists.ambermd.org/mailman/listinfo/amber
> (http://lists.ambermd.org/mailman/listinfo/amber)
>>>>
>>>
>>>
>>>
>>> --
>>> ================================ Current Address =======================
>>> Brian Radak : BioMaPS
>>> Institute for Quantitative Biology
>>> PhD candidate - York Research Group : Rutgers, The State
>>> University of New Jersey
>>> University of Minnesota - Twin Cities : Center for Integrative
>>> Proteomics Room 308
>>> Graduate Program in Chemical Physics : 174 Frelinghuysen Road,
>>> Department of Chemistry : Piscataway, NJ
>>> 08854-8066
>>> radak004.umn.edu :
>>> radakb.biomaps.rutgers.edu
>>> ====================================================================
>>> Sorry for the multiple e-mail addresses, just use the institute
> appropriate
>>> address.
>>> _______________________________________________
>>> AMBER mailing list
>>> AMBER.ambermd.org
>>> http://lists.ambermd.org/mailman/listinfo/amber
> (http://lists.ambermd.org/mailman/listinfo/amber)
>>>
>> _______________________________________________
>> AMBER mailing list
>> AMBER.ambermd.org
>> http://lists.ambermd.org/mailman/listinfo/amber
> (http://lists.ambermd.org/mailman/listinfo/amber)
>
> --
> Dr. Andreas W. Goetz
> Assistant Project Scientist
> San Diego Supercomputer Center
> Tel : +1-858-822-4771
> Email: agoetz.sdsc.edu
> Web : www.awgoetz.de(http://www.awgoetz.de)
>
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
> (http://lists.ambermd.org/mailman/listinfo/amber)"
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
--
Dr. Andreas W. Goetz
Assistant Project Scientist
San Diego Supercomputer Center
Tel : +1-858-822-4771
Email: agoetz.sdsc.edu
Web : www.awgoetz.de
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Thu Nov 15 2012 - 16:00:02 PST
