NBodySimulator.jl
Simulating systems of N interacting bodies.
This project is under development at the moment. The implementation of potential calculations is fairly experimental and has not been extensively verified yet. You can test simulations of different systems now, but be aware of possible changes in the future.
Installation
To install NBodySimulator.jl, use the Julia package manager:
using Pkg
Pkg.add("NBodySimulator")
Basic Components
There are three basic components required for any simulation of systems of N-bodies: bodies
, system
, and simulation
.
Bodies or Particles are the objects that will interact with each other and for which the equations of Newton's 2nd law are solved during the simulation process. Three parameters of a body are necessary, namely: initial location, initial velocity, and mass. MassBody
structure represents such particles:
using StaticArrays
r = SVector(0.0, 0.0, 0.0)
v = SVector(0.1, 0.2, 0.5)
mass = 1.25
body = MassBody(r, v, mass)
For simulation speed, it is advised to use static arrays.
A System covers the bodies and the necessary parameters for the correct simulation of interaction between particles. For example, to create an entity for a system of gravitationally interacting particles, one needs to use GravitationalSystem
constructor:
const G = 6.67e-11 # m^3/kg/s^2
system = GravitationalSystem(bodies, G)
Simulation is an entity determining the parameters of the experiment: time span of simulation, global physical constants, borders of the simulation cell, external magnetic or electric fields, etc. The required arguments for NBodySImulation
constructor are the system to be tested and the time span of the simulation.
tspan = (0.0, 10.0)
simulation = NBodySimulation(system, tspan)
There are different types of bodies, but they are just containers of particle parameters. The interaction and acceleration of particles are defined by the potentials or force fields.
Generating bodies
The package exports quite a useful function for placing similar particles in the nodes of a cubic cell, with their velocities distributed in accordance with the Maxwell–Boltzmann law:
N = 100 # number of bodies/particles
m = 1.0 # mass of each of them
v = 10.0 # mean velocity
L = 21.0 # size of the cell side
bodies = generate_bodies_in_cell_nodes(N, m, v, L)
Molecules for the SPC/Fw water model can be imported from a PDB file:
molecules = load_water_molecules_from_pdb("path_to_pdb_file.pdb")
Potentials
The potentials or force field determines the interaction of particles and, therefore, their acceleration.
There are several structures for basic physical interactions:
g_parameters = GravitationalParameters(G)
m_parameters = MagnetostaticParameters(μ_4π)
el_potential = ElectrostaticParameters(k, cutoff_radius)
jl_parameters = LennardJonesParameters(ϵ, σ, cutoff_radius)
spc_water_paramters = SPCFwParameters(rOH, ∠HOH, k_bond, k_angle)
The Lennard-Jones potential is used in molecular dynamics simulations to approximate interactions between neutral atoms or molecules. The SPC/Fw water model is used in water simulations. The meaning of arguments for SPCFwParameters
constructor will be clarified further in this documentation.
PotentialNBodySystem
structure represents systems with a custom set of potentials. In other words, the user determines the ways in which the particles are allowed to interact. One can pass the bodies and parameters of interaction potentials into that system. If the potential parameters are not set, the particles will move at constant velocities without acceleration during the simulation.
system = PotentialNBodySystem(bodies,
Dict(:gravitational => g_parameters,
:electrostatic => el_potential))
Custom Potential
There exists an example of a simulation of an N-body system at absolutely custom potential.
Here it is shown how to create custom acceleration functions using the tools of NBodySimulator.
First, it is necessary to create a structure for parameters for the custom potential.
struct CustomPotentialParameters <: PotentialParameters
a::AbstractFloat
end
Next, the acceleration function for the potential is required. The custom potential defined here creates a force acting on all the particles, proportionate to their masses. The first argument of the function determines the potential for which the acceleration should be calculated in this method.
import NBodySimulator.get_accelerating_function
function get_accelerating_function(p::CustomPotentialParameters,
simulation::NBodySimulation)
ms = get_masses(simulation.system)
(dv, u, v, t, i) -> begin
custom_accel = SVector(0.0, 0.0, p.a)
dv .= custom_accel * ms[i]
end
end
After the parameters and acceleration function have been created, one can instantiate a system of particles interacting with a set of potentials, which includes the just-created custom potential:
parameters = CustomPotentialParameters(-9.8)
system = PotentialNBodySystem(bodies, Dict(:custom_potential_params => parameters))
Gravitational Interaction
Using NBodySimulator, it is possible to simulate the gravitational interaction of celestial bodies. In fact, any structure for bodies can be used for the simulation of gravitational interaction since all those structures are required to have mass as one of their parameters:
body1 = MassBody(SVector(0.0, 1.0, 0.0), SVector(5.775e-6, 0.0, 0.0), 2.0)
body2 = MassBody(SVector(0.0, -1.0, 0.0), SVector(-5.775e-6, 0.0, 0.0), 2.0)
When solving a gravitational problem, one needs to specify the gravitational constant G.
G = 6.673e-11
Now we have enough parameters to create a GravitationalSystem object:
system = GravitationalSystem([body1, body2], G)
Usually, we solve an N-body problem for a certain period of time:
tspan = (0.0, 1111150.0)
The created objects determine the simulation we want to run:
simulation = NBodySimulation(system, tspan)
sim_result = run_simulation(simulation)
And, finally, we can animate our solution, showing two equal bodies rotating on the same orbit:
using Plots
animate(sim_result, "path_to_animated_particles.gif")
Electrostatic Interaction
Interaction between charged particles obeys Coulomb's law. The movement of such bodies can be simulated using ChargedParticle
and ChargedParticles
structures. The following example shows how to model two oppositely charged particles. If one body is more massive than another, it will be possible to observe the rotation of the light body around the heavy one without adjusting their positions in space. The constructor for the ChargedParticles
system requires bodies and Coulomb's constant k
to be passed as arguments.
r = 100.0 # m
q1 = 1e-3 # C
q2 = -1e-3 # C
m1 = 100.0 # kg
m2 = 0.1 # kg
v2 = sqrt(abs(k * q1 * q2 / m2 / r)) # m/s - using the centrifugal acceleration
t = 2 * pi * r / v2 # s - for one rotation
p1 = ChargedParticle(SVector(0.0, 0.0, 0.0), SVector(0.0, 0, 0.0), m1, q1)
p2 = ChargedParticle(SVector(100.0, 0.0, 0.0), SVector(0.0, v2, 0.0), m2, q2)
system = ChargedParticles([p1, p2], k)
simulation = NBodySimulation(system, (0.0, t))
sim_result = run_simulation(simulation)
Magnetic Interaction
An N-body system consisting of MagneticParticle
s can be used for the simulation of interacting magnetic dipoles, though such dipoles cannot rotate in space. Such a model can represent single domain particles interacting under the influence of a strong external magnetic field. To create a magnetic particle, one specifies its location in space, velocity, and the vector of its magnetic moment. The following code shows how we can construct an iron particle:
iron_dencity = 7800 # kg/m^3
magnetization_saturation = 1.2e6 # A/m
mass = 5e-6 # kg
r = SVector(-0.005, 0.0, 0.0) # m
v = SVector(0.0, 0.0, 0.0) # m/s
magnetic_moment = SVector(0.0, 0.0, magnetization_saturation * mass / iron_dencity) # A*m^2
p1 = MagneticParticle(r, v, mass, magnetic_moment)
For the second particle, we will use a shorter form:
p2 = MagneticParticle(SVector(0.005, 0.0, 0.0), SVector(0.0, 0.0, 0.0), 5e-6,
SVector(0.0, 0.0, 0.00077))
To calculate magnetic interactions properly, one should also specify the value for the constant μ0/4π or its substitute. Having created parameters for the magnetostatic potential, one can now instantiate a system of particles that should interact magnetically. For that purpose, we use PotentialNBodySystem
and pass particles and potential parameters as arguments.
parameters = MagnetostaticParameters(μ_4π)
system = PotentialNBodySystem([p1, p2], Dict(:magnetic => parameters))
simulation = NBodySimulation(system, (t1, t2))
sim_result = run_simulation(simulation, VelocityVerlet(), dt = τ)
Molecular Dynamics (MD)
NBodySimulator allows one to conduct molecular dynamic simulations for the Lennard-Jones liquids, the SPC/Fw model of water, and other molecular systems thanks to implementations of basic interaction potentials between atoms and molecules:
- Lennard-Jones
- electrostatic and magnetostatic
- harmonic bonds
- harmonic valence angle generated by pairs of bonds
The comprehensive examples of liquid argon and water simulations can be found in the examples
folder. Here only the basic principles of the molecular dynamics simulations using NBodySimulator are presented using liquid argon as a classical MD system for beginners. First, one needs to define the parameters of the simulation:
T = 120.0 # °K
T0 = 90.0 # °K
kb = 8.3144598e-3 # kJ/(K*mol)
ϵ = T * kb
σ = 0.34 # nm
ρ = 1374 / 1.6747# Da/nm^3
m = 39.95# Da
N = 216
L = (m * N / ρ)^(1 / 3)#10.229σ
R = 0.5 * L
v_dev = sqrt(kb * T / m)
bodies = generate_bodies_in_cell_nodes(N, m, v_dev, L)
τ = 0.5e-3 # ps or 1e-12 s
t1 = 0.0
t2 = 2000τ
Liquid argon consists of neutral molecules, so the Lennard-Jones potential runs their interaction:
parameters = LennardJonesParameters(ϵ, σ, R)
lj_system = PotentialNBodySystem(bodies, Dict(:lennard_jones => parameters));
Then, a thermostat and boundary conditions should be selected and instantiated:
thermostat = NoseHooverThermostat(T0, 200τ)
pbc = CubicPeriodicBoundaryConditions(L)
simulation = NBodySimulation(lj_system, (t1, t2), pbc, thermostat, kb);
result = run_simulation(simulation, VelocityVerlet(), dt = τ)
It is recommended to use CubicPeriodicBoundaryConditions
since cubic boxes are among the most popular boundary conditions in MD. There are different variants of the NBodySimulation
constructor for MD:
simulation = NBodySimulation(lj_system, (t1, t2));
simulation = NBodySimulation(lj_system, (t1, t2), pbc);
simulation = NBodySimulation(lj_system, (t1, t2), pbc, thermostat);
simulation = NBodySimulation(lj_system, (t1, t2), pbc, thermostat, kb);
The default boundary conditions are InfiniteBox
without any limits, the default thermostat is NullThermostat
(which does no thermostating), and the default Boltzmann constant kb
equals its value in SI, i.e., 1.38e-23 J/K.
Water Simulations
In NBodySImulator the SPC/Fw water model is implemented. For using this model, one has to specify parameters of the Lennard-Jones potential between the oxygen atoms of water molecules, parameters of the electrostatic potential for the corresponding interactions between atoms of different molecules, and parameters for harmonic potentials representing bonds between atoms and the valence angle made from bonds between hydrogen atoms and the oxygen atom.
bodies = generate_bodies_in_cell_nodes(N, mH2O, v, L)
jl_parameters = LennardJonesParameters(ϵOO, σOO, R)
e_parameters = ElectrostaticParameters(k, Rel)
spc_paramters = SPCFwParameters(rOH, ∠HOH, k_bond, k_angle)
water = WaterSPCFw(bodies, mH, mO, qH, qO, jl_parameters, e_parameters, spc_paramters);
For each water molecule here, rOH
is the equilibrium distance between a hydrogen atom and the oxygen atom, ∠HOH
denotes the equilibrium angle made of those two bonds, k_bond
and k_angle
are the elastic coefficients for the corresponding harmonic potentials. Further, one can pass the water system into the NBodySimulation
constructor as a usual system of N-bodies.
simulation = NBodySimulation(water, (t1, t2), pbc, kb);
Thermostats
Usually, during the simulation, a system is required to be at a particular temperature. NBodySimulator contains several thermostats for that purpose. Here the thermostating of liquid argon is presented, for thermostating of water, one can refer to this post.
Andersen Thermostat
τ = 0.5e-3 # timestep of integration and simulation
T0 = 90
ν = 0.05 / τ
thermostat = AndersenThermostat(90, ν)
Berendsen Thermostat
τB = 2000τ
thermostat = BerendsenThermostat(90, τB)
Nosé–Hoover Thermostat
τNH = 200τ
thermostat = NoseHooverThermostat(T0, 200τ)
Langevin Thermostat
γ = 10.0
thermostat = LangevinThermostat(90, γ)
Analyzing the Results of the Simulation
Once the simulation is completed, one can analyze the result and obtain some useful characteristics of the system. The function run_simulation
returns a structure containing the initial parameters of the simulation and the solution of the differential equation (DE) required for the description of the corresponding system of particles. There are different functions that help to interpret the solution of DEs into physical quantities. One of the main characteristics of a system during molecular dynamics simulations is its thermodynamic temperature. The value of the temperature at a particular time t
can be obtained by calling this function:
T = temperature(result, t)
Radial distribution functions
The RDF is another popular and essential characteristic of molecules or similar systems of particles. It shows the reciprocal location of particles averaged by the time of simulation.
(rs, grf) = rdf(result)
The dependence of grf
on rs
shows the radial distribution of particles at different distances from an average particle in a system. Here, the radial distribution function for the classic system of liquid argon is presented:
Mean Squared Displacement (MSD)
The MSD characteristic can be used to estimate the shift of particles from their initial positions.
(ts, dr2) = msd(result)
For a standard liquid argon system, the displacement grows with time:
Energy Functions
Energy is a highly important physical characteristic of a system. The module provides four functions to obtain it, though the total_energy
function just sums potential and kinetic energy:
e_init = initial_energy(simulation)
e_kin = kinetic_energy(result, t)
e_pot = potential_energy(result, t)
e_tot = total_energy(result, t)
Plotting Images
Using the tools of NBodySimulator, one can export the results of a simulation into a Protein Database File. VMD is a well-known tool for visualizing molecular dynamics, which can read data from PDB files.
save_to_pdb(result, "path_to_a_new_pdb_file.pdb")
In the future, it will be possible to export results via the FileIO interface and its save
function. Using Plots.jl, one can draw the positions of particles at any time of simulation or create an animation of moving particles, molecules of water:
using Plots
plot(result)
animate(result, "path_to_file.gif")
Makie.jl also has a recipe for plotting the results of N-body simulations. The example is presented in the documentation.
Contributing
Please refer to the SciML ColPrac: Contributor's Guide on Collaborative Practices for Community Packages for guidance on PRs, issues, and other matters relating to contributing to SciML.
See the SciML Style Guide for common coding practices and other style decisions.
There are a few community forums:
- The #diffeq-bridged and #sciml-bridged channels in the Julia Slack
- The #diffeq-bridged and #sciml-bridged channels in the Julia Zulip
- On the Julia Discourse forums
- See also SciML Community page
Reproducibility
The documentation of this SciML package was built using these direct dependencies,
Status `~/work/NBodySimulator.jl/NBodySimulator.jl/docs/Project.toml`
[e30172f5] Documenter v1.4.1
[0e6f8da7] NBodySimulator v1.10.0 `~/work/NBodySimulator.jl/NBodySimulator.jl`
and using this machine and Julia version.
Julia Version 1.10.3
Commit 0b4590a5507 (2024-04-30 10:59 UTC)
Build Info:
Official https://julialang.org/ release
Platform Info:
OS: Linux (x86_64-linux-gnu)
CPU: 4 × AMD EPYC 7763 64-Core Processor
WORD_SIZE: 64
LIBM: libopenlibm
LLVM: libLLVM-15.0.7 (ORCJIT, znver3)
Threads: 1 default, 0 interactive, 1 GC (on 4 virtual cores)
A more complete overview of all dependencies and their versions is also provided.
Status `~/work/NBodySimulator.jl/NBodySimulator.jl/docs/Manifest.toml`
⌅ [47edcb42] ADTypes v0.2.7
[a4c015fc] ANSIColoredPrinters v0.0.1
[1520ce14] AbstractTrees v0.4.5
[7d9f7c33] Accessors v0.1.36
⌅ [79e6a3ab] Adapt v3.7.2
[ec485272] ArnoldiMethod v0.4.0
⌃ [4fba245c] ArrayInterface v7.7.1
[4c555306] ArrayLayouts v1.9.2
[62783981] BitTwiddlingConvenienceFunctions v0.1.5
[2a0fbf3d] CPUSummary v0.2.4
[fb6a15b2] CloseOpenIntervals v0.1.12
[944b1d66] CodecZlib v0.7.4
[38540f10] CommonSolve v0.2.4
[bbf7d656] CommonSubexpressions v0.3.0
[34da2185] Compat v4.15.0
[a33af91c] CompositionsBase v0.1.2
[2569d6c7] ConcreteStructs v0.2.3
[187b0558] ConstructionBase v1.5.5
[adafc99b] CpuId v0.3.1
[9a962f9c] DataAPI v1.16.0
[864edb3b] DataStructures v0.18.20
[e2d170a0] DataValueInterfaces v1.0.0
[8bb1440f] DelimitedFiles v1.9.1
⌃ [2b5f629d] DiffEqBase v6.145.2
⌅ [459566f4] DiffEqCallbacks v2.36.1
[163ba53b] DiffResults v1.1.0
[b552c78f] DiffRules v1.15.1
[b4f34e82] Distances v0.10.11
[ffbed154] DocStringExtensions v0.9.3
[e30172f5] Documenter v1.4.1
[4e289a0a] EnumX v1.0.4
⌅ [f151be2c] EnzymeCore v0.6.6
[d4d017d3] ExponentialUtilities v1.26.1
[e2ba6199] ExprTools v0.1.10
[7034ab61] FastBroadcast v0.2.8
[9aa1b823] FastClosures v0.3.2
[29a986be] FastLapackInterface v2.0.3
[5789e2e9] FileIO v1.16.3
[1a297f60] FillArrays v1.11.0
[6a86dc24] FiniteDiff v2.23.1
[f6369f11] ForwardDiff v0.10.36
[069b7b12] FunctionWrappers v1.1.3
[77dc65aa] FunctionWrappersWrappers v0.1.3
[d9f16b24] Functors v0.4.10
⌃ [46192b85] GPUArraysCore v0.1.5
[c145ed77] GenericSchur v0.5.4
[d7ba0133] Git v1.3.1
[86223c79] Graphs v1.11.0
[3e5b6fbb] HostCPUFeatures v0.1.16
[b5f81e59] IOCapture v0.2.4
[615f187c] IfElse v0.1.1
[d25df0c9] Inflate v0.1.4
[18e54dd8] IntegerMathUtils v0.1.2
[3587e190] InverseFunctions v0.1.14
[92d709cd] IrrationalConstants v0.2.2
[82899510] IteratorInterfaceExtensions v1.0.0
[692b3bcd] JLLWrappers v1.5.0
[682c06a0] JSON v0.21.4
⌅ [ef3ab10e] KLU v0.4.1
[ba0b0d4f] Krylov v0.9.6
[73f95e8e] LatticeRules v0.0.1
[10f19ff3] LayoutPointers v0.1.15
[0e77f7df] LazilyInitializedFields v1.2.2
[5078a376] LazyArrays v1.10.0
[d3d80556] LineSearches v7.2.0
⌃ [7ed4a6bd] LinearSolve v2.22.1
[2ab3a3ac] LogExpFunctions v0.3.27
[bdcacae8] LoopVectorization v0.12.170
[1914dd2f] MacroTools v0.5.13
[d125e4d3] ManualMemory v0.1.8
[d0879d2d] MarkdownAST v0.1.2
[a3b82374] MatrixFactorizations v2.2.0
[bb5d69b7] MaybeInplace v0.1.2
[e1d29d7a] Missings v1.2.0
[46d2c3a1] MuladdMacro v0.2.4
[0e6f8da7] NBodySimulator v1.10.0 `~/work/NBodySimulator.jl/NBodySimulator.jl`
[d41bc354] NLSolversBase v7.8.3
[2774e3e8] NLsolve v4.5.1
[77ba4419] NaNMath v1.0.2
⌃ [8913a72c] NonlinearSolve v3.1.0
[6fe1bfb0] OffsetArrays v1.14.0
[bac558e1] OrderedCollections v1.6.3
⌃ [1dea7af3] OrdinaryDiffEq v6.66.0
[65ce6f38] PackageExtensionCompat v1.0.2
[d96e819e] Parameters v0.12.3
[69de0a69] Parsers v2.8.1
[f517fe37] Polyester v0.7.13
[1d0040c9] PolyesterWeave v0.2.1
[d236fae5] PreallocationTools v0.4.21
[aea7be01] PrecompileTools v1.2.1
[21216c6a] Preferences v1.4.3
[27ebfcd6] Primes v0.5.6
[8a4e6c94] QuasiMonteCarlo v0.3.3
[3cdcf5f2] RecipesBase v1.3.4
⌅ [731186ca] RecursiveArrayTools v2.38.10
[f2c3362d] RecursiveFactorization v0.2.23
[189a3867] Reexport v1.2.2
[2792f1a3] RegistryInstances v0.1.0
[ae029012] Requires v1.3.0
[7e49a35a] RuntimeGeneratedFunctions v0.5.13
[94e857df] SIMDTypes v0.1.0
[476501e8] SLEEFPirates v0.6.42
⌃ [0bca4576] SciMLBase v2.10.0
[c0aeaf25] SciMLOperators v0.3.8
[efcf1570] Setfield v1.1.1
⌃ [727e6d20] SimpleNonlinearSolve v1.4.0
[699a6c99] SimpleTraits v0.9.4
[ce78b400] SimpleUnPack v1.1.0
[ed01d8cd] Sobol v1.5.0
[a2af1166] SortingAlgorithms v1.2.1
⌃ [47a9eef4] SparseDiffTools v2.18.0
[e56a9233] Sparspak v0.3.9
[276daf66] SpecialFunctions v2.3.1
[aedffcd0] Static v0.8.10
[0d7ed370] StaticArrayInterface v1.5.0
[90137ffa] StaticArrays v1.9.3
[1e83bf80] StaticArraysCore v1.4.2
[82ae8749] StatsAPI v1.7.0
[2913bbd2] StatsBase v0.34.3
[7792a7ef] StrideArraysCore v0.5.6
⌅ [2efcf032] SymbolicIndexingInterface v0.2.2
[3783bdb8] TableTraits v1.0.1
[bd369af6] Tables v1.11.1
[8290d209] ThreadingUtilities v0.5.2
[3bb67fe8] TranscodingStreams v0.10.7
[d5829a12] TriangularSolve v0.2.0
[410a4b4d] Tricks v0.1.8
[781d530d] TruncatedStacktraces v1.4.0
[3a884ed6] UnPack v1.0.2
[3d5dd08c] VectorizationBase v0.21.67
[19fa3120] VertexSafeGraphs v0.2.0
[2e619515] Expat_jll v2.6.2+0
[f8c6e375] Git_jll v2.44.0+2
[1d5cc7b8] IntelOpenMP_jll v2024.1.0+0
[94ce4f54] Libiconv_jll v1.17.0+0
[856f044c] MKL_jll v2024.1.0+0
[458c3c95] OpenSSL_jll v3.0.13+1
[efe28fd5] OpenSpecFun_jll v0.5.5+0
[1317d2d5] oneTBB_jll v2021.12.0+0
[0dad84c5] ArgTools v1.1.1
[56f22d72] Artifacts
[2a0f44e3] Base64
[ade2ca70] Dates
[8ba89e20] Distributed
[f43a241f] Downloads v1.6.0
[7b1f6079] FileWatching
[9fa8497b] Future
[b77e0a4c] InteractiveUtils
[4af54fe1] LazyArtifacts
[b27032c2] LibCURL v0.6.4
[76f85450] LibGit2
[8f399da3] Libdl
[37e2e46d] LinearAlgebra
[56ddb016] Logging
[d6f4376e] Markdown
[a63ad114] Mmap
[ca575930] NetworkOptions v1.2.0
[44cfe95a] Pkg v1.10.0
[de0858da] Printf
[3fa0cd96] REPL
[9a3f8284] Random
[ea8e919c] SHA v0.7.0
[9e88b42a] Serialization
[1a1011a3] SharedArrays
[6462fe0b] Sockets
[2f01184e] SparseArrays v1.10.0
[10745b16] Statistics v1.10.0
[4607b0f0] SuiteSparse
[fa267f1f] TOML v1.0.3
[a4e569a6] Tar v1.10.0
[8dfed614] Test
[cf7118a7] UUIDs
[4ec0a83e] Unicode
[e66e0078] CompilerSupportLibraries_jll v1.1.1+0
[deac9b47] LibCURL_jll v8.4.0+0
[e37daf67] LibGit2_jll v1.6.4+0
[29816b5a] LibSSH2_jll v1.11.0+1
[c8ffd9c3] MbedTLS_jll v2.28.2+1
[14a3606d] MozillaCACerts_jll v2023.1.10
[4536629a] OpenBLAS_jll v0.3.23+4
[05823500] OpenLibm_jll v0.8.1+2
[efcefdf7] PCRE2_jll v10.42.0+1
[bea87d4a] SuiteSparse_jll v7.2.1+1
[83775a58] Zlib_jll v1.2.13+1
[8e850b90] libblastrampoline_jll v5.8.0+1
[8e850ede] nghttp2_jll v1.52.0+1
[3f19e933] p7zip_jll v17.4.0+2
Info Packages marked with ⌃ and ⌅ have new versions available. Those with ⌃ may be upgradable, but those with ⌅ are restricted by compatibility constraints from upgrading. To see why use `status --outdated -m`