Category Archives: NBO

Physical Chemistry Meeting and CCIQS Symposium


Materials Research Institute at Morelia Michoacan (southern Mexico)

Materials Research Institute at Morelia Michoacan (southern Mexico)

So many events going on and so little time to blog about.
Two weeks ago, four members of this group traveled to Morelia in southern Mexico to present their research at the XIII Mexican Physical Chemistry Meeting. The next week after that, they all brought their posters back to Toluca for the internal symposium at CCIQS, where a masters student, María Eugenia, gave a small talk about her research project concerning photosynthesis in bacteria. Below, a short description of their projects is presented in order of seniority.

María Eugenia “Maru” Sandoval
Maru is working in photosynthesis of green sulfur bacteria. Her research deals with the excited states calculations at the Time Dependent DFT level for describing the first stages of photon interaction in antennae complexes of the photosystem II, namely the Fenna-Matthews-Olsen (FMO) complex, which was selected due to its relative structural simplicity over that of more evolved organisms. Maru also gave a talk at the internal Symposium back in Toluca the very next week where she got a positive feedback which will be used in her project.

Howard Diaz
One of the many strategies out there for treatment of HIV-1 infections is to block those proteins used to anchor the virus to a healthy cell. Sort of getting the virus’ hands busy so they can’t attach to a host. 60+ new compounds derived from thiourea were screened and assessed in their interactions with protein GP120, the protein to which the attachment is made, through docking and DFT calculations. Lead compounds are reported. It must be stressed that Howard got an award at CCIQS for having one of the best posters out of 70 in the entire symposium. Kudos and thanks to you, Howard! We now have some MD simulations in order.

Howard Diaz -  GP120 blockers for HIV-1

Howard Diaz – GP120 blockers for HIV-1

Guillermo “Memo” Caballero
His project has some nice philosophical implications if you ask me. Memo started as an experimental chemist and when he ran into a wall trying to obtain a pyridine from the non-aromatic analogue (glutarimide), he came to our group to run some calculations and find out how to force the aromatization process, or at least rationalize if it could be performed at all. Two mechanisms were proposed and now we know that even when the reaction should be quite exothermic, the reaction barriers are too high to be overcome by conventional methods. We now need to find a way to decrease those barriers (cue transition metal simulations). So in a way we are dealing here with the mechanism of a reaction that never happens (at least in an intramolecular way), leading to a reverse reductio ad absurdum reasoning – we assumed the reaction(s) did happen and we found out why is it impossible for them to happen.

No pic. available as of yet

Luis Enrique “Kike” Aguilar
Luis continues to work with calix(n)arenes, this lab’s first love, in drug delivery systems. He is working with two drugs at once: Bosutinib and Sorafenib, second generation drugs for the treatment of Chronic Myeloid Leukemia in cases were resistance to Imatinib has been developed. One of his main goals is to find a calixarene system which is able to discriminate between Bosutinib and pseudo-bosutinib, a commercial isomer which has incorrectly been available for a few years now.

"Kike" Aguilar -  Calix[n]arenes as drug carriers for Bosutinib and Sorafenib

“Kike” Aguilar – Calix[n]arenes as drug carriers for Bosutinib and Sorafenib

Thanks a lot for your efforts; they are paying good results to your ca

The boys from the lab

The boys from the lab

reers and the advancement of our research group. Now back to work, guys!

 

Advertisements

New paper in JACS


Well, I only contributed with the theoretical section by doing electronic structure calculations, so it isn’t really a paper we can ascribe to this particular lab, however it is really nice to see my name in JACS along such a prominent researcher as Prof. Chad Mirkin from Northwestern University, in a work closely related to my area of research interest as macrocyclic recognition agents.

In this manuscript, a calix[4]arene is allosterically opened and closed reversibly by coordinating different kinds of ligands to a platinum center linked to the macrocycle. (This approach has been referred to as the weak link approach.) I recently visited Northwestern and had a great time with José Mendez-Arroyo, the first author, who showed me around and opened the possibility for further work between our research groups.

