%% FINAL
%Define constants and potnetials; units: x in Angstrom, energy in cm-1,
%mass in atomic mass units and time in ps.
clear all;
c = 2.998e10;
mu = 1.6603e-27;
h = 6.626e-34;
hbar = h/2/pi;
A = 1e-10;

w1=200;
w2=200;
mass=1;
k1=1*mu*c*(2*pi*w1)^2/h*A^2;
k2=1*mu*c*(2*pi*w2)^2/h*A^2;
V1= @(x) k1/2*(x+0.5)^2;
V2= @(x) k2/2*(x-0.5)^2 - k2/2;
v12 = 100

e1 = @(x) ( V1(x) + V2(x) - sqrt( ( V1(x) - V2(x))^2 + 4*v12^2)) / 2
e2 = @(x) ( V1(x) + V2(x) + sqrt( ( V1(x) - V2(x))^2 + 4*v12^2)) / 2;
syms x
e11=  diff(e1, x)
syms x
e21= @(x) diff(e2, x)
e12= @(x) diff(e11, x)
e22= @(x) diff(e21, x);

fplot(V1,[-3,3], 'r')
hold on
fplot(V2,[-3,3], 'g')

hold off
figure
fplot(e1,[-3,3],'r')
hold on
fplot(e2,[-3,3],'g')

hold off
figure 
fplot(e11,[-3,3], 'r')
hold on
fplot(e21,[-3,3], 'g')

hold off
figure
fplot(e12, 'r', e22, 'g')
figure