From m10!frisch@uunet.UU.NET Mon Sep 21 18:50:21 1992 Date: Mon, 21 Sep 92 22:50:21 EDT From: m10!frisch@uunet.UU.NET (Michael Frisch) Subject: Re: More details on MOPAC/G90 problem. To: chemistry@ccl.net uunet!helix.nih.gov!milan (Milan Hodoscek) writes Dear netters, I want to put some more details about the problem of different results reported from G90 and MOPAC for the same systems. MNDO MNDO(precise) AM1 MOPAC 21.87985 21.87987 23.06206 G90 21.76903 22.96633 tokcal(mopac) 630.7040 630.1250 8. I changed also position of first atom in benzene for 0.001 A and the difference in energy was only 0.00196 kcal/mole. The difference in porgrams is 0.01082. I consider 0.001 A significant change in geometry in terms of roundoff errors. Conclusion: It seems to me that this isn't roundoff error nor conversion factor problem, but either of the two: 1. parameters for C and/or H atoms are different 2. code is not implemented correctly Regards -- Milan Hodoscek The comparison of the difference in energies vs. the change upon moving an atom is somewhat misleading. I can reproduce the reported energy difference between MOPAC version 6 and Gaussian 92. However, the difference seems to be nearly constant -- optimizing the geometry of benzene with PRECISE in MOPAC and OPT=TIGHT in Gaussian 92 gives agreement in the 2 bond distances to 5 decimal places, namely C-C=1.39504 and C-H=1.09965. The difference in energies at the first point is 0.0942 kcal at the HF/6-31G* structure and 0.0962 at the final optimized AM1 geometry. Calculations on bare C and H atoms (UHF) agree to 5 decimal places printed in MOPAC. So, the source of the 0.1 kcal/mole energy difference isn't clear, but the predicted structures and relative energies are very close between the two programs. Mike Frisch ------- From jwest@ndb.rutgers.edu Tue Sep 22 01:59:53 1992 Date: Tue, 22 Sep 92 05:59:53 -0400 From: jwest@ndb.rutgers.edu (John Westbrook) To: rpm@wag.caltech.edu Subject: Mathematica modules for Brackets Greetings... I would be interested in any information that you obtain on evaluating matrix elements using Mathematica. Thanks.. John Westbrook jwest@ndb.rutgers.edu From m10!frisch@uunet.UU.NET Tue Sep 22 02:24:41 1992 Date: Tue, 22 Sep 92 06:24:41 EDT From: m10!frisch@uunet.UU.NET (Michael Frisch) Subject: Re: More details on MOPAC/G90 problem. To: chemistry@ccl.net Mark Thompson suggest that a lack of rotational invariance might account for the descrepancy. Nice try, but that's not it. The difference between using the standard and z-matrix orientations in Gaussian 92 is 5x10^-7 kcal/mole. Similarly, using two different orientations in MOPAC 6.0 produces no change in the 5 decimal places printed in the energy. Mike Frisch ------- From milan@helix.nih.gov Tue Sep 22 04:28:46 1992 Date: Tue, 22 Sep 92 08:28:46 -0400 From: milan@helix.nih.gov (Milan Hodoscek) To: chemistry@ccl.net Subject: G90/MOPAC: Problem solved -- new problem appeared Dear netters, Thanks to Ross Nobes the problem with differences in energy reported by Gaussian MNDO method and MOPAC is solved, but grep on mopac directory reports the following: grep 0\.52 *.f analyt.f: A0=0.529167D0 block.f: DATA QQPM3 ( 7)/ 0.5293383D0/ block.f: DATA ADPM3 ( 8)/ 0.5299630D0/ block.f: DATA AQPM3 ( 34)/ 0.5283737D0/ block.f: DATA GUESP3( 52,2)/ 0.5242430D0/ delri.f: A0=0.529167D0 ders.f: A0=0.529167D0 diat.f: R=R1/0.529167D0 diat2.f: RAB=R12/0.529167D0 esp.f: BOHR = 0.529167D00 esp.f: DATA BOHR/0.529167D0/ esp.f: BOHR = 0.529167D00 esp.f: DATA BOHR/0.529167D0/ esp.f: DATA BOHR/0.529167D0/ esp.f: DATA BOHR/0.529167D0/ gover.f: R=R/0.529167D0 hcore.f: HTERME = -0.529177D00*DD(NI)*EFIELD(1)*FLDCON hcore.f: HTERME = -0.529177D00*DD(NI)*EFIELD(2)*FLDCON hcore.f: HTERME = -0.529177D00*DD(NI)*EFIELD(3)*FLDCON powsq.f: 350 BMAT(K,ID) = SIG(K)/0.529167D0 repp.f: DATA A0/0.529167D0/ ,EV/27.21D0/, EV1/13.605D0/, EV2/6.8025D0/, HTERME has different constant then the rest of program and this values go into electron core interaction matrix. Regards -- Milan Hodoscek From alper@zohar.bms.com Tue Sep 22 06:45:25 1992 Date: Tue, 22 Sep 92 10:45:25 EDT From: alper@zohar.bms.com (Howard Alper) Subject: dieelctric constants To: chemistry@ccl.net monica(and other interested in the dielectric constant), the computation of the dieelctric constant can be a fairly complicated matter, but for simple dipolar fluids(and polyatomic fluids that to a good approximation can be represented as dipolar fluids) the problem of computing the dielectric constant has been solved in a series of wonderful papers by m. neumann and o. steinhauser in the late 70's and early 80's. a good reference would be: m. neumann, mol. phys. 50, p841(1983), in which a lot of the basic theory is laid out. some other references are: m. neumann, o. steinhauser and g. s. pauley, mol. phys. 52, p97(1984), m. neumann and o. steinhauser, chem. phys. lett. 95(4,5), p417(1983), m. neumann, j. chem. phys. 85, p1567(1986). it is essential in the computation of the dielectric constant to effectively account for the long-range coulombic forces, either through a reaction field or an ewald sum. you may want to check if your program includes these corrections for the long-range electrostatic forces. here's a signpost: when you calculate the so-called GHk factor(which = /Nm**2, where M is the total system diopole moment, N= # molecules, and m=monomer dipole moment), oops that's the Gk factor, look at the value of Gk. if it is somewhat less than 1.0(like Gk=0.15) then your calculation is probably not accounting for the long-range coulombic contributions with a reaction field or ewald sum, and you will have to get or write code to calculate these important interactions. but once you have done that, then the references above will describe how to calculate the dielectric constant. howard alper From d3f012@gator.pnl.gov Tue Sep 22 01:13:30 1992 Date: Tue, 22 Sep 92 08:13:30 PDT From: d3f012@gator.pnl.gov Subject: Re: G90/MOPAC: Problem solved -- new problem appeared To: chemistry@ccl.net > Dear netters, > > Thanks to Ross Nobes the problem with differences in energy reported > by Gaussian MNDO method and MOPAC is solved, but grep on mopac > directory reports the following: > > grep 0\.52 *.f > > analyt.f: A0=0.529167D0 > block.f: DATA QQPM3 ( 7)/ 0.5293383D0/ > block.f: DATA ADPM3 ( 8)/ 0.5299630D0/ > block.f: DATA AQPM3 ( 34)/ 0.5283737D0/ > block.f: DATA GUESP3( 52,2)/ 0.5242430D0/ > delri.f: A0=0.529167D0 > ders.f: A0=0.529167D0 > diat.f: R=R1/0.529167D0 > diat2.f: RAB=R12/0.529167D0 > esp.f: BOHR = 0.529167D00 > esp.f: DATA BOHR/0.529167D0/ > esp.f: BOHR = 0.529167D00 > esp.f: DATA BOHR/0.529167D0/ > esp.f: DATA BOHR/0.529167D0/ > esp.f: DATA BOHR/0.529167D0/ > gover.f: R=R/0.529167D0 > hcore.f: HTERME = -0.529177D00*DD(NI)*EFIELD(1)*FLDCON > hcore.f: HTERME = -0.529177D00*DD(NI)*EFIELD(2)*FLDCON > hcore.f: HTERME = -0.529177D00*DD(NI)*EFIELD(3)*FLDCON > powsq.f: 350 BMAT(K,ID) = SIG(K)/0.529167D0 > repp.f: DATA A0/0.529167D0/ ,EV/27.21D0/, EV1/13.605D0/, EV2/6.8025D0/, > > > HTERME has different constant then the rest of program and this values > go into electron core interaction matrix. > > Regards -- Milan Hodoscek > Am I seeing the same physical constant being defined multiple times, and sometimes with different values? aargh! IMHO one should define ALL physical constants once-and-only-once, usually in an include file or data statement, and refer to the constant throughout the code with that variable or constant's defined name. ie. #define AU_TO_EV 27.212 It looks like somebody did this with the DATA BOHR statements above, but then someone else? used their favorite value of the same physical constant elsewhere in the code. Although the errors incurred from multiple values are probably small, I wonder how this mixture of constants' values has permeated itself into the values of the method's parameters. I assume the optimization procedure used to obtain the parameters employed this same code? A few years ago, we had the same problem with ZINDO. It showed up when the INDO/S method was being reimplemented in Argus and the answers between the two codes did not agree (much like the current case). I believe the design flaws in ZINDO were subsequently fixed. Mark Thompson d3f012@pnlg.pnl.gov From harbowy@CHEMRES.TN.CORNELL.EDU Tue Sep 22 11:01:35 1992 From: Matthew E. Harbowy Subject: Re: More details on MOPAC/G90 problem. To: chemistry@ccl.net Date: Tue, 22 Sep 92 15:01:35 EDT > Mark Thompson suggest that a lack of rotational invariance might > account for the descrepancy. Nice try, but that's not it. The > difference between using the standard and z-matrix orientations in > Gaussian 92 is 5x10^-7 kcal/mole. Similarly, using two different > orientations in MOPAC 6.0 produces no change in the 5 decimal places > printed in the energy. It is my understanding that rotational invariance is the problem which causes the trivial force constants of a force calculation to be non-zero. This should be corrected for during minimization, but the odd force constants are always bothersome. Point to the net, especially JJPS- what is the synergistic effect of lack of rotational invariance with the various modules of MOPAC? Where is it a problem which should be considered? Rotational invariance has always been a 'bug' of MOPAC but is untreated in the manual. matt ----------------------------------------------------------------------------- Matthew E. Harbowy Cornell U Chemistry Department Ithaca NY 14853 harbowy@chemres.tn.cornell.edu (or) ikf@vax5.tn.cornell.edu ----------------------------------------------------------------------------- -- ----------------------------------------------------------------------------- Matthew E. Harbowy Cornell U Chemistry Department Ithaca NY 14853 harbowy@chemres.tn.cornell.edu (or) ikf@vax5.tn.cornell.edu -----------------------------------------------------------------------------