(Ligands: Green = Chloride; Blue = Cyanide)

Closed, semi-open and fully open conformations; selectivity is modulated through cavity size. (Ligands: Green = Chloride; Blue = Cyanide)

Here at UNAM we calculated the interaction energies for the two guests that were successfully inserted into the cavity: N-methyl-pyridinium (Eint = 57.4 kcal/mol) and Pyridine-N-oxide (Eint = +200.0 kcal/mol). Below you can see the electrostatic potential mapped onto the electron density isosurface for one of the adducts. Relative orientation of the hosts within the cavity follows the expected (anti-) alignment of mutual dipole moments. At this level of theory, we could easily be inclined to assert that the most stable interaction is indeed the one from the semi-open compound and that this in turn is due to the fact that host and guest are packed closer together but there is also an orbital issue: Pyridine Oxide is a better electron acceptor than N-Me-pyridinium and when we take a closer look to the (Natural Bonding) orbitals interacting it becomes evident that a closer location does not necessarily yields a stronger interaction when the electron accepting power of the ligand is weaker (which is, in my opinion, both logic and at the same time a bit counterintuitive, yet fascinating, nonetheless).

Electrostatic potential mapped onto the electron density surface of one of the aducts under study

Electrostatic potential mapped onto the electron density surface of one of the adducts under study

All calculations were performed at the B97D/LANL2DZ level of theory with the use of Gaussian09 and NBO3.1 as provided within the former. Computing time at UNAM’s supercomputer known as ‘Miztli‘ is fully acknowledged.

The full citation follows:

A Multi-State, Allosterically-Regulated Molecular Receptor With Switchable Selectivity
Jose Mendez-Arroyo Joaquín Barroso-Flores §,Alejo M. Lifschitz Amy A. Sarjeant Charlotte L. Stern , and Chad A. Mirkin *

J. Am. Chem. Soc., Article ASAP
DOI: 10.1021/ja503506a
Publication Date (Web): July 9, 2014

 Thanks to José Mendez-Arroyo for contacting me and giving me the opportunity to collaborate with his research; I’m sure this is the first of many joint projects that will mutually benefit our groups. 

 

New paper in Computational and Theoretical Chemistry


I always get very happy to have a new paper out there! I find it exciting but most of all liberating since it makes you feel like your work is going somewhere but most of all that it is making its way ‘out there’; there is a strong feeling of validation, I guess.

Two very different families of calix[n]arenes (Fig 1) were tested as drug carriers for a very small molecule with a huge potential as a chemotherapeutic agent against Leukemia, namely, 3-phenyl-1H-[1]benzofuro[3,2-c]pyrazole a.k.a. GTP which has proven to be an effective in vitro Tyrosine Kinase III inhibitor. Having such a low molecular weight it is expected to have a very high excretion rate therefore the use of a carrier could increase its retention time and hence its activity. This time we considered n = 4, 5, 6 and 8 for the size of the cavities and R = -SO3H and -OEt as functional groups on the upper rim as to evaluate only a polar coordinating group and a non-polar non-coordinating one since GTP offers two H-bond acceptor sites and one H-bond donor one along the π electron density that could form π – π stacking interactions between the aromatic groups on GTP and the walls of the calixarene.

Fig 1. Calixarenes under study and their complexes with GTP

Fig 1. Calixarenes under study and their complexes with GTP

Once again calculations were carried out at the B97D/6-31G(d,p) level of theory along with molecular dynamics simulations for over 100 ns of production runs. NBO Deletion interaction energies were computed in order to discern which hosts formed the most stable complexes.

NBO Del interaction energies B97D/6-31G(d,p)

NBO Del interaction energies B97D/6-31G(d,p)

