Package Guide

PhysicalParticles.jl is useful for particle based scientific simulations

Installation

From the Julia REPL, type ] to enter the Pkg REPL mode and run

pkg> add PhysicalParticles

or add from git repository

pkg> add https://github.com/JuliaAstroSim/PhysicalParticles.jl

Test the package by

pkg> test PhysicalParticles

Basic Usage

Vectors

julia> using PhysicalParticles, UnitfulAstro
julia> a = PVector()PVector{Float64}(0.0, 0.0, 0.0)
julia> b = PVector(1.0u"m", 2.0u"m", 3.0u"m")PVector(1.0 m, 2.0 m, 3.0 m)
julia> c = PVector2D(u"m/s")PVector2D(0.0 m s^-1, 0.0 m s^-1)
julia> uconvert(u"m", PVector(1.0, 1.0, 1.0, u"km"))ERROR: UndefVarError: `uconvert` not defined in `Main`
Suggestion: check for spelling errors or missing imports.
Hint: a global variable of this name may be made accessible by importing Unitful in the current active module Main
julia> PVector(BigFloat)PVector{BigFloat}(0.0, 0.0, 0.0)
julia> PVector2D(BigInt, u"m")PVector2D(0 m, 0 m)
julia> PVector(1.0, 1.0) * imPVector2D{ComplexF64}(0.0 + 1.0im, 0.0 + 1.0im)
julia> b * 2.0u"s"PVector(2.0 m s, 4.0 m s, 6.0 m s)
julia> b + PVector(2.0, 2.0, 2.0, u"m") / 2PVector(2.0 m, 3.0 m, 4.0 m)
julia> norm(PVector2D(3.0f0,4.0f0))5.0f0
julia> normalize(PVector(3.0, 4.0))PVector2D{Float64}(0.6, 0.8)
julia> d = PVector(3u"kpc", 4u"kpc")PVector2D(3 kpc, 4 kpc)
julia> norm(d)5.0 kpc
julia> distance(PVector2D(0.0, 0.0), PVector2D(3.0, 4.0))5.0
julia> rotate(PVector(1.0, 0.0), 0.5pi)PVector2D{Float64}(6.123233995736766e-17, 1.0)
julia> rotate(PVector(1.0, 0.0, 0.0), 0.0, 0.0, 0.5pi)PVector{Float64}(6.123233995736766e-17, 1.0, 0.0)
julia> rotate_z(PVector(1.0, 0.0, 0.0), 90.0u"°")PVector{Float64}(0.0, 1.0, 0.0)
julia> rotate(PVector(1.0, 0.0, 0.0), 0.0, 0.0, 90.0u"°", PVector(-1.0, 0.0, 0.0))PVector{Float64}(-1.0, 2.0, 0.0)
julia> rotate(PVector(0.0, 1.0, 0.0), PVector(0.0, 1.0, 1.0), pi)PVector{Float64}(0.0, -2.220446049250313e-16, 0.9999999999999998)
julia> cartesian2cylinderial(PVector(sqrt(2), sqrt(2), 1.0, u"m"))(2.0 m, 0.7853981633974484, 1.0 m)
julia> cylinderial2cartesian(2.0u"m", pi/4, 1.0u"m")PVector(1.4142135623730951 m, 1.414213562373095 m, 1.0 m)
julia> cartesian2spherical(PVector(sqrt(0.5), sqrt(0.5), 1.0, u"m"))(1.4142135623730951 m, 0.7853981633974484, 0.7853981633974483)
julia> spherical2cartesian(sqrt(2)u"m", pi/4, pi/4)PVector(0.7071067811865476 m, 0.7071067811865475 m, 1.0000000000000002 m)
julia> zero(PVector{Float64})PVector{Float64}(0.0, 0.0, 0.0)
julia> iszero(PVector(u"m"))true
julia> isnan(PVector(NaN, NaN))true
julia> PVector2D(1.0, 1.0) ≈ PVector2D(1.0 + 1.0e-8, 1.0 + 1.0e-8)true
julia> ustrip(PVector(1.0, 1.0, 1.0, u"km"))ERROR: UndefVarError: `ustrip` not defined in `Main`
Suggestion: check for spelling errors or missing imports.
Hint: a global variable of this name may be made accessible by importing Unitful in the current active module Main
julia> ustrip(u"m", PVector(1.0, 1.0, 1.0, u"km"))ERROR: UndefVarError: `ustrip` not defined in `Main`
Suggestion: check for spelling errors or missing imports.
Hint: a global variable of this name may be made accessible by importing Unitful in the current active module Main

