Plotting Maxwellian

I have a list of positive and negative values and a single temperature. I am trying to plot the [Maxwell-Boltzmann Distribution][1] using the equation for particles moving in only one direction.

m_e = 9.11E-28 # electron mass [g]
k = 1.38E-16 # boltzmann constant [erg*K^-1]
v = range(1e10, -1e10, step=-1e8) # velocity [cm/s]

T_M = 1e6 # temperature of Maxwellian [K]

function Maxwellian(v_Max, T_Max)
    normal = (m_e/(2*pi*k*T_Max))^1.5
    exp_term = exp(-((m_e).*v_Max.*v_Max)/(3*k*T_Max))
    return normal*exp_term
end

# Initially comparing chosen distribution f_s to Maxwellian F_s
plot(v, Maxwellian.(v, T_M), label= L"F_s" * " (Maxwellian)")
xlabel!("velocity (cm/s)")
ylabel!("probability density")

However, when, plotting this, my whole function is 0:

[![enter image description here][2]][2]

I tested out if I wrote my function correctly by replacing return normal*exp_term with return exp_term (i.e. ignoring any normalization constants) and this seems to produce the distinct of the bell curve:

[![enter image description here][3]][3]

Yet, without the normalization constant, this will not preserve the area under the curve. I was wondering what may I be doing incorrectly with setting up my Maxwellian function and the constant in front of the exponential.
[1]: Maxwell–Boltzmann distribution - Wikipedia
[2]: https://i.stack.imgur.com/tYogq.png
[3]: https://i.stack.imgur.com/wyKeE.png

The data extrema are:

(2.9269998814999053e-123, 1.0769341115495682e-27)

Did you try manually setting the plot() keyword ylim accordingly?

Maybe there’s a choice of units that gives more reasonable size numbers? If so, Unitful and UnitfulRecipes could be of use.

Settings ylims does show the correct plot, as mentioned in my StackOverflow response to this.

Since we have the topic already, I might as well ask this here (I was planning to make a separate post): why are the default ylims of Plots.jl off by several orders of magnitude? Plotting this directly with julia> GR.plot(v, Maxwellian.(v, T_M)) automatically uses appropriate ylims, that lead to a visible plot.

Screenshot_20220214_115901

I think Plots has a lower threshold where it assumes whatever you’re plotting is “supposed” to be zeros. Note that eps(1.0) is about 1e-16, which could be the (somewhat flawed) reason for the specific choice. Either way, I would probably change the units to give a more intuitive scale.

Linking a related Plots.jl issue.

1 Like
using Plots, Unitful, UnitfulRecipes

m_e = 9.11E-28u"g"
k = 1.38E-16u"erg/K" # boltzmann constant [erg*K^-1]
v = range(1e10, -1e10, step=-1e8) .* 1u"cm/s" # velocity [cm/s]

T_M = 1e6u"K" # temperature of Maxwellian [K]

function Maxwellian(v_Max, T_Max)
    normal = (m_e/(2*pi*k*T_Max))^1.5
    exp_term = exp(-((m_e).*v_Max.*v_Max)/(3*k*T_Max))
    return normal*exp_term |> upreferred
end

results in a maximum at 1.077e-21 s^3 m^-3 which is not much more reasonable size number than 1e-26s^3 cm^-3. But in any case Unitful is your friend in doing any physical calculations.

P.S. s3/m3 appears wrong units

for particles moving in only one direction"

Good thing about SI units: 1e-21 s^3/m^3 = 1000 μs^3/m^3

1 Like