Computer Physics Communications

We present H-NESSi (The Hierarchical Non-Equilibrium Systems Simulation package), an open-source software package for solving the Kadanoff-Baym equations (KBE) of nonequilibrium Green’s function (NEGF) theory using hierarchical low-rank compression techniques. The simulation of correlated quantum systems out of equilibrium is severely limited by the cubic scaling in propagation time and quadratic memory growth associated with conventional two-time formulations. H-NESSi overcomes these limitations by combining high-order time-stepping schemes with hierarchical off-diagonal low-rank (HODLR) representations of the retarded and lesser Green’s functions, enabling controllable accuracy at substantially reduced computational cost and memory usage. Imaginary-time quantities are efficiently represented using the discrete Lehmann representation (DLR), enabling compact and accurate treatment of thermal initial states. The implementation supports multiorbital systems, adaptive singular value truncation, and both shared-memory (OpenMP) and distributed-memory (MPI) parallelization strategies suitable for large-scale lattice calculations. The workflow closely mirrors established NEGF frameworks while introducing compression transparently into the propagation procedure. Benchmark applications to driven superconductors within dynamical mean-field theory and to the two-dimensional Hubbard model demonstrate favorable scaling compared to conventional implementations, with asymptotic time complexity significantly below the cubic scaling of uncompressed approaches. H-NESSi thus enables long-time and large-system nonequilibrium simulations of correlated quantum materials, which were previously computationally prohibitive.

A major revision of the ChemNetworks software (originally published in the Journal of Computational Chemistry, 2014, 35, 495–505) is presented. While the original ChemNetworks provided foundational graph construction capabilities for chemical systems, it was limited to simple distance and 3-body angular edge criteria, was not designed for high-performance computing environments or real-time operation alongside running simulations. This release addresses these limitations through three core contributions. First, a recursive Z-matrix-based search algorithm is introduced that enables chemically intuitive, arbitrarily descriptive three-dimensional structure searches, supporting geometric, energetic, and logical criteria. Second, the DataSpaces data staging framework is incorporated as an optional I/O engine, enabling in-memory data exchange between ChemNetworks and running simulations that eliminates persistent storage bottlenecks and supports real-time graph construction and analysis. Third, a modular analysis framework leveraging the igraph library is introduced, providing a straightforward plugin architecture for community-contributed workflows. Benchmark results demonstrate linear scaling with system size and efficient MPI parallelization across up to 64 cores, with total computational complexity of O(N^R T/P), where N is the number of atoms, R is the Z-matrix depth, T is the number of timesteps, and P is the number of MPI processes.

We introduce SoliDualSPHysics, a novel open-source and GPU-accelerated software that extends DualSPHysics to enable the numerical simulation of hyperelastic, finite-strain plastic, and brittle fracture behavior in deformable solids within a unified smoothed particle hydrodynamics (SPH) software framework. The software implements a total Lagrangian formulation for solid mechanics that allows direct application of external loads and boundary conditions, enabling independent solid mechanics simulations. Brittle fracture is modeled through a phase-field approach coupled with SPH, allowing crack initiation, propagation, and branching under dynamic loading without explicit crack tracking, ad hoc crack-path criteria, or local refinement. The framework also supports user-defined mathematical expressions to prescribe time- and space-dependent quantities, complementing the solid and fracture extensions and enhancing flexibility across existing and future DualSPHysics applications. Leveraging DualSPHysics’ native CPU/GPU parallel architecture, the software achieves substantial computational acceleration for large-scale simulations, and the implementation is verified and validated against benchmark numerical problems and experimental data, demonstrating accuracy, robustness, and favorable scaling performance. Comprehensive implementation details and user documentation are provided to ensure reproducibility and to support further development by the community. The framework and source code are freely available through a public GitHub repository.