Particles

We provide 2D version for each type below, for example, the 2D version of Ball is Ball2D:

julia> Massless()Massless 0: Pos = PVector{Float64}(0.0, 0.0, 0.0), Vel = PVector{Float64}(0.0, 0.0, 0.0)
julia> Massless(PVector(0.0, 0.0, 0.0), PVector(), 1)Massless 1: Pos = PVector{Float64}(0.0, 0.0, 0.0), Vel = PVector{Float64}(0.0, 0.0, 0.0)
julia> Massless2D(uCGS)Massless 0: Pos = PVector2D(0.0 cm, 0.0 cm), Vel = PVector2D(0.0 cm s^-1, 0.0 cm s^-1)
julia> Ball()Ball 0: Pos = PVector{Float64}(0.0, 0.0, 0.0), Vel = PVector{Float64}(0.0, 0.0, 0.0), Acc = PVector{Float64}(0.0, 0.0, 0.0), Mass = 0.0
julia> Ball(PVector(0.0u"m", 0.0u"m", 0.0u"m"), PVector(u"m/s"), PVector(u"m/s^2"), 0.0u"kg", 1)Ball 1: Pos = PVector(0.0 m, 0.0 m, 0.0 m), Vel = PVector(0.0 m s^-1, 0.0 m s^-1, 0.0 m s^-1), Acc = PVector(0.0 m s^-2, 0.0 m s^-2, 0.0 m s^-2), Mass = 0.0 kg
julia> Star()Star 0 STAR: Pos = PVector{Float64}(0.0, 0.0, 0.0), Vel = PVector{Float64}(0.0, 0.0, 0.0), Acc = PVector{Float64}(0.0, 0.0, 0.0), Mass = 0.0, Ti_endstep = 0, Ti_begstep = 0, Potential = 0.0, OldAcc = 0.0
julia> SPHGas()ERROR: UndefVarError: `SPHGas` not defined in `Main`
Suggestion: check for spelling errors or missing imports.
julia> a = Star(uAstro)Star 0 STAR: Pos = PVector(0.0 kpc, 0.0 kpc, 0.0 kpc), Vel = PVector(0.0 kpc Gyr^-1, 0.0 kpc Gyr^-1, 0.0 kpc Gyr^-1), Acc = PVector(0.0 kpc Gyr^-2, 0.0 kpc Gyr^-2, 0.0 kpc Gyr^-2), Mass = 0.0 M⊙, Ti_endstep = 0, Ti_begstep = 0, Potential = 0.0 kpc^2 Gyr^-2, OldAcc = 0.0 kpc Gyr^-2
julia> b = SPHGas(uAstro)ERROR: UndefVarError: `SPHGas` not defined in `Main`
Suggestion: check for spelling errors or missing imports.
julia> distance(a,b)ERROR: MethodError: no method matching distance(::Star{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(kpc,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(kpc, Gyr^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-2, Unitful.FreeUnits{(kpc, Gyr^-2), 𝐋 𝐓^-2, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(M⊙,), 𝐌, nothing}}, Unitful.Quantity{Float64, 𝐋^2 𝐓^-2, Unitful.FreeUnits{(kpc^2, Gyr^-2), 𝐋^2 𝐓^-2, nothing}}, Int64}, ::PVector{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}})
The function `distance` exists, but no method is defined for this combination of argument types.

Closest candidates are:
  distance(::AbstractParticle, !Matched::AbstractParticle)
   @ PhysicalParticles ~/work/PhysicalParticles.jl/PhysicalParticles.jl/src/LinearAlgebra.jl:501
  distance(!Matched::AbstractPoint, ::AbstractPoint)
   @ PhysicalParticles ~/work/PhysicalParticles.jl/PhysicalParticles.jl/src/LinearAlgebra.jl:500

Mutate array of immutable particles

