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
Posted on November 11, 2009, in Computational Chemistry, Gaussian, NBO, White papers and tagged Computational Chemistry, Gaussian, NBO. Bookmark the permalink. 35 Comments.
I am currently investigating transition state structures of Diels-Alder and 1,3-dipolar cycloadditions. I have performed an NBO anaysis on a few of these calculations (QST2 or QST3 in G03). I have some problems interpreting the output but I am better after reading some posts here and elsewhere.
However, in my readings I have come to realize that one of the most important criteria used to identify bond development i.e. asynchronousity in the TS, is via analysis of Wiberg bond indices. Is there any literature, examples, tutorials, etc that allows one to teach themselves to read a Gaussian output and obtain the relevant information?
In general, I am concerned TS calculations and their structures. What insight can an NBO anlysis furnish oneself? Specifically, what keywords would one want to include in their calculation that afford insight into transition states?
References and tutorials in any capacity are welcome.
Best regards,
Glenn Abbott
Hello!
First let me appologize for my late response, I was on vacations and i completely neglected real life for a while.
In order to obtain Wiberg Indexes in your input you have to use the BNDIDX keyword in the $NBO section at the end of your input file.I don’t know of any tutorial about reading Wiberg bond indexes in Gaussian outputs, except for the NBO manual which I think you can consult directly at:
http://www.chem.wisc.edu/~nbo5/
Wiberg’s original reference is:
Wiberg, K.B. Tetrahedron, 24, (1968) 1083.
Now, reading wiberg’s indexes is very easy, you just have to search for the string ‘Wiberg’ in your input and that will take you to the Wiberg bond index (WBI matrix section in which you can read the WBI for any par of atoms by just crossing them in the matrix. The next place where the string is found is the “Wiberg Bond Index Total by atom” This section adds ALL of the wbi’s displayed above for each atom, which roughly resembles (or should resemble) the number of total bonds formed by each atom.
About the keywords suitable for TS analysis I think it depends on what kind of information you want to get. Since NBO is a Lewis reinterpretation of the wavefunction I think its use for TS structures would allow you to understand which electron pairs are being broken and possibly the Second Order Perturbation Analysis (obtained with the E2PERT keyword) would also yield the corresponding energy to such breaking process. I don’t have any experience with NBO and TS structures, I think it’s interesting and worth the research so please keep me posted about your results; I will sure do some trials of my own and will post them here.
Best regards and thanks for reading!
I forgot something…
If you are performing Geometry Optimizations you have to be careful what NBO analysis you are looking at, since it gets performed both at the starting geometry and at the optimized one; so of course you want to take a look at the final NBO analysis, or in other words you will find the ‘Wiberg’ string four times in your file (the firs two corresponding to the analysis of the starting geometry)
I hope this helps and as usual I’m available for any other thing you need
Best wishes!
Hi,
I have the following case and I’ll be grateful for any suggestion(s).
- I have crystal structure, in which there is interesting connection between two atoms and I want to study this bond by studying the NBO and to check the atomic orbitals that make this connection.
The problem is I tried many optimization methods and changed the bassis sets and played with the charge and multiplicity and in every time after optimization I lost the new bond. So, how to solve the optimization problem and is there any way to make single point calculation, followed by NBO?
- The best optimization method gave me 0.2A difference with the crystal structure however there is some modification in the bond conformation and that made me worried that the new structure is different. So, is that divination accepted?
- I’ll be happy to listen to any idea or suggestion.
Thank you very much
Best regards
Kem
Hi Kem!
) try calculating its properties without optimizing. Do not play with the charge and multiplicity, your molecule or ion has the number of electrons it has and that is that!
When you take a crystal structure and you optimize it you do it in a vacuum, that means you are taking any other interaction or force field out of the equation. Out of the crystal field your molecule might look quite different so maybe you should try to use a single point calculation unless you find too many imaginary frequencies. You may play with the basis set and method (a.k.a. level of theory) but playing with the charge and multiplicity makes no sense, you are removing/adding electrons to your structure! you can’t expect the wavefunction to behave the same way when you are in fact calculating different electronic configurations. If you still insist in optimizing it then try the option MaxStep=N where N relates to the size of the step as 0.01N bohr or radians (default N=30) so the optimization is more constrained. Another approach you might use is to freeze the bond; use opt=modredundant and after the molecule specification and before the NBO section type atomnumber1 atomnumber2 valueinangstroms F
The “F” will freeze that bond.
So in a nutshell, differences between optimized and crystal structures are to be expected. Maybe you are finding out that such a bond is an artifact produced within the context of the crystal field! (cant wait to read your paper
Have a nice day and thanks for reading
dear sir,
i am karthick research scholar in the field of vibrational spectroscopy. our team tried to get potential energy scan from Gaussian 03. but scan option in 03 requires input variables. but our doubt is “where the input variables to be attached?”. from your website we studied a lot in our research basics. we hope your favorable reply will help us a lot and best coordination between ours.
thanking you,
Hello Karthick!
I’m glad to know the site has been helpful to you. About your question, take a look at this other post of mine under this link.
The variables are placed at the end of the file after the molecule specification (coordinates), their syntax depends on the kind of scan you want to perform. I hope this answers your question but if it doesn’t then let me know and I will try to do it again.
Have a nice day!
Hello.
I want to study the nucleophilic attack of bromide to some complexes and oxidation of bromide to HOBr. I think maybe investigation about HOMO of Br- and LUMO of related complexes are useful. Lower gap energy between these two orbitals indicate the faster bromide oxidation reaction. Maybe survey the orbital size and spreading is helpful. Actually, I am not sure which kind of caclucaltions with Gaussian program are suitable for this situation.
Your guidance is my pleasure.
Best regards,
Mohammad Chahkandi
Inorganic and Analytical Chemistry II
Friedrich-Schiller-Universität Jena, Germany
Hi Mohammad!
I think you are on the right path, analyzing the gap between frontier orbitals is a good approach. What you need to do then is to run a population analysis, NBO is a robust and suitable method for your needs so give it a try! use pop=NBO in your route section and then proceed as described in this post. I also encourage you to run a Fukui indices analysis (find the post in this same blog)
Have a nice day!
dear sir,
I am saran. Im working in metal -organic molecule interaction studies. how to define the K-points for a particular system for eg.,Al(100) plane.
Hello Saran
So sorry for such a delayed response. I assume you have already constructed the 100 plane, if not check out the tutorials available at this site:
http://www.gaussian.com/g_tech/gv5ref/tutorials.htm
Then use in your input the keyword PBC=(NKPoint=Z) where Z is the number of k points you want to calculate.
I hope this helps! Thanks for reading
Концерты, песни, хиты, треки, видеоролики на Concertus.Info
Dear sir,
Iam Mrs.Rajamanikumaran research scholar in the field of vibrational spectroscopy.i wish to get your suggestions for my work.i wish to calculate hyperpolaraizability and NBO.i got 10 components for hyperpolarizability. output i got from gauss03w was Final packed hyperpolarizability:
K= 1 block:
1
1 0.170195D+03
K= 2 block:
1 2
1 0.150704D+03
2 -0.150383D+02 -0.229488D+03
K= 3 block:
1 2 3
1 -0.514043D-01
2 -0.215172D-01 -0.114820D-01
3 0.660429D+02 -0.636988D+02 -0.457934D-01 how to take beta components from this and how to convert it to esu.To calculate NBO i gave #rhf/6-311+G(d,p) pop=nbo in the route section.i didnt get the nbo output.please give me your suggestion regarding this.
your suggestions are very valuable for me and this site is very useful to me.
Thanking you
Mrs.Rajamanikumaran
Dear sir,
I am devi reseach scholar in the field of computational chemistry. I want to know how to calculate wibergs index, and how to give input file for finding wibergs index in guassian 03,and also i request you to give the details beta components of hyperpolarisibility in guassian 03 output file, for example beta value is 0.660429D+02 like that what is the value of D here.
Hi Devi
You can calculate the wiberg bond index by incluiding pop=NBORead in your route section and then at the end of the file leave a blank line followed by:
$NBO BNDIDX $END
This will calculate the Wiberg bond index.
About the D in 0.6604D+02 it is a way to write scientific notation, i.e., it means 10^2 hence 0.6604D+02 = 0.6604×10^2 = 66.04
A negative sign after the D would mean a negative exponent for example 12.34D-02 = 12.34×10^-2 = 0.1234
Hope this helps
Dear sir,
I am Devi, Thank you sir, your valuable suggestions helped me a lot. I have one more query in wibergsindex, how to write the input file in route section for finding wiberg index in guassian 03. I does not know about the route section basics. I request you to give the suggestions for me.
Thaking you,
Yours sincerely,
Devi..
Thank you so much for this helpful blog. I have run the NBO analysis, and have the data in hand. However it seems that the output from the NBO reorders my MOs from what I view them in with gaussview. So instead of NBO 1 being the first surface I view, it may be something like #17 or something like that.
I’m talking about output like this:
Natural Bond Orbitals (Summary):
Principal Delocalizations
NBO Occupancy Energy (geminal,vicinal,remote)
====================================================================================
Molecular unit 1 (C6H4)
1. BD ( 1) C 1 – C 2 1.99006 -0.89691 124(g),125(g),50(v),92(v)
128(v),134(v)
2. BD ( 2) C 1 – C 2 1.93227 -0.36735 128(v),134(v),124(g),125(g)
126(v),131(v),49(v),91(v)
130(r),133(r)
These two orbitals when I view them in gaussview actually correspond to #7 and #17. Is there any easy way to correlate these two without either comparing the energy values or the output from the “sorting of NBOs”?
Thanks again for this great blog!
Hi,
Thanks for your kind words about the blog.
It would seem like the orbitals were reordered but remember the basis was re-orthonormalized. Plus the NBO sorting uses a different criteria than the canonical MOs. I wouldn’t worry about this, for they only resemble in shape but not in energy.
Hope this helps!
Best wishes
Respected Sir
Natural Bond Orbitals (Summary) consist numerous values. From which how we interpret the calculation results. Sir any common methods is used for sorting these results.
Senior research fellow
Dept of chemistry
University of calicut
Great blog!
Quick question. We are trying to determine the relative interaction with a protein and a inhibitor. We have crystal structure data (from a .pdb) for the protein and are attempting to put different inhibitors in the pocket (which we assume we know where it is). We have optimized ONIOM calculations for several inhibitors and now want to look at NBO’s to distinguish what interactions are important.
However when i run a NBO with the following (obtained from your site, many thanks!)
#P opt rhf/6-31G(d,p) pop=(full,nboread) gfoldprint geom=connectivity
stuff
blank iine
$NBO PLOT $END
i am seeing the structure of our cut out piece of protein and inhibitor move with the different cycles. This seems alarming to me as we spent much time trying to get their optimal structure already. Is this the correct way to go about gleaming interactions between a protein and a inhibitor? And should i be worried about the movement of the molecules?
Hi Steve!
Thanks for your kind words about the blog!
This seems to me like a very thorough way of performing a docking analysis. It could get very expensive if you try to probe many different substrates, however if you can manage the calculations in reasonable amounts of time then go for it!
From what I see in your input you are asking gaussian to optimize the structure again. If you already did it with the ONIOM method then you can drop the ‘opt’ keyword and just run the NBO analysis of the intermolecular interactions and the obtained results will correspond to the ONIOM structure.
So in short. This is a correct, although probably unorthodox, way to probe the interactions between a substrate and a protein; I hope I can take a look at your paper when its done! About the movement of the molecules, this is not a problem but it will for sure take a long time to re-optimize the full system (which already took a long time optimizing with the ONIOM method, right?)
I hope this helps; if it doesn’t, please don’t hesitate to contact me again . Thanks for reading!
Hi,
This blog really is very helpful. I’m a grad student. I have been recently trying to use NBO analysis to try and get an indication of which orbitals are interacting the strongest. I do a lot of CASSCF and CASPT2 work to try and find transition states for electrocyclic isomerization reactions where I break sigma bonds. A lot of them are very bi-radical and I know that the HF wave function fails to describe them correctly. My question is do you think that the info I get from the perturbation theory analysis of the Fock matrix is relevant if the HF wave function is doing a poor job? I was hoping that this method would still have some qualitative use to help rationalize why several activations barriers were lower in energy than I expected. Any advice would be greatly appreciated.
Hi Jeff,
This is a very interesting question! I think this is a very good way to at least probe the space of those orbitals that might interact more strongly than others. I think that at least qualitatively you are safe and I see no problem in doing so if you save a lot of time on the CASSCF and CASPT2 calculations. Of course a systematic study would ensure this approach is sound and theoretically consistent but in the mean time I think you are ok. Nice work!
Thanks a lot for your kind words about the blog. I’m glad you’ve found it useful
Have a nice day!
Sir
I am having Gaussian 03, and i want to know how to calculate NBO, internal coordinates, and symmetry coordinates with model calculation
Hi Jeff,
This is a very interesting question! I think this is a very good way to at least probe the space of those orbitals that might interact more strongly than others. I think that at least qualitatively you are safe and I see no problem in doing so if you save a lot of time on the CASSCF and CASPT2 calculations. Of course a systematic study would ensure this approach is sound and theoretically consistent but in the mean time I think you are ok. Nice work!
Thanks a lot for your kind words about the blog. I’m glad you’ve found it useful
Have a nice day!
Hi sir
I performed a NBO analysis (with gaussian program) of an
open shell molecule employing the MP2 density (density=mp2) and
the following message appeared:
A bond orbital with an occupancy of – 0.00002 electrons was found!
Please check you input data.
Hi Zahra,
I would need to see your input data to be able to help you. An electron occupancy below zero is just wrong, I don’t know what is going on.
Have a nice day
Hi Joaquin
Thanks very much for your attention.my input is:# mp2/6-311++g(d,P) pop=(saveNBO) Gen 6D scf=(maxcyc=1000) density
and gave this error:A bond orbital with an occupancy of -0.00009 electrons was found!
Best regard for you
Hi Dr.joaquin
I performed an optimization on Cu with this input
# opt=(maxcycle=300,calcALL,steep,gdiis) b3lyp/6-31g(d,p) scrf=(solven
t=water,read,pcm) scf=(maxcyc=2000,vshift=800) integral(grid=ultrafine
)
and give this error
GPN1st: Unable to complete.
Error termination via Lnk1e in /share/apps/g09/l701.exe at Sat Apr 21 14:26
How I can solve this problem?
Have a nice day
Not to say that i know what i am doing. But it does appear that you have two max cycles at different numbers, 300 following opt and 200 following scf?
sir i run my compound using gaussian03,i got one negative frequency and 77 positive frequency.so how to modefy negative frequency.
Hi
negative frequency means that there is a transition state
my input #p b3lyp/6-31G(d,p) opt freq=raman
got one nagative freq and 77 positive freq
i did’t understand
Pingback: The Chemical blogspace (Cb) has ranked this blog! « Dr. Joaquin Barroso's Blog