You may find a link to the ScienceDirect website for downloading the paper from this link. Last, but certainly not least, I’d like to thank all coauthors for their contributions and patience in getting this study published: Dr. Rodrigo Galindo-Murillo; Alberto Olmedo-Romero; Eduardo Cruz-Flores; Dr. Petronela M. Petrar and Prof. Dr. Kunsági-Máté Sándor. Thanks a lot for everything!

fig8

Donor and acceptor H-bond sites increases the probability of keeping the drug in place for a higher retention rate

Donor and acceptor H-bond sites increases the probability of keeping the drug in place for a higher retention rate

New paper in Journal of Chemical Theory and Computation


Happy new year to all my readers!

Having a new paper published is always a matter of happiness for this computational chemist but this time I’m excedingly excited about anouncing the publishing of a paper in the Journal of Chemical Theory and Computation, which is my highest ranked publication so far! It also establishes the consolidation of our research group at CCIQS as a solid and competitive group within the field of theoretical and computational chemistry. The title of our paper is “In Silico design of monomolecular drug carriers for the tyrosine kinase inhibitor drug Imatinib based on calix- and thiacalix[n]arene host molecules. A DFT and Molecular Dynamics study“.

In this article we aimed towards finding a suitable (thia-) calix[n]arene based drug delivery agent for the drug Imatinib (Gleevec by Novartis), which is a broadly used powerful Tyrosine Kinase III inhibitor used in the treatment of Chronic Myeloid Leukaemia and, to a lesser extent, Gastrointestinal Stromal Tumors; although Imatinib (IMB) exhibits a bioavailability close to 90% most of it is excreted, becomes bound to serum proteins or gets accumulated in other tissues such as the heart causing several undesired side effects which ultimately limit its use. By using a molecular capsule we can increase the molecular weight of the drug thus increasing its retention, and at the same time we can prevent Imatinib to bind, in its active form, to undesired proteins.

We suggested 36 different calix and thia-calix[n]arenes (CX) as possible candidates; IMB-CX complexes were manually docked and then optimized at the B97D/6-31G(d,p) level of theory; Stephan Grimme’s B97D functional was selected for its inclusion of dispersion terms, so important in describing π-π interactions. Intermolecular interaction energies were calculated under the Natural Bond Order approximation; a stable complex was needed but a too stable complex would never deliver its drug payload! This brings us to the next part of the study. A monomolecular drug delivery agent must be able to form a stable complex with the drug but it must also be able to release it. Molecular Dynamics simulations (+100 ns) and umbrella sampling methods were used to analyse the release of the drug into the aqueous media.

Optimized geometries for all complexes under study (B97D/6-31G*)

Optimized geometries for the 20 most stable complexes under study (B97D/6-31G*)

Potential Mean Force profiles for the four most stable complexes for position N1 and  N2 from the QM simulations are shown below (Red, complexes in the N1 position, blue, N2 position). These plots, derived from the MD simulations  give us an idea of the final destination of the drug respect of the calixarene carrier. In the next image, the three preferred structures (rotaxane-like; inside; released) for the final outcome of the delivery process are shown. The stability of the complexes was also assessed by calculating the values of ΔG binding through the use of the Poisson equations.

PMF for the most stable compounds

PMF for the most stable compounds

General MD simulation final structures

General MD simulation final structures

Thanks to my co-authors Maria Eugenia Sandoval-Salinas and Dr. Rodrigo Galindo-Murillo for their enormous contributions to this work; without their hard work and commitment to the project this paper wouldn’t have been possible.

First graduated student!


Today is truly a landmark in our lab because on this day, María Eugenia “Maru” Sandoval-Salinas has defended her thesis and has thus obtained her B. Sc. in Chemistry. She is the first student under my supervision to achieve this goal, and I hope it won’t be long until we get some more, although now the bar has been set quite high. For the time being, Maru is pursuing a career in the pharmaceutical industry but has every intention of coming back to the lab for her Masters degree; she has a reserved spot here with us at CCIQS.

