ArgonMD.jl: Molecular Dynamics Simulation of Liquid Argon

ArgonMD.jl is a Julia package for molecular dynamics (MD) simulations of liquid argon in a cubic box with periodic boundary conditions (PBCs). Interactions between argon atoms are modeled using the Lennard-Jones (LJ) potential.

The package serves as a case study in performance engineering with Julia for scientific computing on supercomputers.

• Serial code optimization: Through improved type stability, elimination of unnecessary heap memory allocations, optimized memory layouts, exploitation of pairwise symmetry, cache blocking, SIMD vectorization, and minimized expensive operations such as division, square root, and exponential functions, ArgonMD.jl achieves 400x speedup over a baseline implementation on a single CPU core of an AMD EPYC 9655 processor. Benchmarks were performed on the Otus system at the Paderborn Center for Parallel Computing (PC2).
• Shared-memory parallelization: Using Julia multithreading, dynamic scheduling, and [cache blocking optimization](Details from Alex), ArgonMD.jl provides lock-free workload-balanced parallelization across all 192 CPU cores of a single compute node.
• Distributed-memory parallelization. Using MPI.jl, ArgonMD.jl combines a static balanced scheduler for inter-node workload distribution with dynamic scheduling for intra-node workload sharing. For an MD simulation of 8.4 million argon atoms, ArgonMD.jl achieves 99% parallel efficiency, demonstrating near-ideal strong scaling on 8 compute nodes (1536 CPU cores) of the Otus system at PC2.

ArgonMD.jl follows a modular design, where each implementation is organized as a Julia module within the package. To benchmark a specific implementation, run run_benchmark.jl with the module name, followed by the number of argon atoms, and optionally the number of MD steps for benchmarking (default = 10).

Examples for benchmarking serial, multithreaded, and MPI-parallel codes are shown below.

• To benchmark the serial code:

julia --project=<path_to_ArgonMD.jl> -- <path_to_ArgonMD.jl>/benchmark/run_benchmark.jl <module_name> <number_of_argon_atoms> [number_of_MD_steps]

See list of module names and short descriptions.

• To benchmark the multithreaded code:

export JULIA_NUM_THREADS=<number_of_worker_threads>,0 # number of Julia worker threads

julia --project=<path_to_ArgonMD.jl> -- <path_to_ArgonMD.jl>/benchmark/run_benchmark.jl <module_name> <number_of_argon_atoms> [number_of_MD_steps]

See list of module names and short descriptions.

• To benchmark the MPI code:

export JULIA_NUM_THREADS=<number_of_worker_threads>,0 # number of Julia worker threads per MPI process
# use `srun` as the MPI launcher with 8 MPI processes, each on ${JULIA_NUM_THREADS} CPU cores
srun --ntasks=8 --cpus-per-task=${JULIA_NUM_THREADS} julia --project=<path_to_ArgonMD.jl> -- <path_to_ArgonMD.jl>/benchmark/run_benchmark.jl MPI_01 <number_of_argon_atoms> [number_of_MD_steps]

The source code for serial implementations can be found in the "src/1-Serial" directory.

1. Module name: Serial_01
   • Source code: "01_baseline.jl"
   • Description: Baseline implementation for serial code
2. Module name: Serial_02
   • Source code: "02_type_stable.jl"
   • Description: Julia code with full type stability
3. Module name: Serial_03
   • Source code: "03_fixedsize_vectors.jl"
   • Description: Use fixed-size vectors to avoid frequent heap allocations
4. Module name: Serial_04
   • Source code: "04_SoA_memory.jl"
   • Description: Store data with SoA (Structure of Arrays) memory layout
5. Module name: Serial_05
   • Source code: "05_unique_pairs.jl"
   • Description: Exploit pairwise symmetry for computation
6. Module name: Serial_06
   • Source code: "06_cache_blocking.jl"
   • Description: Utilize cache blocking for memory accesses
7. Module name: Serial_07
   • Source code: "07_turbo_if_error.jl"
   • Description: Use @turbo to vectorize a loop containing if-condition
   • ❗Serial_07 is not included in the benchmark due to compile-time errors.
8. Module name: Serial_08
   • Source code: "08_simd_if.jl"
   • Description: Use @simd to vectorize a loop containing if-condition
9. Module name: Serial_09
   • Source code: "09_turbo_masked.jl"
   • Description: Use masked SIMD vectorization via @turbo
   • ❗Serial_09 is not included in the benchmark results due to low performance.