To simulate particle dynamics on GPU (bit types), and for better performance, particle data are stored in immutable sturcts. One of the shortcomings is that you have to allocate a new object when changing member data, however, memory operation is optimized when doing this in array:

julia> using BangBangERROR: ArgumentError: Package BangBang not found in current path.
- Run `import Pkg; Pkg.add("BangBang")` to install the BangBang package.
julia> p = Star()Star 0 STAR: Pos = PVector{Float64}(0.0, 0.0, 0.0), Vel = PVector{Float64}(0.0, 0.0, 0.0), Acc = PVector{Float64}(0.0, 0.0, 0.0), Mass = 0.0, Ti_endstep = 0, Ti_begstep = 0, Potential = 0.0, OldAcc = 0.0
julia> p = setproperties!!(p, Mass = 321.0)ERROR: UndefVarError: `setproperties!!` not defined in `Main`
Suggestion: check for spelling errors or missing imports.
Hint: a global variable of this name may be made accessible by importing BangBang in the current active module Main
julia> a = [Star()]1-element Vector{Star{Float64, Float64, Float64, Float64, Float64, Int64}}:
 Star 0 STAR: Pos = PVector{Float64}(0.0, 0.0, 0.0), Vel = PVector{Float64}(0.0, 0.0, 0.0), Acc = PVector{Float64}(0.0, 0.0, 0.0), Mass = 0.0, Ti_endstep = 0, Ti_begstep = 0, Potential = 0.0, OldAcc = 0.0
julia> a[1] = setproperties!!(a[1], Mass = 321.0) # In general, `Mass` is mutated right at its original memory address.ERROR: UndefVarError: `setproperties!!` not defined in `Main`
Suggestion: check for spelling errors or missing imports.
Hint: a global variable of this name may be made accessible by importing BangBang in the current active module Main

Random and Conversion

julia> p = rand(PVector{Float64}, 3)3-element Vector{PVector{Float64}}:
 PVector{Float64}(0.6868243848717273, 0.937955939458905, 0.677173937577704)
 PVector{Float64}(0.5057976095264569, 0.0724289577092524, 0.33455899778369014)
 PVector{Float64}(0.9178148576899107, 0.9929176467480426, 0.36329190817422574)
julia> pu = rand(PVector{Float64}, 3) * u"m"3-element Vector{PVector{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}}}:
 PVector(0.830344417615241 m, 0.17958472778512902 m, 0.11422169494351941 m)
 PVector(0.9999039584083748 m, 0.6175235661219828 m, 0.6277062707271623 m)
 PVector(0.6193178708409812 m, 0.8633908038453204 m, 0.9420866365531839 m)
julia> p_Ball = [Ball(uSI) for i=1:3]3-element Vector{Ball{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(m, s^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-2, Unitful.FreeUnits{(m, s^-2), 𝐋 𝐓^-2, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(kg,), 𝐌, nothing}}, Int64}}:
 Ball 0: Pos = PVector(0.0 m, 0.0 m, 0.0 m), Vel = PVector(0.0 m s^-1, 0.0 m s^-1, 0.0 m s^-1), Acc = PVector(0.0 m s^-2, 0.0 m s^-2, 0.0 m s^-2), Mass = 0.0 kg
 Ball 0: Pos = PVector(0.0 m, 0.0 m, 0.0 m), Vel = PVector(0.0 m s^-1, 0.0 m s^-1, 0.0 m s^-1), Acc = PVector(0.0 m s^-2, 0.0 m s^-2, 0.0 m s^-2), Mass = 0.0 kg
 Ball 0: Pos = PVector(0.0 m, 0.0 m, 0.0 m), Vel = PVector(0.0 m s^-1, 0.0 m s^-1, 0.0 m s^-1), Acc = PVector(0.0 m s^-2, 0.0 m s^-2, 0.0 m s^-2), Mass = 0.0 kg