The MCPlas toolbox represents a collection of MATLAB functions for the automated generation of an equation-based fluid-Poisson model for non-thermal plasmas in the multiphysics simulation software COMSOL. Following the development of the new generation of the LXCat platform, all input data are prepared in a structured and interoperable JSON format and can be supplied and validated using existing JSON schemas. The toolbox includes fully transparent, editable MATLAB source code and offers an advanced description of electron transport in addition to commonly used approaches in the plasma modelling community. It supports one-dimensional and two-dimensional modelling geometries employing Cartesian, polar and cylindrical coordinate systems. MCPlas is tested on two reference cases: DC- and RF-driven low-pressure glow discharges in argon. Comparison of MCPlas results with results obtained by employing COMSOL’s Plasma Module verifies the reliability of the plasma model implemented by MCPlas and demonstrates the significance of electron transport treatment and boundary conditions applied in the toolbox. Using the same examples, the easy handling of complex reaction kinetic models in MCPlas and the reusability of its JSON input data across different modelling platforms are illustrated. This demonstrates that MCPlas provides a transparent and reproducible workflow for the simulation of non-thermal plasmas using COMSOL.

We present C-BerryTrans, a C++ code designed for first-principles calculations of Berry-curvature-driven transverse transport properties, namely the anomalous Hall conductivity (AHC) i.e., σ_{μν}^{AHC} and anomalous Nernst conductivity (ANC) i.e., α_{μν}^{ANC}. The code directly extracts eigenvalues and momentum-matrix elements from WIEN2k calculations and evaluates the Berry curvature (Ω) using a Kubo-like formalism, thereby avoiding the interpolation errors inherent in Wannier-based approaches. To ensure computational efficiency, C-BerryTrans parallelizes Ω evaluation over k-points using OpenMP and stores band-resolved curvature data in binary format, significantly reducing memory usage. This design enables rapid post-processing of AHC and ANC over a wide range of temperature (T) and chemical potential (ω) values in a single run. The code has been benchmarked on well-studied ferromagnetic materials- Fe, Fe3Ge, Pd, Fe3Al, and Co2FeAl. For Fe, the σ_{xy}^{AHC} is obtained to be  ∼ 775 ( ∼ 744) S/cm at 0 (300) K. In case of Fe3Ge, the calculated value of σ_{xy}^{AHC} is found to be 311 S/cm at 300 K. Nextly, for Co2FeAl, the magnitude of computed value of σ_{xy}^{AHC} at 2 K is found to be  ∼ 56 S/cm. Moving further, the room temperature magnitude of α_{xy}^{ANC} for Pd is obtained to be  ∼ 0.97 AK^{-1}m^{-1}. In case of Fe3Al, the maximum magnitude of α_{xy}^{ANC} for T ≤ 500 K is computed as  ∼ 2.83 AK^{-1}m^{-1}. Lastly, for Co2FeAl, the value of α_{xy}^{ANC} is obtained to be  ∼ 0.10 AK^{-1}m^{-1} at 300 K. These results show excellent agreement with previously reported data. With its accuracy, scalability, and user-friendly workflow, C-BerryTrans provides a powerful tool for exploring Ω-driven transport phenomena and is well suited for high-throughput materials discovery. The code further enables the evaluation of Ω-derived AHC/ANC contributions along user-defined high-symmetry k-point paths. This provides valuable microscopic insight into how specific band-structure features contribute to Ω-driven AHC/ANC. Additionally, the code is equipped with a visualization module that allows analysis of k-point contributions to AHC or ANC in any material. This further enhances its capability for exploring topological materials.

We present a Monte Carlo event generator specifically developed for the study of hadronic two-photon scattering events in two-photon scattering at electron-positron colliders. The code enables the generation of events with exact leading-order QED coupling and a flat phase space decay of the hadronic state into an arbitrary number of final state particles as selected by the user. Thus, this generator is well-suited for the use of partial wave analyses tools to study the two-photon production of higher-multiplicity final states across a wide range of energies and photon virtualities. Furthermore, the code integrates both experimental and theoretical inputs on the two-photon couplings of hadrons to simulate two-photon production processes. Motivated by the investigations of the BESIII collaboration, the final states π^{+}π^{-}, π^{0}π^{0}, π^{0}η, K^{+}K^{-}, K^{0}_{S}K^{0}_{S}, ηη, and f_1(1285) → ηπ^{+}π^{-} via a^{±}_{0}(980)π^{∓} and f_{0}(500)η are currently included. The code is sufficiently flexible to easily add additional final states as well as quickly change the already included channels.