10. Module name: Serial_10

• Source code: "10_simd_masked.jl"
• Description: Use masked SIMD vectorization via @simd

11. Module name: Serial_11

• Source code: "11_simd_if_min_ops.jl"
• Description: Based on Serial_08 with minimized arithmetic operations

12. Module name: Serial_12

• Source code: "12_simd_masked_min_ops.jl"
• Description: Based on Serial_10 with minimized arithmetic operations

Benchmark calculations were performed for two simulation systems of different sizes:

• 50000 argon atoms in a cubic box (length = 133 Å) for 10 MD steps (Δt = 1 fs)
• 100000 argon atoms in a cubic box (length = 168 Å) for 3 MD steps (Δt = 1 fs)

This enables evaluation of the effect of system size, e.g., the number of atoms, while keeping the overall runtime approximately comparable.

• One single CPU core of AMD EPYC 9655 96-Core Processor
• Julia version 1.12.6 x86-64
• Results in the "PC2-results/Otus/1-Serial" directory

Note: Serial_07 and Serial_09 are not included due to compile-time errors and low performance.

Note: Serial_07 and Serial_09 are not included due to compile-time errors and low performance.

The source code for multithreaded implementations can be found in the "src/2-Multithreading" directory.

1. Module name: Multithreaded_01
   • Source code: "01_spawn_atomic_lock.jl"
   • Description: Multithreaded code using @spawn with atomic operations and lock to avoid data race conditions
2. Module name: Multithreaded_02
   • Source code: "02_threads_lock_free.jl"
   • Description: Multithreaded code using @threads without lock serialization
3. Module name: Multithreaded_03
   • Source code: "03_cyclic_scheduler.jl"
   • Description: Multithreaded code using @threads with a cyclic scheduler to balance workload
4. Module name: Multithreaded_04
   • Source code: "04_static_balanced_scheduler.jl"
   • Description: Multithreaded code using @threads with a static balanced scheduler
5. Module name: Multithreaded_05
   • Source code: "05_dynamic_balanced_scheduler.jl"
   • Description: Multithreaded code using @threads with a dynamic balanced scheduler
6. Module name: Multithreaded_06
   • Source code: "06_cyclic_scheduler_no_cache_blocking.jl"
   • Description: Based on 03 without cache blocking to demonstrate the adverse impact on performance and parallelization
7. Module name: Multithreaded_07
   • Source code: "07_static_balanced_scheduler_no_cache_blocking.jl"
   • Description: Based on 04 without cache blocking to demonstrate the adverse impact on performance and parallelization
8. Module name: Multithreaded_08
   • Source code: "08_dynamic_balanced_scheduler_no_cache_blocking.jl"
   • Description: Based on 05 without cache blocking to demonstrate the adverse impact on performance and parallelization

Benchmark calculations were performed for a simulation system of 524288 argon atoms in a cubic box (length = 292 Å) for 6 MD steps (Δt = 1 fs).

• Two AMD EPYC 9655 96-Core CPU sockets (192 CPU cores)
• Julia version 1.12.6 x86-64
• Results in the "PC2-results/Otus/2-Multithreading" directory

Dotted lines indicate the parallel efficiency normalized to a single-core CPU frequency of 4.5 GHz.

The source code for MPI + multithreaded parallelization can be found in the "src/3-MPI" directory.

1. Module name: MPI_01
   • Source code: "01_internode-static_intranode-dynamic.jl"
   • Description: Hybrid MPI + multithreaded parallelization, using the static scheduler for inter-node workload distribution and the dynamic scheduler for intra-node workload sharing

Benchmark calculations were performed for two simulation systems of different sizes:

• 4194304 argon atoms in a cubic box (length = 584 Å) for 16 MD steps (Δt = 1 fs)
• 8388608 argon atoms in a cubic box (length = 736 Å) for 4 MD steps (Δt = 1 fs)

This enables evaluation of the effect of system size, e.g., the number of atoms, while keeping the overall runtime approximately comparable.

• From 1 to 8 compute nodes (1536 CPU cores) of the Otus system
• Julia version 1.12.6 x86-64
• Results in the "PC2-results/Otus/3-MPI" directory

Molecular Modelling: Principles and Applications (2nd Edition) Andrew R. Leach, Pearson Education Limited 1996, 2001.