Read Chapter 5 (readbook) and the handouts on Computational Determination of the Molecular Constants of HCl


Goals:

1 Construct the HCl  Morese potential from experimental measurements 
2 Construct the HCl  Morese potential from theoretical calculations

Why: Because this is chemistry nad molecules do not stay in one state, they REACT and CHANGE in the process of interaction/reaction!!

THIS IS THE MOST VALUABE INFORMATION FOR CHEMISTRY and chemists(this is all you did so far in chemistry - drawing reactions, but very  approximately and descriptively - vith "arrows" and "hills"...)

Note also that the energy you plot here is directly related to the electronic and nuclear structure and changes in the chemical reaction!!)

1

Start from your values from measurements ( and literature)
Morse potential needs ( you can do all in SI unit system or in cm-1, BUT BE CONSISTENT!!!!!   :)

re (equilibrium distance in m also used as l), use yours from measurement
De (binding energy or isociation energy or energy of the bond - this one can be obtained from spectroscopy or thermodynamic measurements--bomb, use the one from Herzberg table in eV and convert in J) ,

 and k,  ββ  (or β in your book) see p 165  k== 2 D ((ββ)^2)  , so you can use your k and De
from Herzberg

(see page 165 book)

Hertzber has eV for De (bond energy or binding energy)

See also Fig 1 p 437 textbook

De= 4.4336 eV

1 hartree = 4.359748 (10^(-18))  J    (check)

1 eV = 1.6  (10^(-19))  J    (check)

Dehclvib=7.10647  (10^(-19))  kg (m^2)  (s^(-2))

You must also correct for the zero vibrational energy = v/2

Evib=2.8535      (10^(-20))  kg (m^2)  (s^(-2))

Total "Denergy is now the sum, that is the bottom of the HO and Mores curves

Dehcl=Dehclvib+Evib

De=7.39182`*^-19

Back in eV that is

Dehcl/(1.6  (10^(-19)))

4.61988749 eV

Dev=4.61988

and if you use cm-1 (textbook) then in cm-1
1eV is 8065.5 cm-1

Dcm=4.61988749 *8065.5

CHECK PREVIOUS NUMBERS _ THEY ARE ONLY HERE TO SHOW YOU CONVERSIONS!!!

HERE  NOW YOU HAVE TO USE YOUR EXPERIMENTAL DATA FOR ALL AND YOUR CHOICE OF UNITS

c=2.99799 (10^8)  m/s

m1=1.00794*(1.660539*10^-27)

m2=35.453*(1.660539*10^-27)

μ=(m1 m2)/(m1+m2)

μ=1.62747 (10^(-27)) kg

νe =2.886 (10^3) /cm

Compare β with other values

Note νe=   (units constant) (Sqrt[(k/μ)])

THIS  c 2 Pi or c h indicates cm-1 or m-1 units

k=μ(2 π( c )( νe )  )^2

For example:  (USE YOUR VALUE HERE!!!!!!!)

νe =2.886 (10^5) /m

c=2.99799 10^8

νe=2.886 (10^5)

μ

h=6.626069*10^-34

Your lab textbook expression for beta (p 439 (10):

c=2.999 10^10

β1=(k/(2h c Dcm))^.5

(same result if De is in m^-1 and c is in m/s)

Units in SI system
:
J m^(-2)/(2  J s m s^(-1) m^(-1))

Clear[J,m,s]

Evaluate[(J m^(-2)/(2  J  s  m s^(-1) m^(-1)))^.5]

THIS  beta is in (m^-1) that comes from k (De and c must be in the same units )

1/β1

(AS SHOWN BEFORE)

PowerExpand[((1.8118576344890297`*^10 (Sqrt[((s^2)/(m^2))]))/s)]

USE YOUR re

Hertzberg re =

re =1.27455 10^-10

khcl=478 N (m^(-1))

khcl=478

ββ

De

Morse potential (see book),-it is a simple interpolation of the experimental values for diffferent vibrational levels of the ground excited state, alternative would be any analytical function such as polynom

 Use the expression that does not have additional constant, relative to zero, so that minimum energy is zero

So this is what book has.......

De

Plot[De (1-  Exp[-β1 (r-re)])^2,{r,0.1 (10^(-10)),4.1 (10^(-10))}]

Plot[Dcm  (1-  Exp[-β1  (r-re)])^2,{r,0.1 (10^(-10)),4.1 (10^(-10))}]

Plot[Dev  (1-  Exp[-β1  (r-re)])^2,{r,0.1 (10^(-10)),4.1 (10^(-10))}]

but, you can eliminate the constant (part that does not depend on variables
other than relative energy reference)- do it by calculating  
the square..^2

Then you will have two parallel curves shifted by a constant.....

Plot[{De  (- 2  Exp[-ββ  (r-re)]+ Exp[-ββ   2 (r-re)]),De  (1- 2  Exp[-ββ  (r-re)]+ Exp[-ββ   2 (r-re)])},{r,0.1 (10^(-10)),4.1 (10^(-10))}]

