Category Archives: Software
Some atomic properties such as an atomic charge are isotropic, but every now and then some derivations of them become anisotropic, for example the plotting of the Molecular Electrostatic Potential (MEP) on the electron density surface can exhibit some anisotropic behavior; quantifying it can be a bit challenging.
It is well known that halogen atoms such as Chlorine can form so-called halogen-bonds of the type R-Cl-R in crystals with a near perfect 180° angle. This finding has lead to the discovery of σ-holes in halogens. σ-Holes are electrophilic portions of the anisotropic electrostatic potential in an otherwise nucleophilic atom. Recently, Guillermo “Memo” Caballero and I calculated the MEP for a series of trichloromethyl-containing compounds at the MP2/cc-pVQZ level of theory and the mapping shows evidence of such σ-holes as seen in Figure1. Those small blue portions on an otherwise red atom indicate that some electron density is missing in that position, which by the way is located at 180° away from the carbon atom.
But having the picture is not enough. We want to quantify just how strong are those σ-holes to effectively attract a nucleophile and perhaps perform some chemistry on the C-Cl bond. That’s when we resorted to MultiWFN, a Multifunctional Wavefunction Analyzer developed by Tian Lu (卢天) at the Beijing Kein Research Center for Natural Sciences. You can check the project leader list of publications here. Among many other capabilities, MultiWFN is able to print details about properties along a surface.
In order to work with MultiWFN you need to generate a *.wfn file, if you have a previous Gaussian calculation for which you want to analyze their surface you can run a guess=only calculation in order to extract the wavefunction from the checkpoint file. Here is a dummy of the input for such calculation
%chk=oldfile.chk # output=wfn geom=check guess=(only,read) density=current --blank-- Title Card --blank-- 0 1 --blank-- filename.wfn --blank--
In our case, having a post-Hartree-Fock calculation, the use of density=current is mandatory to get the MP2 density matrix and not just the HF one. Running this calculation will generate the file filename.wfn which is now used with MultiWFN. When starting MultiWFN you get to see a terminal window like the one below in which you are asked to input the path of your wfn file:
After loading it you will get the following window with the various options available. Type 12// (these two slashes are mandatory) to get the quantitative analysis of molecular surface option.
Then you will be asked to define some elements of that surface (we used the default options 0)
The following screenshot shows the results section in which several maxima and minima of electrostatic potential were found (7 and 11 in our case); a star is placed on the side of the global maximum. The value of the MEP at those points is given in Hartrees, eV and kcal/mol which I personally hate because there isn’t such a thing as a mole of ‘potentials’ (same argument as giving an orbital’s energy in kcal/mol, moles of what? orbitals? Personally, I don’t like it even if its valid).
Their visualizer is activated through the option 0 and although it is far from pretty it is quite good enough to find the numbers corresponding to maxima and minima of the MEP on the isodensity surface. If we look for the maxima then we find for our example (CHCl3) that a maximum is located in front of each Cl atom in a straight line from the C atom. Now we get to put a number on the mapped isosurface provided by Gaussian or even import the file into Chimera.
We are still working our way around MultiWFN, I hope we can find the batch option, it would be most useful. In the mean time, Guillermo and I will continue to search for σ-holes in chlorinated reagents. Thanks to Guillermo for his ongoing work in this and other topics within the realm of organic reactivity.
Have you any suggestions or ideas to work with MultiWFN? Please share them!
A bit outside the scope of this blog (maybe), but just too cool to overlook. Augmented reality in chemistry education.
This is a guest post from Samantha Morra of EdTechTeacher.org, an advertiser on FreeTech4Teachers.com.
Augmented Reality (AR) blurs the line between the physical and digital world. Using cues or triggers, apps and websites can “augment” the physical experience with digital content such as audio, video and simulations. There are many benefits to using AR in education such as giving students opportunities to interact with items in ways that spark inquiry, experimentation, and creativity. There are a quite a few apps and sites working on AR and its application in education.
There are 6 physical paper cubes printed with different symbols from the periodic table. It takes a while to cut out and put together the cubes, but it…
View original post 475 more words
It’s been a long time since I last posted something and so many things have happened in our research group! I should catch up with them in short but times have just been quite hectic.
Here is a contribution from Igor Marques at the University of Aveiro in Portugal (Group Website); the original text can be found as a comment in the original NBO Visualization post but it is pretty much the same thing you can find in this post. Here is a link to Chemcraft’s website. Thanks for sharing this, Igor!
=> Examples provided by Igor Marques used Chemcraft Version 1.7, build 365 <=
In the Gaussian input, with the NBORead option included under the population keyword, we should include the PLOT option as illustrated below. The gfoldprint keyword will print the basis set to the output file in the old G03 format. Some visualization programs require a certain format of the basis set to be printed to the output file in order to plot orbitals and other surfaces like the electron density; therefore, if you want to play safe, use gfoldprint, gfprint and gfinput in the same line. gfprint will print the basis set as a list but in the new G09 format, whereas gfinput will print the basis set using Gaussian’s own input format. (The used level of theory and number of shared processors are shown as illustrations only; also the Opt keyword is not fundamental to the visualization of the NBO’s)
%chk=filename.chk %nprocshared=8 #P b3lyp/6-311++g** Opt pop=(full,nboread) gfoldprint filename 0 1 molecular coordinates $NBO BNDIDX PLOT $END
this will generate files from *.31 to *.41
For the visualization of NBOs, you’ll need FILE.31 and FILE.37. Open FILE.31 from chemcraft. It will automatically detect FILE.37 (if in the same directory).
Tools > Orbitals > Render molecular orbitals
select the NBOs of interest (whcih are in the same order of the output),
Adjust settings > OK
On the left side of the window, select the NBO of interest and then click on ‘show isosurface’. Adjust the remaining settings. To represent another orbital, click on ‘keep this surface’ and then select another orbital from the rendered set and follow the previous steps.
> It’s possible to open a formated checkpoint file, containing the NBOs, in chemcraft.
%Chk=filename.chk %nprocshared=4 #P b3lyp/6-311++g** Opt pop=(full,nboread,savenbo) gfoldprint filename 0 1 molecular coordinates $NBO BNDIDX $END
the procedure is identical, but it is only necessary to read the *fchk file and then render the desired orbitals.
However, two problems might arise:
a) Orbitals in the checkpoint are reordered, thus requiring some careful inspection of the output.
b) Sometimes, for a larger molecule, the checkpoint might not be properly saved and the Gaussian job (as previously reported – http://goo.gl/DrSgA ) will end with:
Failed in SchOr1 in NBStor.
Error termination via Lnk1e in /data/programs/g09/l607.exe at Wed Mar 6 15:27:33 2013.
As usual, thanks to all for reading/commenting/rating this and other posts in this blog!
This is the first time I reblog a post from a fellow computational chemist and the reason why I do it is because of its beautiful simplicity and usefulness. Given the scope this blog has taken I think this post becomes most appropriate. This post will show you how to create an energy level diagram using nothing but MS Excel.
Kudos to ‘Eutactic’, from Australia, for coming up with a nice solution to this problem. Check out his blog at eutactic.wordpress.com.
Thanks for letting me repost it 🙂
I worked out a very quick and easy way to generate level schemes in Excel, based on a query from one of the other students in the group. Normally I would resort to something like the astonishing TikZ for this sort of task, however our group is very much a Microsoft Office ‘What You See Is A Metaphor For Cosmic Horror‘ group and recommending that a colleague learns two new markup languages to produce a figure is probably not helpful in the short term. One of the issues with charting energy levels in Excel is that levels are typically represented by horizontal bars connected at their vertices with lines representing transitions. Whilst Excel does have a horizontal bar as a marker, it possesses two show-stopping limitations:
- It is only uniformly scalable, and can only be scaled so far – we cannot make it anywhere near wide and…
View original post 222 more words
Due to extensive popular demand, I hereby make available the necessary files to run Molekel in its old 4.3 version. The program has been compiled to work under Windows 32 bit architecture. Just extract it and place the main folder (as provided here) in any location and run the .exe file located inside. You can generate a direct access to it from your desktop and it even includes a small icon to be used for this purpose.
Also, the manual in pdf format is included; in it you can find the proper citation which must be included in any publication that makes use of Molekel. Just in case you can’t find it, here it is:
MOLEKEL 4.3, P. Flükiger, H.P. Lüthi, S. Portmann, J. Weber, Swiss Center for Scientific
Computing, Manno (Switzerland), 2000-2002.
Stefan Portmann & Hans Peter Lüthi. MOLEKEL: An Interactive Molecular Graphics Tool.
CHIMIA (2000) 54 766-770.
Now some considerations:
- This is an old program. It was generated back in the WindowsXP days, so it wouldn’t be a surprise if it doesn’t work in more recent platforms or under any other 32 bit OS.
- The manual is included. Please read it. This blog is not Molekel’s help desk; I may help but I can’t solve everything, I just don’t have the time for it.
- I didn’t participate/collaborate/helped or got involved in any way in the development of this program, i.e., don’t shoot the messenger. The Molekel homepage is: http://molekel.cscs.ch/wiki/pmwiki.php
- I strongly recommend to make use of the NEW version. A bit more obscure but also great once you figure it out, plus there is available support for it from the actual developers.
- I also strongly recommend to look over the internet for other visualization softwares. I don’t recall having reviewed any in this blog. Perhaps some other time.
- Since this is not my development I will remove it from the server upon the request of the rightful owners A.S.A.P! My guess is they wont mind all that much since its an old version and it was given away for free from their server anyway.
- I can’t think of anything else to put on this list right now but I reserve the right to come back to it and add something more. I just don’t want any trouble.
So, here it is! Right click on the link and download it; Use it to generate nice plots of your orbitals, densities, electrostatic potentials, etc. Consider this a Happy New Year’s gift!
Rate and comment this and all the other posts you find interesting in this blog. Please!
UPDATE: Thanks to Yuekui Wang for the following information.
This copy doen’t work on some WinXP machine with ATI monitor card. The original copy is still available on the cscs web site. Download link is as follows:
It works fine on many macjines, I am sure.
What if you could convince people to help you doing your research on their spare time? What if you could convince a million people to contribute to a specific scientific effort without the need of recruiting them yourself? Even better if you can get all these exo-collaborators without making a huge dent in your budget, which sometimes is just impossible even if you are willing to do it. The world is an interconnected global one these days; millions of voices sound through the web which makes it hard for yours to stand out. As an interconnected entity, communicating to a large mass has become feasible but it comes with a price: you need to be attractive! That’s right, you can get a million people to help your scientific efforts, its called crowdsourcing (think about wikipedia for instance), but you have to make it rewarding in some way, and since you are trying to convince them to work for you in their free time you have to make it look like something they’d do on that free time; So why not making it a video game?
There are nowadays some serious games which are nothing more than an internet-based platform in which tons of data are loaded and accessed by many users who analyze them while being scored on their achievements in many different ways according to the rules of each game. Remember that famous NASA screensaver (SETI@home) which used the idle time on your computer to crunch data from the SETI (Search for Extraterrestrial Intelligence) project; that is a more passive example of exo-collaboration simply called distributed computing, since the user has to do nothing but allowing the entrance of data into their computers for further processing.
Fold.it is a videogame (that actually began as a distributed computing screensaver called Rosetta@home) that allows you to play around with a protein and fold it in many different ways while you score points according to the conformer’s plausibility. Obtaining the native tertiary (and even the secondary) structure of a protein from no other information than the primary structure is extremely difficult given the enormous amount of available degrees of freedom. Molecular dynamics alone is unable to predict the native tertiary structure of the protein; the number p of possible disulfide bonds present in a protein is p = n!/[(n/2)!2^n/2] where n is the number of cysteine residues available, plus computers know nothing about proteins or enzymatic catalysis so a hand from us fellow humans and our chemical insight is widely needed. Therefore our previous knowledge of chemistry, biochemistry and the nature of related proteins can help us help those programs in finding the best possible answer to ‘how does this protein look like in 3D space?’ but since this human-helped process is slow and cumbersome you need thousands of people working on it a great deal of time; almost as if every person playing with the same structure was a single core in your computer. Fold.it thus, is a sort of protein self docking, if you will, in which players are ranked according to their skills and rewarded according to how well your structure complies with three simple rules: 1) lack of voids (packing) 2) keeping the orange hydrophobic chains unexposed to the aqueous exterior and 3) avoiding clashes. Scoring functions for these three concepts are calculated and then yield a score for the player which is then ranked to other players folding the same protein (or to other players in their overall performance).
Fold.it has already collected some major success stories such as the one published on Nature Structural & Molecular Biology by David Baker (founder of Fold.it) et al. on September 2011 (doi:10.1038/nsmb.2119) in which players helped in solving the crystal structure of a protease from a retrovirus which causes AIDS in monkeys. The determination of this structure had already taken 15 years of work with only partial success; but the data was available in Fold.it for only three weeks when the appropriate match to the diffraction experiments was found! This case alone has stirred too much attention and for a beautifully written piece about it, you can check this article at the Discover Magazine by Ed Yong.
Other such examples of crowdsourcing in scince, more specifically in astronomy, are Galaxy Zoo and Moon Zoo in which thousands of images from the Hubble telescope and numerous moon probes are made available for users to sort and classify. The aim of Moon Zoo is to study the amount, shape and occurrence of craters, which basically never erode unlike those on Earth. This analysis will let us know more about the origin of our natural satellite and ultimately about the origins of our solar system.
To the participants in this specific kind of scientific crowdsourcing the term Citizen Science is applied and even publications such as the Scientific American magazine host a section where you can call out for volunteers in your projects. Some sort of classified ads for the lonely scientists in their labs in search for “idle” hands that can make a significant contribution to science. Some Citizen Science projects are intended for kids and teenagers as a way to get more people interested in scientific disciplines by engaging them directly in activities with a measurable progress of their own contributions. It is worth mentioning that projects like Fold.it, Moon Zoo and Galaxy Zoo are developed in a way that can be used by people with no expertise in the field in order to recruit as many people as possible just to perform a very specific task, proving thus that the human brain is a powerful and beautiful machine whose insight isn’t equaled by any artificial system, yet.
Well, it is now time to go back to work before I’m deemed a permanent exo-collaborator by my bosses. Just a final thought: What were our mothers saying about us playing too much with our video games?
2011, International Year of Chemistry
As usual please share your thoughts in the comments section, rate this post and let me know that you are out there reading this.
I am frequently asked how to include an extra set of basis functions in a calculation or how to use an entirely external basis set. Sometimes this question also implies the explicit declaration of an external pseudopotential or Effective Core Potential (ECP).
New basis sets and ECPs are published continuously in specialized journals all the time. The same happens with functionals for DFT calculations. The format in which they are published is free and usually only a list of coefficients and exponents are shown and one has to figure out how to introduce it in ones calculation. The EMSL Basis Set Exchange site helps you get it right! It has a clickable periodic table and a list of many (not all) different basis sets at the left side. Below the periodic table there is a menu from which one can select which program we want our basis set for; finally we click on “get basis set” and a pop-up window shows the result in the selected format along with the corresponding references for citation. A multiple query can be performed by selecting more than one element on the table, which generates a list that almost sure can be used as input without further manipulations. Dr. David Feller is to be thanked for leading the creation of this repository. More on the history and mission of the EMSL can be found on their About page. Because of my experience, the rest of the post addresses the inclusion of external basis sets in Gaussian, other programs such as NwChem will be addressed in a different post soon.
The correct format for inclusion of an external basis set is exemplified below with the inclusion of the 3-21G basis set for Carbon as obtained from the EMSL Basis Set Exchange site (blank lines are marked explicitly just to emphasize their location:
spin multiplicity Molecular coordinates - blank line - C 0 S 3 1.00 172.2560000 0.0617669 25.9109000 0.3587940 5.5333500 0.7007130 SP 2 1.00 3.6649800 -0.3958970 0.2364600 0.7705450 1.2158400 0.8606190 SP 1 1.00 0.1958570 1.0000000 1.0000000 **** - blank line -
The use of four stars ‘****’ is mandatory to indicate the end of the basis set specification for any given atom. If a basis set is to be declared for a second atom, it should be included after the **** line without any blank line in between.
WARNING! Sometimes we can find more than one basis set in a single file this is due to different representations, spherical or cartesian basis sets. Gaussian by default uses cartesian (5D,7F) functions. Pure gaussian use 6 functions for d-type orbitals and 10 for f-type orbitals (6D, 10F). Calculations must be consistent throughout, hence all basis functions should be either cartesian or pure.
Inclusion of a pseudopotential allows for more computational resources to be used for calculation of the electronic structure of the valence shell by replacing the inner electrons for a set of functions which simulate the presence of these and their effect (such as shielding) on the valence electrons. There are full core pseudopotentialas, which replace the entire core (kernel). There are also medium core pseudopotentials which only replace the previous kernel to the full one, allowing for the outermost core electrons to be explicitly calculated. The correct inclusion of a pseudopotential is shown below exemplified by the LANL2DZ ECP by Hay and Wadt for the Chlorine atom.
spin multiplicity Molecular coordinates - blank line - basis set for atom1 **** basis set for atom2 (if there is any) **** - blank line - CL 0 CL-ECP 2 10 d potential 5 1 94.8130000 -10.0000000 2 165.6440000 66.2729170 2 30.8317000 -28.9685950 2 10.5841000 -12.8663370 2 3.7704000 -1.7102170 s-d potential 5 0 128.8391000 3.0000000 1 120.3786000 12.8528510 2 63.5622000 275.6723980 2 18.0695000 115.6777120 2 3.8142000 35.0606090 p-d potential 6 0 216.5263000 5.0000000 1 46.5723000 7.4794860 2 147.4685000 613.0320000 2 48.9869000 280.8006850 2 13.2096000 107.8788240 2 3.1831000 15.3439560
If a second ECP is to be introduced, it should be placed right after the first one without any blank line! If a blank line is detected then the program will assume it’s done reading all ECPs and Basis Sets.
Finally, here is an example of a combination of both keywords. If a second ECP was needed then we’d place it at the end of the first one without a blank line. The molecule is any given chlorinated hydrocarbon (H, C and Cl atoms exclusively)
#P B3LYP/gen pseudo=read ADDITIONAL-KEYWORDS - blank line - 0 1 Molecular Coordinates - blank line - H 0 S 3 1.00 19.2384000 0.0328280 2.8987000 0.2312040 0.6535000 0.8172260 S 1 1.00 0.1776000 1.0000000 **** C 0 S 7 1.00 4233.0000000 0.0012200 634.9000000 0.0093420 146.1000000 0.0454520 42.5000000 0.1546570 14.1900000 0.3588660 5.1480000 0.4386320 1.9670000 0.1459180 S 2 1.00 5.1480000 -0.1683670 0.4962000 1.0600910 S 1 1.00 0.1533000 1.0000000 P 4 1.00 18.1600000 0.0185390 3.9860000 0.1154360 1.1430000 0.3861880 0.3594000 0.6401140 P 1 1.00 0.1146000 1.0000000 **** Cl 0 S 2 1.00 2.2310000 -0.4900589 0.4720000 1.2542684 S 1 1.00 0.1631000 1.0000000 P 2 1.00 6.2960000 -0.0635641 0.6333000 1.0141355 P 1 1.00 0.1819000 1.0000000 **** - blank line - CL 0 CL-ECP 2 10 d potential 5 1 94.8130000 -10.0000000 2 165.6440000 66.2729170 2 30.8317000 -28.9685950 2 10.5841000 -12.8663370 2 3.7704000 -1.7102170 s-d potential 5 0 128.8391000 3.0000000 1 120.3786000 12.8528510 2 63.5622000 275.6723980 2 18.0695000 115.6777120 2 3.8142000 35.0606090 p-d potential 6 0 216.5263000 5.0000000 1 46.5723000 7.4794860 2 147.4685000 613.0320000 2 48.9869000 280.8006850 2 13.2096000 107.8788240 2 3.1831000 15.3439560 - blank line -
If you like this post or found it useful please leave a comment, share it or just give it a like. It is as much fun to find out people is reading as it is finding the answer to ones questions in someone else’s blog 🙂
One of the most successful posts this blog has ever published was on certain nuances of the solvation calculations on PCM in G03. However there are some differences in the SCRF modules between G09 and G03 and I here present some of them as well as some tips to work with the new version.
The SCFVAC keyword used to calculate the Gibbs Solvation Energy change is no longer available. It is now replaced by DoVacuum which should be included in the SCRF options as SCRF=(DoVacuum,etc.). However, the absolute solvation energy now requires a gas-phase optimization along with a frequency calculation followed by the same calculations with the SCRF=SMD option in the desired solvent and with the appropriate variables.
Gaussian 03 used to calculate and report non-electrostatic contributions to the solvation energy, however they were not included in the total energy nor during optimization procedures. These non-electrosatic interactions are no longer calculated in the default. In order to include these terms during the SCF procedure, and to have them reported separately, the SCRF=SMD option should be used.
My previous post on PCM mentioned the usage of the options OFac=0.8 & RMin=0.5 as part of the additional input. These ‘magic numbers’ (I hate the term) were used to modify the way by which the overlapping spheres were treated in order to create the surface which in turn defined the cavity. G09 uses a new algorithm to make the overlaps generate a smoother surface. I recommend to use the default values before including ‘magic parameters’.
All the default values which G03 used can be retrieved with the G03Defaults keyword, but it is strongly suggested to use it only for comparison with calculations previously done with the older version.
As with some other so-called ‘white papers’ this post will be further improved as more information arises during my own calculations. Thanks for reading! Please comment/like/share this post, as well as others in the blog, if you found useful the information within.
Last week I had a presentation at the Chemistry Institute in which I talked about the research I’ve been doing during the last year. Here I insert a link to my prezi presentation in which I make an outline of the project’s scope and goals, so not many results available yet.
Unfortunately wordpress.com doesn’t allow embedding of flash files so here is the link to my presentation (in Spanish, sorry) directly to the prezi website.
I’ve been using Prezi lately (at www.prezi.com) and although I still can’t say I got the hang of it, I like the way this presentations flow way better than Powerpoint slides. With prezi you make a single slide in which all the information is contained and then you zoom back and forth (or just forth) through the topics you want to review. The presentation is created online, making it available for public/private use or online viewing; it may also be saved as a flash file which allows you to play it in almost any computer. If you know your way around mental maps, prezi is definitely for you!
This post will become updated continuously in order to include as much of these useful tips as I find along the way, so if you are interested please subscribe so you don’t have to visit often. Some are old or even a bit obsolete in terms of software versions but I still include them so they can all be gathered in one place; this is basically a summary of some common errors in the use of Gaussian (G03 and G09). Please feel free to comment on this post, sharing your tricks with various programs in whatever branch of computational chemistry you might be working on.
Vibrations calculated with Gaussian 09 can’t be visualized with GaussView 3.x
*Change the lines “Atom AN” with “Atom AN”. They look the same but there are TWO spaces in the first (this is the way G09 prints the output files) and in the second one there is only ONE space. The credit is entirely due to Jean Poully who posted this trick on the CCL a few months ago.
*Reader John Keller from Alaska suggests to use the formatted chk file (.fchk) to avoid this problem. As usual with formatted chk files make sure you format them in the same computer where you calculated it, or at least make sure both versions are compiled under the same architecture (32 or 64 bits).
Molecules appear flat on Molekel versions previous to 5.x
Molekel 4.x reads Gaussian98 files only. Gaussian03 and Gaussian09 files need to have the header corrected in the following part, right at the beginning of the file:
Gaussian 09: EM64L-G09RevA.02 11-Jun-2009
Gaussian 98: EM64L-G09RevA.02 11-Jun-2009
A rather old one and kind of a weird one too but it works and some of us (myself very much included) like the previous Molekel far more than the new 5.x version.
Electrostatic potentials can not be visualized with Molekel 4.x
This is an error that stems from the fact that atomic charges are not read by Molekel, so all electrostatic properties derived from them are not calculated nor visualized. The problem is that Gaussian also changed the header for the charges when going from Gaussian98 to Gaussian03. In Gaussian 98 charges are labeled as “Total Atomic Charges”. In G03 and G09, atomic charges are labeled “Mulliken Atomic Charges”. Therefore find and change this label for the latter. PLEASE NOTE that Molekel will look for the title “Total Atomic Charges” so if your calculation was a geometry optimization then you have to find the set of charges that correspond to the last step of it. Look for the last string and change that one; you may also change them all with the right edition option in your editor, Molekel will use the last one since that should be the one corresponding to the optimized structure.
I’ll post something else when I find/remember it. Thanks for reading!