[Physics] Plotting the maxwell-Boltzmann velocity distribution in matlab

thermodynamics

I'm trying to plot a maxwell-boltzman velocity distribution in matlab.

I have also asked this question at cross validated without much luck

The PDF is

f(v)=sqrt(m/2*pi*k*T) * exp(-m*v^2/2*k*T)

See this link under Distribution for the velocity vector

I have defined the following function:

(This simple function is split over so many lines for debugging purposes.)

function [ output_args ] = mypdf( V_i )
    m=2.18e-25;
    k=1.38e-23;
    T=500;

    intm = m/(2*pi*k*T);
    intmsqr = sqrt(intm);
    exponent = -m*V_i .*V_i;
    exponent = exponent/(2*k*T);
    exponent = exp(exponent);
    output_args = intmsqr*exponent;
end

I then ran this function for a variety of input velocities, saved them and plotted them like this:

>> velocities=-1000:1000;
>> results=mypdf(velocities);
>> plot(results)

To my horror I simply get a perfectly symmetrical distribution rather than classic maxwell-boltzman shape.

I also tried simply plotting the speeds of particles using another matlab function:

function [ output_args ] = mb_speed( V_i )
%UNTITLED Summary of this function goes here
%   Detailed explanation goes here

    m=2.18e-25;
    k=1.38e-23;
    T=500;

    term1 = m/(2*pi*k*T);

    term1 = term1* term1 * term1;

    term1 = sqrt(term1);

    term2 = 4*pi*V_i .* V_i;

    term3 = -m*V_i.*V_i;

    term3= term3/2*k*T;

    term3 = exp(term3);

    output_args = term1 * term2 .* term3;

end

I then plotted it using

>> speeds=0:1000;
>> r2=mb_speed(speeds);
>> plot(r2)

This just produces an exponential curve and not the classic shape.

What am I doing wrong in both of these cases?

Thanks.

Best Answer

First, this question has nothing to do with Matlab. Besides, if you think there is an error in your code implementation, this is most certainly the wrong stackexchange.

That said, the error is a physical one. You are confusing the distribution of velocity vectors $\vec{v}$, which is a 3D Gaussian (by that I mean a function from $\mathbb{R}^3$ to $\mathbb{R}$ that is the product of Gaussians in three Cartesian coordinates) with the 1D distribution of those vectors' magnitudes. The former had better be symmetric about 0 - do you expect that gas in thermal equilibrium is moving in one particular direction? This is in fact exactly what you plot (well, actually you plot a 1D slice of the 3D domain, essentially the distribution of $x$-velocities).

You need a $v^2$ factor in the distribution (and a few modified constants) to get the asymmetric thing you are looking for. This is the same $v^2$ (or $r^2$) that arises, for instance, when converting a spherically symmetric integrand from Cartesian to spherical coordinates. When you do that, by the way, don't try to plug in negative values - they don't make physical sense, since again it is a distribution of magnitudes.

In terms of the notation in that Wikipedia article, the distribution of vectors is $f_\mathbf{v}$, which is $f_v^3$, and you are plotting $f_v$. What you want to plot is what the article calls unsubscripted $f$ (well, one of the many distinct things it calls $f$): $$ f(v) = \sqrt{\left(\frac{m}{2\pi kT}\right)^3} 4\pi v^2 \exp\left(\frac{-mv^2}{2kT}\right). $$

Related Question