use the one that bottoms at zero- now we have Morse relative to the zero

func=De  (- 2  Exp[-ββ  (r-re)]+ Exp[-ββ   2 (r-re)])

Plot[De  (1- 2  Exp[-ββ  (r-re)]+ Exp[-ββ   2 (r-re)]),{r,0.1 (10^(-10)),4.1 (10^(-10))}]

THIS IS THE RESULT - YOUR EXPERIMENTAL/MEASURED energy dependence on the distance!!!!!

Q:Do you know what is the derivative of energy??

A: .........................................

You can also investigate properties of that Morse function as a function of variable k, De or ββ.

Example function of k

f[r_,k_]:=De  (1-  Exp[-ββ  (r-k re)])^2

Plot3D[f[r,k],{k,0.1,1.2},{r,0.1 (10^(-10)),4.1 (10^(-10))}]

You can also solve the QM solutions of a vibrator (harmonic oscillator) -see book

f[n_,x_]:= Abs[((1/Pi)^(1/4)HermiteH[n,x])/(E^(x^2/2)Sqrt[2^n n!])]^2

Plot[Evaluate@
  Append[Table[f[n,x]+n+1/2,{n,0,7}],x^2/2],{x,-4,4},Filling->Table[n->n-1/2,{n,1,8}]]

2 part 2 
Become familiar with the molecular electronic calculations (input , output of Gaussian program)

SetDirectory["C:\PC07f"]

hcl=Import["hcl.txt", "Table"]

datahcl={{1.,-460.262},{1.1,-460.316},{1.2,-460.339},{1.3,-460.343},{1.4,-460.336},{1.5,-460.324},{1.6,-460.308},{1.7,-460.291},{1.8,-460.275}}

ListPlot[datahcl]

<<LinearRegression`;

<<Statistics`LinearRegression`

dplothcl=ListPlot[datahcl,Ticks->{Automatic,Automatic},PlotRange->All,AxesLabel->{"r","energy in hartree"},PlotJoined->False]

approach 1, p 41-5 handouts

regressc=Regress[datahcl,{1,m,m^2,m^3},m];
Chop[regressc,10^(-6)]

regressc = Regress[datahcl, {1,m,m^2,m^3}, m,
RegressionReport ->
 {FitResiduals, SinglePredictionCITable,
 ParameterConfidenceRegion}]

or go to higher order (note ....max 8, for 9  points) note SE is large!!!!

regressc=Regress[datahcl,{1,m,m^2,m^3,m^4,m^5,m^6},m];
Chop[regressc,10^(-6)]

Conclusion - needs more points, subtraction of a large constant, units SI,and /  or better analytical form

BACK TO MORSE (now from the calculations)

datahcl={{1.,-460.262},{1.1,-460.316},{1.2,-460.339},{1.3,-460.343},{1.4,-460.336},{1.5,-460.324},{1.6,-460.308},{1.7,-460.291},{1.8,-460.275}}

Clear[d,b],const

r

d

b

const

re=1.275`*^-10

You can try directly.....

FindFit[datahcl,const+d(1-Exp[- b (r-re)])^2,{d,b,const},r]

not that good... 
or you can use one known parameter such as De or β and then find the other

Say you want top find calculated De :

De experiment (Herzberg)

4.4336*1.6 10^-19

.. using β ( to get to the correct region of fitting)

Start with all SI units ...

datahclsi={{1. 10 ^-10,-2.00306 10 ^-15},{1.1 10 ^-10,-2.0033 10 ^-15},{1.2 10 ^-10,-2.0034 10 ^-15},{1.3 10 ^-10,-2.00341  10 ^-15},{1.4 10 ^-10,-2.00338 10 ^-15},{1.5 10 ^-10,-2.00333 10^-15},{1.6 10 ^-10,-2.00326 10 ^-15},{1.7 10 ^-10,-2.00319 10^-15},{1.8 10 ^-10,-2.00312 10 ^-15}}

... and fit for De.... (use the const so that you do not depend on the relative zero, and also take -460 hartrees in Jouls so that you dpo not have  a large constant...

FindFit[datahclsi,const-460*4.359 10^(-18)+(d)(1-Exp[- (1.8118576344890297`*^10) (r-re)])^2,{d,const},r]

(8.143453757165877`*^-19)/(1.6 10^-19)

that is 5.1 eV  to be compared to De= 4.4336 eV

or in J

d->8.143453757165877`*^-19                           versus                    De=7.39182`*^-19

errordpercent = 100*(7.39182`*^-19-8.143453757165877`*^-19   )/(8.143453757165877`*^-19  )

CONCLUSION. more points better fit can bring us close to the experimental values for  potential !!!!

you can try also with :

{De  (- 2  Exp[-ββ  (r-re)]+ Exp[-ββ   2 (r-re)])