How to plot the Free particle wave function

Hello!

I have the function:

And a plot from my book is :

As you see in the picture, as a function of x/a is ploted. I tried to program this in Matlab with optional parameters but I dont get the correct picture

My code is:

clc
clear all


a=10;
k=-100:1:100;

c1=pi*sqrt(2*a);

h=(6.626/(2*pi))*10^-34;
m=9.11*10^-31;
N=length(k);
xb=60; 
h=2*xb/(N-1);
x=-xb:h:xb;
t=m*a*a/h;
V=zeros(N,1);
S=0;
f=@(k,x)(1/c1)*((sin(k*a)./k)*exp(i*(k.*x-(h/(2*m))*t.*(k.^2))));

for i=1:length(k)-1
    c=f(k(i),x(i))+f(k(i+1),x(i));
    S=S+c*0.5;
    V(i)=S;
end
Ps=abs(V).^2;
plot(x/a,a*Ps)
grid on

The plot is :

How can I get the correct plot?!

Answer 1

[David Hill]

a=10;
c1=pi*sqrt(2*a);
h=(6.626/(2*pi))*10^-34;
m=9.11*10^-31;
X=-60:.1:60;
t=m*a^2/h;
V=zeros(1,201);
for j=1:length(X)
    f=@(k)((sin(k*a)./k).*exp(1i*(k.*X(j)-(h/(2*m))*t.*(k.^2))));
    V(j)=integral(f,-100,100);
end
plot(X/a,a*abs(V).^2);

[Pouyan Msgn]

This is the only message I get :

The integral may not exist, or it may be difficult to approximate numerically to the requested accuracy.

> In integralCalc/iterateScalarValued (line 372)

In integralCalc/vadapt (line 132)

In integralCalc (line 75)

In integral (line 88)

In QM3 (line 10)

[Pouyan Msgn]

It worked for me in intervall -10 to 10! thank you