Hard work pays off - We wait for you to come back for your Masters Degree!

Hard work pays off – We wait for you to come back for your Masters Degree!

Maru’s thesis deals mainly, but not exclusively, with calculating the interaction energies of calix- and thia-calix[n]arenes with the tyrosine kinase inhibitor Imatinib, which is widely used in the treatment of Chronic Myeloid Leukemia (CML), in order to rationally design a drug delivery agent for this drug. Her work is (a huge) part of an article currently under revision that I only wish had been published before her defense. Still, we await for that paper to be published in the next few weeks.

Throughout her stay at our lab, Maru was a dedicated student willing to learn new skills every time. As she replied today to one of the questions: “it’s not so much how many calculations I got right, but how many I got wrong!“. I find deep meaning in this sentence, perhaps deep enough as to consider it an aphorism, because indeed the more we try the more we fail, and the more we fail the more we learn and the closer we get to success.

Congratulations, Maru! I personally thank you for all the hard work invested in your thesis, all the long hours in front of the computer and your disposition to learn and work during the last 1.5 years. I’m certain you’ll find success in any venture you undertake; and I’m certain of it because you never stop trying.

Taking the Oath after being unanimously approved

Taking the Oath after being unanimously approved

Natural Bond Orbitals (NBO) Visualization with Chemcraft


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.

I’m glad to publicly thank Prof. Frank Weinhold’s gesture to include this blog in the bibliography section of the new NBO6.0 website under the NBO-Related Websites tab.

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.

Some Considerations:

> It’s possible to open a formated checkpoint file, containing the NBOs, in chemcraft.
Gaussian input:

%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!

Paper published in Canadian Journal of Chemistry


A new paper has been published and that is always good news. The paper entitled “Synthesis of new γ-lactones from preactivated monosubstituted pyrazines and TMS–ketene acetals” coauthored by Azucena Garduño-Alva, M. Carmen Ortega-Alfaro, José G. López-Cortés, Isabel Chávez, Joaquin Barroso-Flores, Rubén A. Toscano, Henri Rudler and Cecilio Álvarez-Toledano was a fruitful collaboration between several researchers from within the UNAM Institute of Chemistry and from other labs.

Therein, the lactone formation from pyrazines is analyzed, with some resulting orientations not quite in accordance with the common orientation patterns yield by electrondonor and electronwithrdawing groups. In order to assess the electronic structure of the intermediates and its influence on the resulting orientations, I performed a Fukui analysis based on the Natural Population formalism.

I will come back to this post and expand on the information once I get some more free time, thanks for your understanding. As usual the link to the paper can be located below and it is also available as a pdf file upon request to this author.
Thanks for reading!

Azucena Garduño-Alva, M. Carmen Ortega-Alfaro, José G. López-Cortés, Isabel Chávez, Joaquin Barroso-Flores, Rubén A. Toscano, Henri Rudler, Cecilio Álvarez-Toledano

Canadian Journal of Chemistry, 2012, 90(5): 469-482, 10.1139/v2012-016

Natural Bond Orbitals Deletion analysis (NBODel) in Gaussian03 & Gaussian09


The Natural Bond Orbitals Deletion analysis provides an excellent approach to the assessment of bonding energy within a single molecular fragment or between many. It deletes specific elements of the Fock matrix (this means it sets their values to 0.000) and then re-diagonalizes it in order to find the difference in energy respect to the original matrix. About nine different kinds of deletions are available, which will be briefly summarized in the following section.

One of the main strengths of the NBO derived methods is their almost complete basis set independence, which allows us to obtain comparable numbers under different levels of theory.

