Chemical bonding/code

Matlab® code to plot figure 2 of chemical bonding (Heitler-London curves). See J. C. Slater, Quantum Theory of Molecules and Solids, vol. 1, McGraw-Hill, New York (1963), Table 3.1 and L. Pauling and E. Bright Wilson, Introduction to Quantum Mechanics, McGraw-Hill, New York (1935), pp. 342, 343.

function HL close all; R = [0.5 :0.1  :4.5]; [Eplus, Emin, E0] = VB(R);

plot(R, [Eplus; Emin; E0], 'linew', 2) line([R(1) R(end)], [0, 0], 'color', 'k', 'linew', 1.5, 'linest', '--')

xlabel('R/a_0', 'fonts', 16) ylabel('E/E_H', 'fonts', 16) ht = text(0.92,   2, 'E_-', 'fonts', 14 ); ht = text(0.85, 0.5, 'E_0', 'fonts', 14 ); hs = text(2.0, -0.36, 'E_+', 'fonts', 14 );

function [Eplus, Emin, H0] = VB(R) % Heitler-London equation for H2 % Slater, Molecules and Solids, vol. 1, p. 50 % Slater's equation for K' has a typo: w^2/3 --> w^3/3 a = 1.0; g = -psi(1);     % Euler's constant w = a*R;

S = exp(-w).*(1+w+w.^2/3); Sp = exp( w).*(1-w+w.^2/3);

J = 2*(-1./w + exp(-2*w).*(1+1./w)); K = -2*(exp(-w).*(1+w)); Jp = 2*(1./w - exp(-2*w).*(1./w + 11/8 + 3*w/4 + w.^2/6) ); Kp = (2/5)*(-exp(-2*w).*(-(25/8)+ (23*w)/4 +3*(w.^2) + (w.^3)/3) ...    +6./w.*(S.^2.*(g+log(w))+Sp.^2.*Ei(-4*w) -2*S.*Sp.*Ei(-2*w) ));

H0 = 2*J+Jp+2./R; H1 = 2*K.*S + Kp + 2*(S.^2)./R;

Eplus = (H0+H1)./(1+S.^2); Emin =  (H0-H1)./(1-S.^2); return function E = Ei(x) E = real(-expint(-x));