From giovanni <-at-> sg2.csrsrc.mi.cnr.it Tue Jan 21 06:17:47 1997 Received: from sg2.csrsrc.mi.cnr.it for giovanni "-at-" sg2.csrsrc.mi.cnr.it by www.ccl.net (8.8.3/950822.1) id FAA04166; Tue, 21 Jan 1997 05:21:31 -0500 (EST) Received: from localhost by sg2.csrsrc.mi.cnr.it via SMTP (951211.SGI.8.6.12.PATCH1042/951211.SGI.AUTO) for id CAA08222; Tue, 21 Jan 1997 02:22:01 -0800 Date: Tue, 21 Jan 1997 02:21:58 -0800 (PST) From: Giovanni Scalmani To: Computational Chemistry List Subject: Summary: Electrostatic effects in molecular crystals Message-ID: MIME-Version: 1.0 Content-Type: TEXT/PLAIN; charset=US-ASCII Dear friends, some days ago I submitted a question about "Electrostatic effects in molecular crystals": > >I would like to know if someone has already faced >and (possibly) solved the following problem: > >I would like to model - within an ab-initio calculation - >the crystalline enviroment effect on a single molecule >in a molecular crystal, i.e. I would like to perform a >(possibly high level) calculation of a single molecule's >electronic structure under the influence of a >neighbouring molecules' charge distribution. >An iterative procedure should be necessary to attain >self-consistency between the ab-initio charge distribution >of the reference molecule and the neighbouring molecules' >charge distribution, which will be in turn derived from the >former. > >Since I have no deep knowledge of crystallography, I will >appreciate advises about two issues: >1) Widely speaking: Is such an approach correct? Can you > give me some references in the literature? If I have more > than one molecule (or even molecular species) in the cell > a similar approach on more than one molecule is feasible, > but is it worth while? How should I model the neighbouring > molecules' charge distribution (Mulliken charges, charges > fitted to the molecular electrostatic potential, > Atom-In-Molecule derived charges, orientable atomic > dipoles, ... ) ? >2) More technically: do you know a program which reads the > crystalline structure (for example from the Cambridge > Database) and performs all coordinates manipulation I > will need ? In particular to obtain cartesian coordinates > of any specified molecule in the cell surrounded by all the > neighbouring molecules whose centres of mass are within > a certain radius from the reference molecule's centre of > mass ? > Here I collect the answers I received. Thank you very much to people who have contributed. Giovanni. =Answer 1===================================================================== From: R29CLOSE %-% at %-% ETSU.ETSU-TN.EDU Organization: East Tennessee State University --- Dear Giovanni: I have worked on both these problems. In particular I wrote a program to generate a unit cell (or more) from crystallographic data. I sub- mitted it to CCL last year as DRAWCRYS.FOR. I can help you with all this. In particular, if you send me your xyz coordinates and the space group, I can do the calculations. Then with BABEL I can generate an ALCHEMY file to view the added molecules. The first part is more of a problem. In MOPAC or AMPAC there is a description of "sparkles" These are pseudo charges added to a calculation to simulate the charge of neighbors. But first I need to know what you want to calculate. Often times the easiest thing to do is simulate H-bands between neighboring molecules with water molecules in the right place. ... =Answer 2===================================================================== From: Bart Rousseau Organization: University of Antwerp --- Here are some references on solid state ab initio calculations: Ab-Initio Studies of Crystal Field Effects in Acetylene. P.Popelier, A.T.H. Lenstra, C. Van Alsenoy en H.J. Geise Acta Chemica Scandinavica, A42, 539-543 (1988). An ab-initio study of crystal field effects : Solid and gasphase geometry of acetamide. P. Popelier, A.T.H. Lenstra, C. Van Alsenoy and H.J. Geise Journal of the American Chemical Society, 111, 5658-5660 (1989). An Ab-Initio Study of Crystal Field Effects. Part 3 : Solid and Gas Phase Geometry of Formamide, Modeling the changes in a peptide group due to hydrogen bonds. P. Popelier, A.T.H. Lenstra, C. Van Alsenoy, H.J. Geise Structural Chemistry, 2, 3-9 (1991). Solids Modelled by crystal Fiels ab-initio methods. Part 4 : Thermal Vibrational Parameters and Lattice Expansion Coefficient of the Cubic Phase of Acetylene. K. Verhulst, A.T.H. Lenstra, C. Van Alsenoy, P. Popelier, H.J. Geise. Acta Crystallographica, submitted. Solids Modelled by crystal field ab-initio methods. Part 5 : The phase transitions in biphenyl from a molecular point of view. A.T.H. Lenstra, C. Van Alsenoy, K. Verhulst, H.J. Geise Acta Crystallographica, B50, 96-106 (1994). Solid State Moddeling Part VI : 2,3-diketopiperazine. On the Integration of crystallographic and spectroscopic evidence. A.T.H. Lenstra, B. Bracke, B. Van Dijk, S. Maes, C. Van Alsenoy, Herman O. Desseyn, Spiros P. Perlepes. Acta Crystallographica, submitted. Ab-initio Studies of Crystal Field Effects. Part VII : Structure of 2,3-Diketopiperazine Using a 13-Molecule Cluster, a Calculation involving 1092 Basis Functions. Anik Peeters, C. Van Alsenoy, A.T.H. Lenstra, H.J. Geise International Journal of Quantum Chemistry, 46, 73-80 (1993). Ab-initio studies of crystal field effects. Part 8 : Structure of formamide oxime using a 15-molecule cluster. Anik Peeters, C. Van Alsenoy, A.T.H. Lenstra, H.J. Geise Journal of Molecular Structure (THEOCHEM), 304, 101-107 (1994). Solids Modelled by Crystal Field Ab-Initio Methods. Part 9 : Stereoselective Order-disorder in Tri-ortho-thymotide- 3-Buten-2-ol (2/1) Clathrate. K. Verhulst, A.T.H. Lenstra, C. Van Alsenoy Acta Crystallographica, B51, 1016-1020 (1995). Solids Modelled by Crystal Field Ab-initio Methods. Part 10 : Structure of alpha-glycine, beta-glycine and gamma-glycine using a 15-molecule cluster. A. Peeters, C. Van Alsenoy, A.T.H. Lenstra, H.J. Geise Journal of Chemical Physics, 103, 6608-6616 (1995). Solids Modelled by crystal Fiels ab-initio methods. Part 11 : Integration of Chemical Substitution and Packing Schemes Exploiting Crystallographic and Spectroscopic Evidence Illustrated via Acetamide and Thioacetamide. Koen Verhulst, Stefan Maes, Christian Van Alsenoy, Albert T.H. Lenstra Journal of Molecular Structure (THEOCHEM), submitted. =Answer 3===================================================================== From: Gabriele VALERIO Organization: Ecole Nationale Superieure de Chimie de Montpellier --- Caro Giovanni, I can give you an answer for your second (technical) question. MOLDRAW (http://www.ch.unito.it/ch/DipIFM/Software/MOLDRAW/moldraw.html) is useful to obtain cartesian coordinates or construct z-matrix from crystalline structural data. =Answer 4===================================================================== From: "Gustavo A. Mercier Jr" --- Hi! I can probably shed some light into your problem. Your problem is analogous to performing a QM computation of an active site of a protein including the electrostatic effects of the surrounding protein "environment". A lot of work went into this area during the 80's when I was doing my Ph.D. work. Check out names like Orlando Tapia, Peter vanDuijnen, Hand Dijkman, Arieh Warshel, Michaeil Levitt,Harel Weinstein, and obviously myself. Warshel and Levitt started the whole thing with a paper in J. Mol. Biol. in the 70's. There is also an early paper by Umeyama where he simply included the effect of the "environment" as a set of point charges in the hamiltonian. Today you are in luck. Many ab initio packages include the ability to incorporate "solvent" effects and you can use this feature for your work. Obviously there are several ways of doing this. - Reaction Field methods (I believe G9x) as originally implemented by O. Tapia (Although he did semiempirical methods). - Multipole Representation. Here you have a set of discrete point charges centered on the atoms, or even higher multipole (some even centered on lone pairs etc.). G9x, Hondo, and even ADF allow for computations with point charges. There is a version of Hondo developed by the the people at the National Bureau of Standards (USA) that incorporates higher multipoles. (Check out the work of Morris Krauss). - Multipole Representation plus polarization. This is probably the most complete treatment. You can check my work and Hans Dijkman's. We never implemented this in a very user friendly manner. Indeed, I was doing serial computations with manual steps in between to do so. I then left the field. Hans and I developed complementary mathematical treatments. All of it is predicated in the group function method of McWeeney. I must point out several caveats. As Hans Dijkman (P. vanDuijnen's grad student at the time) pointed out, you cannot have a direct interaction between the "environment" and the SCF density at each SCF step. You need to have an interaction between the converged density of your molecule and the converged polarized "environment". For example: Run 1 : Qm computation in the presence of point charges. effect 1: polarization of the Qm molecule. Run 2: Use the converge density of 1 to polarize the environment. For example, with atom centered polarizabilities. Run 3: Fix the polarized environment and repeat Run 1 using the new environment, but starting the SCF from the output of Run 1. effect 3: Second order polarization of Qm molecule. Run 4: Repeat 2 with the output of Run 3. Run 5: Repeat 3 with the polarized environment of Run 4. ... The above is a true implementation of the method of group functions by McWeeney when one of the groups is "classical" and the other is Qm. Early work in the area was done by H. Weinstein (my Ph.D. advisor). This is tedious and unfortunately I never automated it because In my case it converged extremely quickly. The polarization of a charged environment is very small. The reason I describe all of this is because people took the wrong turn many times. They simply included the environment during the SCF process and this led to wrong result in thre presence of extended basis sets (no convergence!), but "worked" with minimal basis sets. The effect when you include the environment during the SCF process is to incorporate instantaneous fluctuations in the density of the Qm and environment groups. This resulted in spurious term like a the dispersion term that vanDuijnen described. The densities at each SCF step are not physical observables that can be used to compute interactions between groups, as McWeeney shows. Finally, don't forget the factor of (1/2) that occurs in the energy term because the cost of polarizing the systems is (1/2) the polarization energy computed from the converged polarized system (i.e. Qm + environment). This was also a source of trouble then. ... =Answer 5===================================================================== From: "Donald E. Williams" --- Dear Dr. Scalmani: Your questions raise many different issues. In most of my work I have considered that a given molecule is not polarized by surrounding molecules in the crystal. My procedure is to get a good wavefunction for the gas molecule (I use gaussian-94 for that), calculate the MEP on a geodesic grid around the molecule. Then I fit the MEP, using my program PDM97, with net atomic charges, if possible, or adding additional sites (e.g. lone pairs) if necessary. Then to pursue my interest in ab initio crystal structure prediction I put one or more molecules in an arbitrary cell and minimize the total energy including exchange repulsion and dispersion attraction. This mimization is carried out with my program mpa/mpg. If everything works the observed crystal structure (including space group symmetry) is predicted. For your specific interest in locating molecules around a reference molecule in the crystal, mpa/mpg would of course do that. A recent reference can be found in Acta Cryst. 1996, A52, 326-328. -Donald Williams =Answer 6===================================================================== From: iguana -x- at -x- one.net (Ray Crawford) --- Giovanni, Could you please forward the results of this query to me? We are currently trying to do something similar in that we are trying to estimate the energetic differences between a small molecule crystal conformation and a bound macromolecule/ligand crystal structures. The programs we are using to do this are Gaussian and Spartan. Just a few considerations that you might want to keep in mind are: 1.) crystal structure atomic placements are different than force field/ab initio atom placements. What this means is that if you try to do a Single Point calculation on the single molecule crystal structure, you will get an extremely high energy as compared to a nearby low energy conformation. We found that the majority of this energetic difference in due to the placement of the H's (both bond lengths (which are VERY influential) and bond angles/torsions). Other energetic problems result from the placements of the heavy atoms (although these are not as influential). 2.) not many programs allow you to get around these problems. We are currently using Gaussian because of it's ablility to use very explicative z-matrices. Another option are to use Charmm because of its "restrant" option. I hope these two considerations are useful although they may not be exactly related/relavent to what you are doing... But I do realize that we are examining different aspects of similar problems... Any info you could pass along would be greatly appreciated. =Answer 7===================================================================== From: "Francois GILARDONI - Gp. J.WEBER - Univ. Geneve" (The original message was in Italian. I am sorry to Francois for the errors I will put in this translation!) --- I think that there could be many solutions to your problems: 1. Make use of point charges. Gaussian94 can perform that type of calculations. I used instead the DFT-code deMon and I had good results (J.Chem.Phys., 104 (19), p.7624). The set of charges can be reduced by using Ewald summation technique. In that case, the code I used wasn't able to fit the electrostatic field with the desired precision. 2. Make use of Lennard-Jones, Morse or Buckingham functions to fit the potential energy curve obtained from ab-initio calculations (involving the nearest neighbors). According to my experience, Buckingham functions are the most suited ... but it will strongly depend upon the particular situation. I tried to model Ru(bz)2 and [Ru(cp)2]2+ crystals with this method. 3. The Car-Parinello method and the program MARVIN (pseudo-potential) could also be good choice, but I don't know to much about them. When you are using point charges, don't forget to use a distribution that reflects the symmetry of the crystal. ---------------------------------------------------------------------- ^^^ | SCALMANI Giovanni giovanni()at()sg2.csrsrc.mi.cnr.it o o | Universita' degli Studi di Milano | | Dipartimento di Chimica Fisica ed Elettrochimica \_/ | via C.Golgi, 19 Phone: ++39-2-26603254 | 20133 Milano (Italy) Fax : ++39-2-70638129 ----------------------------------------------------------------------