I asked on Math.SE and was advised to try here instead.

I need to draw from a Maxwell-Boltzmann velocity distribution to initialise a molecular dynamics simulation. I have the PDF but I'm having difficulty finding the correct CDF so that I can make random draws from it.

The PDF I am using using is:

$$f(v)=\sqrt \frac{m}{2\pi kT} \times exp \left( \frac{-mv^2}{2kT}\right) $$

I am told that to find the CDF from the PDF we perform:

$$CDF(x)= \int_{-\infty}^x PDF(x) dx $$

After integrating $ f(v) $ I get:

$$ CDF(v)= \sqrt \frac{m}{2\pi kT} \times \left( \frac{\sqrt\pi\times erf \left( \frac{mv}{2\pi kT} \right) }{2\times \left( \frac{m}{2kT} \right)} \right) $$

$$CDF(v)= _{-\infty} ^{x} \left[ {\sqrt \frac{m}{2\pi kT} \times \left( \frac{\sqrt\pi\times erf \left( \frac{mv}{2\pi kT} \right) }{2\times \left( \frac{m}{2kT} \right)} \right)} \right] $$

  1. After I reach this point I am unable to proceed as I do not know how to evaluate something between $x$ and ${-\infty}$.

  2. I am also concerned that I have not done the integration correctly.

  3. I want to implement the CDF in C++ in the end so I can draw from it. Does anyone know if there will be a problem with doing this because of the erf, or will I be alright with this GSL implimentation ?

Thanks for your time.

@bryansis2010 on Math.SE says that I can evaluate in the range $x$ to $0$ instead of $-\infty$ as we do not drop below 0 Kelvin.

Would this then make the CDF:

$$ CDF(v)= \sqrt \frac{m}{2\pi kT} \times \left( \frac{\sqrt\pi\times erf \left( \frac{mv}{2\pi kT} \right) }{2\times \left( \frac{m}{2kT} \right)} \right) $$

as $erf(0)=0$

Is this correct?


closed as off topic by David Z Feb 17 '13 at 20:15

Questions on Physics Stack Exchange are expected to relate to physics within the scope defined by the community. Consider editing the question or leaving comments for improvement if you believe the question can be reworded to fit within the scope. Read more about reopening questions here. If this question can be reworded to fit the rules in the help center, please edit the question.

  • 2
    $\begingroup$ Hi RRs_Ghost - I suspect that whoever told to bring this question here was incorrect. I don't think it's on topic for us. Even though the function you're integrating comes from physics, you're still just asking how to do an integration, which itself is a pure math problem. I won't close this immediately so that other people have a chance to object if they would like to, but we may just be sending you back to Mathematics. (Also, for future reference: you shouldn't cross-post a question to multiple SE sites. Ask for it to be migrated if it's off topic on the first place you put it.) $\endgroup$ – David Z Feb 17 '13 at 16:34

The solution is to realise that that function is merely a Gaussian. In fact Each component of the velocity vector has a normal distribution with mean =0 and st-dev $\sqrt {kT/m}$.

All that is left to do at that point is to get the Gaussian CDF (well known) and sample from it, making sure to plug in our mas and temperature.

$$CDF(x)=\frac{1}{2}\times \left[ 1+ erf\left( \frac{x-\mu}{\sqrt{2{\sigma}^2}} \right) \right]$$

GSL implements gsl_cdf_ugaussian_P (double x) here


Not the answer you're looking for? Browse other questions tagged or ask your own question.