Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix Maxwell-Boltzmann Velocity function #1

Merged
merged 1 commit into from Oct 29, 2018

Conversation

Ellipse0934
Copy link
Contributor

Replacing by the Normal distribution function in Distributions.jl for performance and accuracy improvement.

Replacing by the Normal distribution function in Distributions.jl for performance and accuracy improvement.
@Ellipse0934
Copy link
Contributor Author

function maxwellboltzmann(mass::Real, T::Real) # Currently in Molly
    norm_dist = sum(rand(12)) - 6
    return abs(norm_dist) * sqrt(T / mass)
end
function NEWmaxwellboltzmann(mass::Real, t::Real) # Changed function
    return rand(Normal(0,sqrt(8.3144598*t/mass)))
end

Benchmarks: using @benchmark in BenchmarkTools

BenchmarkTools.Trial:  # Old Function 
  memory estimate:  176 bytes
  allocs estimate:  1
  --------------
  minimum time:     49.129 ns (0.00% GC)
  median time:      50.660 ns (0.00% GC)
  mean time:        61.214 ns (15.26% GC)
  maximum time:     44.944 μs (99.81% GC)
  --------------
  samples:          10000
  evals/sample:     991

BenchmarkTools.Trial: # New Function 
  memory estimate:  0 bytes
  allocs estimate:  0
  --------------
  minimum time:     4.364 ns (0.00% GC)
  median time:      4.789 ns (0.00% GC)
  mean time:        4.806 ns (0.00% GC)
  maximum time:     8.139 ns (0.00% GC)
  --------------
  samples:          10000
  evals/sample:     1000

Also, I think the current function's implementation might be wrong

return abs(norm_dist) * sqrt(T / mass)

By taking the absolute value of the normal distribution the resultant of the normalization of the three velocities does not end up in a Maxwell Distribution. (I began with studying MD with no Chemistry background 2 days ago so please forgive me if I made a huge blunder!)

screenshot from 2018-10-30 00-44-08
Blue: Current function and Green: Chi(3) scaled by a factor of sqrt(RT/M)
I also plotted against my function which seems to agree with the scaled Chi(3) and various other sources on the Internet.

@jgreener64
Copy link
Collaborator

Great work, thanks for fixing this and improving the implementation.

@jgreener64 jgreener64 merged commit dc6d731 into JuliaMolSim:master Oct 29, 2018
jgreener64 pushed a commit that referenced this pull request Feb 25, 2022
jgreener64 pushed a commit that referenced this pull request Apr 27, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

2 participants