Help plotting the Lennard Jones Potential

I have the following code to plot the Lennard-Jones potential for Xe However, when I plot the code instead of getting the expected curve, I almost just get an "L" I tried changing the axis, but it still doesn't yield the proper curve.

%%1.1: Plotting the Leonard Jones Potential
eps = 1.77; %kJ/mol
sig = 4.10; %A
r = linspace(0.01*sig,6*sig,10000);

for i = 1:length(r)
    V(i) = 4*eps*((sig/r(i))^12-(sig/r(i))^6);
end

plot(V,r)
ylim([-1 1])

Answer 1

[Image Analyst]

Try it this way:

% 1.1: Plotting the Leonard Jones Potential
eps = 1.77; %kJ/mol
sig = 4.10; %A
r = linspace(0.01*sig,6*sig,10000);
for i = 1:length(r)
    V(i) = 4*eps*((sig/r(i))^12-(sig/r(i))^6);
end
semilogy(r, V, 'b-')
grid on;
% ylim([-1 1])

[Amanda Chun]

I'm not sure this is right, ideally my plot would look similar to this

[Image Analyst]

Look over the formula. I know nothing of that formula or whether you entered it correctly. I just plotted what you had so you could see it better.

Answer 2

[David Goodmanson]

Hi Amanda, try

plot(r,V)  instead of plot(V.r)
and
ylim([-2,1])

For a more Matlablike approach you could calculate the array V all at once as a function of the array r, which you have already defined as a vector.

V = 4*eps*((sig./r).^12-(sig./r).^6);

Using ./ and .^ means all the division and power calculations are done element-by-element with the array r.