Both G03 and G09 use the NBO3.1 program. The 5.0 version is sold separately by their creators, namely prof. Frank Weinhold, who can be contacted through their website. It’s not available for geometry optimizations (gradients), some people insist on trying to get a different geometry by eliminating a certain interaction and that is just not possible directly with the NBODel method It is indeed possible to perform NBODeletions along with optimizations in G09 (Thanks to prof. Weinhold for his clarifying message) but there are some restrictions: molecular coordinates should be in Z-Matrix format and the number of variables to be optimized should not exceed 50; prof. Weinhold also recommends to use print=0 in the $NBO keylist in order to prevent the output files to become too big. Be sure to start with a proper geometry (close to the desired minimum) since given the nature of this analysis some apretiable geometry effects are to be expected.

The general syntax for its usage includes the string pop=NBODel in the route section of the GaussianX input file. Then, at the end of the file, the following is required:

--End of Input File--
--blank line--
$NBO $END
$DEL
Interactions to be deleted
$END

ENTIRE BLOCKS OF ATOMS

In this kind of deletion one is able to delete all the elements between specific groups of atoms, as if their orbitals (and hence their common Fock elements) did not overlap.

--End of Input File--
--blank line--
$NBO $END
$DEL
ZERO 2 ATOM BLOCKS
2 BY 3
1 2
3 4 5
3 BY 2
3 4 5
1 2
$END
--blank line--

The first line after $DEL indicates how many groups of atoms will be set to zero and the following lines indicate how many atoms belong to each group (i.e. the size of each block which in this case are 2 and 3, respectively). After this line the groups of atoms are listed, in this example all elements from atoms 1 and 2 with those of atoms 3, 4 and 5 will become zero. The next three lines are used for symmetry, so all the interactions from (1,2)->(3,4,5) are deleted along with (3,4,5)->(1,2)

DELETIONS BETWEEN ENTIRE MOLECULAR FRAGMENTS (Intermolecular deletions)

If we want to assess the interaction energy between two molecules, the previous method would consume a lot of time in declaring the size of each block with every atom of each molecule in it, plus there seems to be a limit to the size of the block. In this kind of deletion one is able to delete all the elements between two or more molecular fragments.

--End of Input File--
--blank line--
$NBO $END
$DEL
ZERO 2 DELOC FROM 1 TO 2 FROM 2 TO 1
$END
--blank line--

The delocalizations can also be calculated only in one direction (FROM 1 to 2), in the case above both interactions 1->2 and 2->1 have been deleted. The input for a trimer in which all three fragments interacted with each other would look like this:

ZERO 6 DELOC FROM 1 TO 2 FROM 2 TO 1 FROM 2 TO 3 FROM 3 TO 2 FROM 1 TO 3 FROM 3 TO 1

In short, the number of bilateral delocalizations to be deleted is equal to twice the number of edges in a graph depicting the intermolecular interactions (A post on topology in chemistry is now due).

Reading the output file

Almost at the very end of the output file the following section can be found:

>>>>>>>>>> Convergence criterion not met.
SCF Done:  E(RHF) =  -4728.57245403     A.U. after    2 cycles
Convg  =    0.2354D-03             -V/T =  2.0012
——————————————————————————
Energy of deletion :      -4728.572454034
Total SCF energy :      -4728.604640956
——————-
Energy change :          0.032187 a.u.,          20.198 kcal/mol
——————————————————————————

The warning about the convergence can be disregarded without any concern about the accuracy of the outcome and it will show in every $DEL calculation. The SCF energy displayed in the second line is the energy corresponding to the modified Fock Matrix, which is the same as the one labeled as Energy of deletion. The Total SCF energy corresponds to the original Fock Matrix; the difference between them is labeled as Energy change and the value is reported in both atomic units as well as kcal/mol.

Some common errors and possible solutions

–> Sometimes you get the following error message at the beginning of the calculation making it crash:

************************************************
** ERROR IN INITNF. NUMBER OF VARIABLES ( 57) **
**  INCORRECT (SHOULD BE BETWEEN 1 AND 50)  **
************************************************