Studying the interactions between energetic particles (EPs) and Alfvén eigenmodes (AEs) is essential for understanding alpha particle transport in burning plasma conditions within future fusion reactors. To numerically simulate EP-driven Alfvén instabilities, we have developed GMEC[1,2] (Gyrokinetic-Magnetohydrodynamics Energetic-particle Code) on the central processing unit (CPU) platform. However, long-time and large-scale simulations require higher grid resolution and more particles, demanding substantial computational resources. Therefore, applying more efficient numerical algorithms and optimization strategies on advanced computing architectures becomes critical. Computer scientists have designed the graphics processing unit (GPU) to implement single instruction multiple thread (SIMT) with fewer control units and higher throughput than CPUs. The GPU is optimized for data-parallel tasks and has emerged as a powerful platform for scientific computing. Based on NVIDIA’s Compute Unified Device Architecture (CUDA), this work presents cuGMEC, a high-performance gyrokinetic-MHD hybrid code developed from scratch for GPU platforms, incorporating newly added equations and physical terms. After testing on NVIDIA GPUs, the results show favorable scaling and acceleration, achieving up to 20 times speedup compared to the Intel Xeon Gold 6348 processor. Several benchmarks with theories and other codes have also been conducted to verify the gyrokinetic-MHD hybrid code.

Molecular dynamics is widely used to study various phenomena, such as diffusion, shock wave propagation, and plasma dynamics. A wide range of software packages supports the expanding scope of molecular dynamics applications. However, the quality of simulations depends on force field approximations, ranging from simple models to direct quantum solutions. Recently, machine learning approaches for constructing accurate interatomic potentials have received significant attention. In MDcraft, we integrate these advances into a scalable, physically accurate framework. MDcraft is a comprehensive, modern molecular dynamics platform. It offers a high-level Python API (Application Programming Interface) with a user-friendly, script-based interface. The core simulation algorithms are implemented in C++ to ensure robustness and computational efficiency. MDcraft is built for high-performance computing on modern clusters and supports dynamic domain decomposition and load balancing via the Message Passing Interface (MPI) for scalable parallelization. Additionally, MDcraft leverages multithreading within nodes through standard C++ parallelism, enabling efficient use of heterogeneous architectures. We demonstrate the code's capabilities through several examples, including the shock response in aluminum, the shock Hugoniot in argon, and the cold curve of copper.

Optical Bloch equations and rate equations serve as powerful tools to model light-matter interactions from textbook-like two-level atoms to the complex internal dynamics of molecules. A particular challenge in this context is posed by molecular laser cooling, where many dozens or hundreds of levels need to be taken into account for a comprehensive modeling. Here, we present MoleCool, a numerically efficient Python toolbox to implement and solve the corresponding differential equation systems. We illustrate both the capabilities of the toolbox and some of the intricacies of molecular laser cooling by educational examples, which range from simple Rabi oscillations to spontaneous and coherent cooling schemes for various currently studied or considered molecular species. This includes, in particular, a comprehensive modeling of laser cooling dynamics with full hyperfine structure resolution in radioactive radium monofluoride (RaF), as well as studies of other complex species such as barium monofluoride (BaF) and ytterbium monohydroxide (YbOH).

Density functional theory (DFT) remains the most widely used electronic structure method. Although exact in principle, in practice, it relies on approximations to the exchange-correlation (XC) functional, which is known to be a unique functional of the electron density. Despite 50 years of active research, existing XC approximations remain far from general purpose chemical accuracy of various thermochemical and materials properties. In that light, the inverse DFT problem, of finding the exact XC potential corresponding to an accurate groundstate density, offers an insightful tool to understand the nature of the XC functional as well as aid in the development of more accurate functionals. However, solving the inverse DFT problem is fraught with several numerical challenges, such as non-uniqueness or spurious oscillations in the solution and non-convergence. We present invDFT as an open-source framework to address the outstanding challenges in inverse DFT and computed XC potentials solely from a target density. We do so by use of a systematically convergent differential finite-element basis—higher-order finite-elements for the Kohn-Sham orbitals and linear finite-elements for the XC potential—which together render the inverse DFT problem well-posed. We also employ necessary asymptotic corrections to the target density to avoid any unphysical oscillations in the resulting XC potential. We also employ several numerical and high-performance computing (HPC) advances that affords both efficiency and parallel scalability, on CPU-GPU hybrid architectures. We demonstrate the accuracy and scalability of invDFT using accurate full-configuration interaction (FCI) densities as well as model densities, ranging up to 100 electrons and spanning both weakly and strongly correlated molecules.