julia> assign_particles(p_Ball, :Pos, pu)
julia> p_Ball3-element Vector{Ball{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(m, s^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-2, Unitful.FreeUnits{(m, s^-2), 𝐋 𝐓^-2, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(kg,), 𝐌, nothing}}, Int64}}:
 Ball 0: Pos = PVector(0.830344417615241 m, 0.17958472778512902 m, 0.11422169494351941 m), Vel = PVector(0.0 m s^-1, 0.0 m s^-1, 0.0 m s^-1), Acc = PVector(0.0 m s^-2, 0.0 m s^-2, 0.0 m s^-2), Mass = 0.0 kg
 Ball 0: Pos = PVector(0.9999039584083748 m, 0.6175235661219828 m, 0.6277062707271623 m), Vel = PVector(0.0 m s^-1, 0.0 m s^-1, 0.0 m s^-1), Acc = PVector(0.0 m s^-2, 0.0 m s^-2, 0.0 m s^-2), Mass = 0.0 kg
 Ball 0: Pos = PVector(0.6193178708409812 m, 0.8633908038453204 m, 0.9420866365531839 m), Vel = PVector(0.0 m s^-1, 0.0 m s^-1, 0.0 m s^-1), Acc = PVector(0.0 m s^-2, 0.0 m s^-2, 0.0 m s^-2), Mass = 0.0 kg
julia> pconvert([1.0, 2.0, 3.0])PVector{Float64}(1.0, 2.0, 3.0)
julia> pconvert([1.0u"m" 4.0u"m";
                 2.0u"m" 5.0u"m";
                 3.0u"m" 6.0u"m"])2-element StructArray(::Vector{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}}, ::Vector{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}}, ::Vector{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}}) with eltype PVector:
 PVector(1.0 m, 2.0 m, 3.0 m)
 PVector(4.0 m, 5.0 m, 6.0 m)

Extent

julia> p = [Ball(PVector(-1.0u"m", 1.0u"m", 1.0u"m"), PVector(u"m/s"), PVector(u"m/s^2"), 1.0u"kg", 1),
            Ball(PVector(3.0u"m", -3.0u"m", -3.0u"m"), PVector(u"m/s"), PVector(u"m/s^2"), 3000.0u"g", 2)]2-element Vector{Ball{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(m, s^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-2, Unitful.FreeUnits{(m, s^-2), 𝐋 𝐓^-2, nothing}}, M, Int64} where M}:
 Ball 1: Pos = PVector(-1.0 m, 1.0 m, 1.0 m), Vel = PVector(0.0 m s^-1, 0.0 m s^-1, 0.0 m s^-1), Acc = PVector(0.0 m s^-2, 0.0 m s^-2, 0.0 m s^-2), Mass = 1.0 kg
 Ball 2: Pos = PVector(3.0 m, -3.0 m, -3.0 m), Vel = PVector(0.0 m s^-1, 0.0 m s^-1, 0.0 m s^-1), Acc = PVector(0.0 m s^-2, 0.0 m s^-2, 0.0 m s^-2), Mass = 3000.0 g
julia> minimum_x(p)-1.0 m
julia> maximum_x(p)3.0 m
julia> center(p)PVector(1.0 m, -1.0 m, -1.0 m)
julia> pos_center(p)PVector(1.0 m, -1.0 m, -1.0 m)
julia> mass_center(p)PVector(2.0 m, -2.0 m, -2.0 m)
julia> median(p)ERROR: MethodError: no method matching isless(::Ball{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(m, s^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-2, Unitful.FreeUnits{(m, s^-2), 𝐋 𝐓^-2, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(g,), 𝐌, nothing}}, Int64}, ::Ball{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(m, s^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-2, Unitful.FreeUnits{(m, s^-2), 𝐋 𝐓^-2, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(kg,), 𝐌, nothing}}, Int64})
The function `isless` exists, but no method is defined for this combination of argument types.

Closest candidates are:
  isless(!Matched::Missing, ::Any)
   @ Base missing.jl:87
  isless(::Any, !Matched::Missing)
   @ Base missing.jl:88
  isless(!Matched::VersionNumber, !Matched::VersionNumber)
   @ Base version.jl:202
  ...
julia> extent(p)Extent: xMin = -1.0 m, xMax = 3.0 m, yMin = -3.0 m, yMax = 1.0 m, zMin = -3.0 m, zMax = 1.0 m, SideLength = 4.0 m, Center = PVector(1.0 m, -1.0 m, -1.0 m)