I have found that changing the molecule specification section from Z-matrix to Cartesian coordinates, or vice versa, overcomes this difficulty. Also, if the Opt keyword appears in the route section the previous message will be shown. Opt is not available under the NBODel method (read the first paragraph for the proper correction).

–> Possible conflicts between NBODel and the usage of DFT methods:

In some revisions of Gaussian 03 there is a conflict when using NBODel and DFT methods. The IOp(5/48=10000) should be included to repair such conflict. This issue was solved in some revision of Gaussian 03 but I don’t know which, so try this if you have problems. Gaussian09 has taken care of the issue although still the usage of DFT to obtain NBODel calculations is not advised.

–> The following error is not self-explanatory:

NtrOpn-Old failed
Error termination via Lnk1e in ‘/../../path’

This particular error arises from the absence of the ‘$NBO $END’ line before the $DEL instruction. The previous line may or may not include additional keywords. If you are interested in computing some kind of deletion energy just leave the line as presented above in all previous examples. My guess is that the $DEL instruction does not calculate the corresponding NBO’s from which to make the deletion but it rather takes all the results from the $NBO instruction and works from there. Bottom line: don’t forget this line!

As with other posts tagged as ‘white papers’, this one will be updated and expanded every time new information is found. In the mean time, thanks to everyone for reading, commenting and rating, this keeps me going with the blog. Have you encountered problems with NBODel methods? share your experiences and solutions with the rest in the comments section.

Have a nice day!

Analyzing Natural Bond Orbitals (NBO) results from Gaussian ’03 & ’09


I’m posting this white paper “as is” but I will keep on coming back to update it; the thing is that I haven’t had the time to post anything else lately and if I don’t do it like this then I will never get it out. Please be kind on your comments and ratings!

For a more thorough discussion of NBO analysis please check out the references.

Requesting an NBO population analysis on the route section of a Gaussian job (via the simplest option pop=nbo) generates a lot of information on the output file which sometimes is not too straightforward to read. Please refer to previous posts in this blog (all under the NBO category) for more info on NBO calculations. The first part of this post deals with the analysis of the default information obtained from an NBO calculation; the second part deals with some useful/popular options of the NBO analysis which are controlled through the pop=nboread keyword available in NBO3.1 which is the version incorporated in Gaussian.

pop=nbo

Perhaps the most common issue when reading NBO results is making sense of the hybridization shown. Take the example below as taken from the NBO manual for Methylamine:

(Occupancy)   Bond orbital/ Coefficients/ Hybrids

————————————————————-

1. (1.99858)  BD  (1)  C 1 – N 2

(40.07%)    0.6330*  C1 s(21.71%)p 3.61(78.29%)

-0.0003  -0.4653  -0.0238  -0.8808  -0.0291

-0.0786  -0.0110  -0.0000  -0.0000

(50.93%)     0.7742* N  2  s(30.88%)p 2.24(69.12%)

-0.0001  -0.5557  0.0011  0.8302  0.0004

0.0443  -0.0098  0.0000 0.0000

The first line shows the occupancy (between 0.00 and 2.00 electrons), then the label of the NBO (BD = Bonding (2 centers); CR = Core (1 center); LP = Lone pair; RY = Rydberg; BD* = Antibond), after that a ‘serial number‘ which corresponds to the connectivity (or idealized bond index between the atoms, i.e., single double triple bond), and finally the atom(s) to which the NBO belongs.

The next line describes the natural atomic hybrids of which the NBO is composed, giving the percentage (100|cA|^2) of the NBO on each hybrid (40% the C hybrid orbital and 60% the N hybrid orbital), the polarization coefficient cA (Nitrogen is more electronegative than Carbon, hence the polarization coefficients look like 0.633 for C and 0.7742 for N, i.e. this NBO is more polarized towards N), the atom label (self explanatory) , and a hybrid label showing the sp^x composition: for the C hybrid orbital it is sp^3.61 while for the one on N is sp^2.24. The number in parenthesis is just the percentage in composition of the hybrid (75%p and 25%s corresponds to an ideal sp^3 hybridization) Now, for the tricky part: Sometimes you’ll get numbers larger than 3 for an sp^x hybrid or even s^yp^x hybrids (for y>1 and x>3), this has to do with the basis set employed and the number of functions used to describe each atomic orbital. (EXAMPLE SOON, IM LAZY)

