x=1:1:1000;
b=180;a=1000;
for i=1:length(x)
y(i)=sqrt(b/(2*pi*x(i)^3))*exp(-b*(x(i)-a)^2/(2*a^2*x(i)));
end
plot(x,y)