Below the hybridization we find the set of coefficients that specify how the NHO is written as a linear combination of Natural Atomic Orbitals. The previous NHO for C has the largest coefficients from the second (-0.4653) and fourth (-0.8808) natural atomic orbitals, corresponding to a rough description like: NHOc = -0.4653(2s)c -0.8808(2px)c

pop=nboread

To control the options in the Natural Bond Orbital population analysis the line $NBO keywords $END must be included at the end of the input file (after the infamous blank line located after the molecule specification)

Wiberg bond indexes

This is a feature I use a lot when dealing with adducts with closed shell molecules. It is achieved through the BNDIDX keyword. Two indexes are obtained: The Wiberg bond index which is presented as a matrix; and the Wiberg bond index total which is the summation of all Wiberg indexes for every atom. The latter index roughly resembles the number of covalent bonds each atom forms. When analyzing this total index on a given atom, it should be compared to a well defined one in another part of the molecule just to make sure that all numbers are reflecting the same trends. Wiberg’s is not the only bond index provided by NBO analysis: the BOAO keyword generates bond index on the basis of the natural atomic orbitals.

Some common errors and possible solutions

->> Concerning NBO Deletion analysis

***** WARNING ***** The variational principle has been violated and the above deletion energy is invalid!!

This usually caused by one of the following:

1) Wavefunction symmetry breaking: Use the Nosymm option on the route section of your input file.

2) The use of DFT methods: Strong deletions often lead to densities that derail the density functional.  Safest is to use HF, where the variational theorem can be counted on, and
DFT artifacts are averted.  The FIXDM keyword (available only in NBO 5.0) may also correct some numerical problems associated with large basis sets

3) Deletion of some degenerated orbitals which breaks the symmetry of the wavefunction. Sometimes it could be that although there are no formal degenerate states accidental degeneracies might occur, deletion of which will cause this error.

->> Also on NBO Deletion Analysis

************************************************
** ERROR IN INITNF. NUMBER OF VARIABLES ( 57) **
**  INCORRECT (SHOULD BE BETWEEN 1 AND 50)  **
************************************************

This is one weird error. I have solved it by changing the molecule specification from Z-Matrix to Cartesian coordinates. Also, if you are running a DFT calculation in an older revision of Gaussian 03, try using the IOp(5/48=10000). Later revisions of Gaussian 03 and the more recent Gaussian 09 claim to have fixed this problem.

->> NBO Analysis

WARNING: Population inversion found on atom X#

This is not a problem! NBO lists orbitals within each atom according to their energy and their population afterwards. This warning only states that some low energy lying orbitals are less populated than others with higher energy. It cannot be “fixed” since this is a natural consequence of the orthonormalization process of the NBO generation.

References

NBO5.0 Home. Site of the creator of NBO and NBOView (Dr. Frank Weinhold)

FAQ‘s on NBO Home. I can’t recommend it enough!!!

Population Analysis Keyword in Gaussian

Natural Bond Orbitals (NBO) visualization


This is the second post on a series which will try to address common technical questions in computational chemistry that recursively appear on the CCL.

The Natural Bond Orbitals analysis is a powerful tool in population analysis calculations which is more robust than the traditional Mulliken approach, if for no better reason because its almost insensitive to the change of basis set while Mulliken’s P.A. is highly sensitive to basis set effects. Another advantage of the NBO analysis is that it provides a localized depiction of the electron density over a molecule, making it more related to chemists intuition. So far I have only worked with Gaussian 98, Gaussian 03 and only recently with Gaussian 09 in calculating NBO’s although it is possible also to perform them with GAMESS and the standalone NBO5.0 program created by Frank Weinhold. Visualizing them, however is never a straightforward process, and quite often we see more questions on the CCL than answers trying to address the matter. Most of the answers are concerned with what visualization programs to use but they seldom provide step by step instructions, furthermore most manuals are a bit cryptic about the procedure to plot this orbitals.

In Molekel 4.3

http://molekel.cscs.ch/wiki/pmwiki.php

Make sure that the route section in your input file includes the following options: pop=nboread; gfoldprint (in case you are using G03 or even G09) or gfprint (in case you are still running G98), e.g.:

#P opt rhf/6-31G(d,p) pop=(full,nboread) gfoldprint geom=connectivity

This calculation requests a geometry optimization followed by a Natural Bond Orbital population analysis (with keywords to be read) using the Restricted Hartree-Fock method with the split valence basis set 6-31G(d,p)

NOTE: I have my own issues and concerns about the use of DFT along with NBO but maybe that will be treated in another post.

At the end of the input file, after the infamous blank line type the following:

$NBO PLOT $END

You may include other keywords such as BNDIDX which generates a Wiberg bond index (order) matrix; or BOAO which generates the same matrix but in the Natural Atomic Orbitals basis. The PLOT option will generate a series of files with numerical extensions. The one you want to pay attention to is filename.47

If you are using Gaussian 03 or Gaussian 09 and STILL want to use Molekel 4.3 then you probably already know you have to change 03 for 98 on the header of the output file:

****************…              ****************
Gaussian 03:  x86…  –>       Gaussian 98:  x86-
2-Feb-2009                               2-Feb-2009
***************…                ****************

That aside, load your output file (filename.out or filename.log) on Molekel 4.3 as usual. Then go to Load -> nbo orb and load filename.47. Now, go to Compute -> Orbital and now select the Natural Bond Orbital you are interested in. This should do it! As an interesting exercise try computing the same orbitals (lets say HOMO and HOMO-1) with and without loading filename.47 in order to observe the difference between the shape of the MO’s and NBO’s. Molekel 4.3, though, is filled with bugs that will make it close unexpectedly, specially when running under Windows Vista. Sometimes the window closes because some sort of resolution problem, specially when taking snapshots (interestingly enough this happened to me when the background color was set to white) decrease the resolution of your monitor before taking the snapshot to prevent this problem. Some people complain about the look and feel of the latest molekel version so they stick to this old bugged one, so that is why I’m posting this method.

In Gasusview 3.0

This is the gaussview version I currently work with. When performing the NBO analysis on Gaussian (by the way, Gaussian 09 cannot be visualized with gaussview 3.x) use the savenbo option in the rout section, for example:

#P opt HF/6-31G(d) freq=noraman pop=(ful,nboread,savenbo) geom=connectivity

This will save the NBO coefficients into the checkpoint file. Load the output file with gaussview normally and then click on the Molecular Orbital icon (or go to Edit -> MO’s). This will open a new window with four tabs at the bottom of the molecule image. Click on the ‘New’ tab and load your checkpoint file. Gaussview will automatically format the chk file (this can cause some troubles when crossing architectures or platforms, so be careful to generate the proper formatted chk file!) Once load select the orbital you need to visualize and go to the ‘Visualize’ tab and click update; the orbital displayed is the Natural Bond Orbital. In this tab you can also adjust certain parameters like the isovalue (which is set to 0.02 by default) or the cube grid which controls how smooth the surface looks. UPDATE: According to John Keller from Alaska you can use this methodology when calculating with G09 and using Gview3.x, this also allows to visualize vibrations when this same software combination is used.

As usual, this post will be updated whenever I find some more useful information about the matter. Rate this post or leave a comment, just to know if you found it useful. Thanks!

%d bloggers like this: