<div><div>The Discrete Element Method (DEM) is widely used to simulate the mechanical behavior of granular materials across a broad range of applications and industrial domains. Particle shape is a key feature playing a crucial role for physics-fidelity of DEM simulations. However, accurately representing complex particle shapes within DEM frameworks presents significant challenges such as defining unambiguous contact normals or managing geometric singularities. Rigid particles are often modeled as convex polyhedra, which inherently suffer from ill-defined outward normal vectors at sharp edges and vertices. To represent non-convex geometries, these polyhedra must typically be combined, further increasing the computational and geometric complexity. In this work, we adopt an efficient and robust strategy to overcome these limitations by using <em>R</em>-shapes, defined as rounded-edge shapes, also known as sphero-polyhedra, obtained by sweeping a sphere of radius <em>R</em> along the edges and faces of a base polyhedral shape. This construction results in smooth surface transitions and circumvents common issues associated with traditional polygonal representations. This paper provides a detailed presentation of the implementation, structure, and advantages of <em>R</em>-shapes in DEM simulations. The proposed solutions are implemented in a fully open-source software package called <span>Rockable</span>, developed in <span>C++</span>, which integrates state-of-the-art numerical techniques and shared-memory parallelization for enhanced performance. Beyond the geometric modeling aspects, we also address several methodological challenges, including the treatment of contact elasticity and the numerical integration scheme. The combined contributions of this work offer a practical and efficient framework for simulating complex particle shapes in DEM with high physics fidelity and computational efficiency.<strong>Program summary</strong></div><div><em>Program Title:</em><span>Rockable</span></div><div><em>CPC Library link to program files:</em> (to be added by Technical Editor)</div><div><em>Developer’s repository link:</em></div><div><span><span>https://github.com/richefeu/rockable</span><svg><path></path></svg></span></div><div><em>Licensing provisions:</em> CeCILL-B</div><div><em>Programming language:</em>C<span>++</span>11</div><div><em>Supplementary material:</em></div><div><em>Nature of problem(approx. 50–250 words):</em></div><div>The open-source software <span>Rockable</span> addresses key challenges in simulating the mechanical behavior of granular materials using the Discrete Element Method (DEM), widely applied in both industrial applications and academic studies particularly where particle shape plays a critical role. Accurate modeling of the diversity of particle shapes in DEM remains non-trivial, due in part to ambiguities in defining contact normals. Rigid particles are often represented as convex polyhedra, which suffer from poorly defined
{"title":"Advanced strategies for discrete simulations with three-dimensional R-shapes in rockable framework","authors":"Vincent Richefeu , Gaël Combe , Lhassan Amarsid , Raphaël Prat , Jean-Mathieu Vanson , Saied Nezamabadi , Patrick Mutabaruka , Jean-Yves Delenne , Farhang Radjaï","doi":"10.1016/j.cpc.2025.109997","DOIUrl":"10.1016/j.cpc.2025.109997","url":null,"abstract":"<div><div>The Discrete Element Method (DEM) is widely used to simulate the mechanical behavior of granular materials across a broad range of applications and industrial domains. Particle shape is a key feature playing a crucial role for physics-fidelity of DEM simulations. However, accurately representing complex particle shapes within DEM frameworks presents significant challenges such as defining unambiguous contact normals or managing geometric singularities. Rigid particles are often modeled as convex polyhedra, which inherently suffer from ill-defined outward normal vectors at sharp edges and vertices. To represent non-convex geometries, these polyhedra must typically be combined, further increasing the computational and geometric complexity. In this work, we adopt an efficient and robust strategy to overcome these limitations by using <em>R</em>-shapes, defined as rounded-edge shapes, also known as sphero-polyhedra, obtained by sweeping a sphere of radius <em>R</em> along the edges and faces of a base polyhedral shape. This construction results in smooth surface transitions and circumvents common issues associated with traditional polygonal representations. This paper provides a detailed presentation of the implementation, structure, and advantages of <em>R</em>-shapes in DEM simulations. The proposed solutions are implemented in a fully open-source software package called <span>Rockable</span>, developed in <span>C++</span>, which integrates state-of-the-art numerical techniques and shared-memory parallelization for enhanced performance. Beyond the geometric modeling aspects, we also address several methodological challenges, including the treatment of contact elasticity and the numerical integration scheme. The combined contributions of this work offer a practical and efficient framework for simulating complex particle shapes in DEM with high physics fidelity and computational efficiency.<strong>Program summary</strong></div><div><em>Program Title:</em><span>Rockable</span></div><div><em>CPC Library link to program files:</em> (to be added by Technical Editor)</div><div><em>Developer’s repository link:</em></div><div><span><span>https://github.com/richefeu/rockable</span><svg><path></path></svg></span></div><div><em>Licensing provisions:</em> CeCILL-B</div><div><em>Programming language:</em>C<span>++</span>11</div><div><em>Supplementary material:</em></div><div><em>Nature of problem(approx. 50–250 words):</em></div><div>The open-source software <span>Rockable</span> addresses key challenges in simulating the mechanical behavior of granular materials using the Discrete Element Method (DEM), widely applied in both industrial applications and academic studies particularly where particle shape plays a critical role. Accurate modeling of the diversity of particle shapes in DEM remains non-trivial, due in part to ambiguities in defining contact normals. Rigid particles are often represented as convex polyhedra, which suffer from poorly defined ","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 109997"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145836598","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-03-01Epub Date: 2025-11-27DOI: 10.1016/j.cpc.2025.109957
Nicolás F. Barrera , Javiera Cabezas-Escares , Mònica Calatayud , Francisco Munoz , Tatiana Gómez , Carlos Cárdenas
FukuiGrid is a Python-based code that calculates Fukui functions and Fukui potentials in systems with periodic boundary conditions, making it a valuable tool for solid-state chemistry. It focuses on chemical reactivity descriptors from Conceptual Density-Functional Theory (c-DFT) and enables the calculation of Fukui functions through methods such as finite differences and interpolation. FukuiGrid addresses the challenges associated with periodic boundary conditions when calculating the electrostatic potential of a Fukui function (known as the Fukui potential) by integrating various corrections to alleviate the compensating background of charge. These corrections include the electrode approach and self-consistent potential correction as post-processing techniques. This package is compatible with VASP outputs and specifically designed to study the reactivity of surfaces and adsorbates. It generates surface reactivity maps and provides insights into adsorption site preferences, as well as regions prone to electron donation or withdrawal. FukuiGrid has been designed to make c-DFT easier for the surface chemistry community.
{"title":"FukuiGrid: A Python code for c-DFT in solid-state chemistry","authors":"Nicolás F. Barrera , Javiera Cabezas-Escares , Mònica Calatayud , Francisco Munoz , Tatiana Gómez , Carlos Cárdenas","doi":"10.1016/j.cpc.2025.109957","DOIUrl":"10.1016/j.cpc.2025.109957","url":null,"abstract":"<div><div>FukuiGrid is a Python-based code that calculates Fukui functions and Fukui potentials in systems with periodic boundary conditions, making it a valuable tool for solid-state chemistry. It focuses on chemical reactivity descriptors from Conceptual Density-Functional Theory (c-DFT) and enables the calculation of Fukui functions through methods such as finite differences and interpolation. FukuiGrid addresses the challenges associated with periodic boundary conditions when calculating the electrostatic potential of a Fukui function (known as the Fukui potential) by integrating various corrections to alleviate the compensating background of charge. These corrections include the electrode approach and self-consistent potential correction as post-processing techniques. This package is compatible with VASP outputs and specifically designed to study the reactivity of surfaces and adsorbates. It generates surface reactivity maps and provides insights into adsorption site preferences, as well as regions prone to electron donation or withdrawal. FukuiGrid has been designed to make c-DFT easier for the surface chemistry community.</div></div>","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 109957"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145681707","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-03-01Epub Date: 2025-11-05DOI: 10.1016/j.cpc.2025.109929
Isam Ben Soltane , Mahé Roy , Rémi Andre , Nicolas Bonod
The Singularity Expansion Method Parameter Optimizer - SEMPO - is a toolbox to extract the complex poles, zeros and residues of an arbitrary response function acquired along the real frequency axis. SEMPO allows to determine this full set of complex parameters of linear physical systems from their spectral responses only, without prior information about the system. The method leverages on the Singularity Expansion Method of the physical signal. This analytical expansion of the meromorphic function in the complex frequency plane motivates the use of an accuracy-driven improved version of the Cauchy method constrained by properties of physical systems, as well as an auto-differentiation-based optimization approach. Both approaches can be sequentially associated to provide highly accurate reconstructions of physical signals in large spectral windows. The performances of SEMPO are assessed and analysed in several configurations that include the dielectric permittivity of materials and the optical response spectra of various optical metasurfaces. SEMPO’s performances are thoroughly analyzed and benchmarked with other state-of-the-art methods to highlight its capability to retrieve the natural poles of a physical system. Program summary • Program title: Singularity Expansion Method Parameter Optimizer (SEMPO) • Library link to program files: • Developer’s repository link: https://doi.org/10.5281/zenodo.15210008 • Licensing provisions: Creative Commons Attribution 4.0 International • Programming language: Python • Nature of problem: Spectral functions of scattering coefficients of physical systems can be rigorously expanded through the Singularity Expansion Method. This expansion method permits to cast the spectral function through two equivalent forms: a sum of complex rational functions depending on complex poles and residues, and a factorized expression that depends on complex poles and zeros. Spectral functions are usually acquired along the real frequency axis, the challenge is thus to identify the set of complex poles, zeros and residues from these acquisitions at real frequencies. • Solution method: SEMPO is a numerical toolbox that aims at extracting the complex poles, zeros and residues from arbitrary spectral functions. SEMPO relies on both an improved algebraic Cauchy method associated with the factorized expression, and an open-source auto-differentiation method associated with the expansion in terms of poles and residues. Both methods are based on the Singularity expansion method and leverage the properties of physical systems and signals. The accuracy of SEMPO is thoroughly assessed and analysed in different configurations including the spectral responses of optical metasurfaces and the dielectric permittivity of optical materials.
{"title":"SEMPO - Retrieving Complex poles, residues and zeros from arbitrary real spectral responses","authors":"Isam Ben Soltane , Mahé Roy , Rémi Andre , Nicolas Bonod","doi":"10.1016/j.cpc.2025.109929","DOIUrl":"10.1016/j.cpc.2025.109929","url":null,"abstract":"<div><div>The Singularity Expansion Method Parameter Optimizer - SEMPO - is a toolbox to extract the complex poles, zeros and residues of an arbitrary response function acquired along the real frequency axis. SEMPO allows to determine this full set of complex parameters of linear physical systems from their spectral responses only, without prior information about the system. The method leverages on the Singularity Expansion Method of the physical signal. This analytical expansion of the meromorphic function in the complex frequency plane motivates the use of an accuracy-driven improved version of the Cauchy method constrained by properties of physical systems, as well as an auto-differentiation-based optimization approach. Both approaches can be sequentially associated to provide highly accurate reconstructions of physical signals in large spectral windows. The performances of SEMPO are assessed and analysed in several configurations that include the dielectric permittivity of materials and the optical response spectra of various optical metasurfaces. SEMPO’s performances are thoroughly analyzed and benchmarked with other state-of-the-art methods to highlight its capability to retrieve the natural poles of a physical system. <strong>Program summary</strong> • Program title: Singularity Expansion Method Parameter Optimizer (SEMPO) • Library link to program files: • Developer’s repository link: <span><span>https://doi.org/10.5281/zenodo.15210008</span><svg><path></path></svg></span> • Licensing provisions: Creative Commons Attribution 4.0 International • Programming language: Python • Nature of problem: Spectral functions of scattering coefficients of physical systems can be rigorously expanded through the Singularity Expansion Method. This expansion method permits to cast the spectral function through two equivalent forms: a sum of complex rational functions depending on complex poles and residues, and a factorized expression that depends on complex poles and zeros. Spectral functions are usually acquired along the real frequency axis, the challenge is thus to identify the set of complex poles, zeros and residues from these acquisitions at real frequencies. • Solution method: SEMPO is a numerical toolbox that aims at extracting the complex poles, zeros and residues from arbitrary spectral functions. SEMPO relies on both an improved algebraic Cauchy method associated with the factorized expression, and an open-source auto-differentiation method associated with the expansion in terms of poles and residues. Both methods are based on the Singularity expansion method and leverage the properties of physical systems and signals. The accuracy of SEMPO is thoroughly assessed and analysed in different configurations including the spectral responses of optical metasurfaces and the dielectric permittivity of optical materials.</div></div>","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 109929"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145681711","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-03-01Epub Date: 2025-11-17DOI: 10.1016/j.cpc.2025.109961
Jian Zhang , Ping Zhu , Chris C. Hegna
A hybrid spectral/finite-element code is developed to numerically solve the resistive finite-pressure magnetohydrodynamic equilibria without the necessity of postulating nested magnetic flux surfaces in the non-axisymmetric toroidal systems. The adopted approach integrates a hyperbolic parallel damping equation for pressure updating, along with a dynamic resistive relaxation for magnetic field. To address the non-axisymmetry in toroidal geometry, a pseudo flux mapping is employed to relate the axisymmetric computational domain to the physical domain. On the computational mesh, an isoparametric C1-continuous triangular element is utilized to discretize the poloidal plane, which is complemented with a Fourier decomposition in the toroidal direction. The versatility of the code is demonstrated through its application to several different non-axisymmetric toroidal systems, including the inherently three-dimensional equilibria in stellarators, the helical-core equilibrium states in tokamak plasmas, and the quasi-single-helicity states in a reversed-field pinch.
{"title":"Numerical solutions of resistive finite-pressure magnetohydrodynamic equilibria for stellarator and non-axisymmetric toroidal plasmas","authors":"Jian Zhang , Ping Zhu , Chris C. Hegna","doi":"10.1016/j.cpc.2025.109961","DOIUrl":"10.1016/j.cpc.2025.109961","url":null,"abstract":"<div><div>A hybrid spectral/finite-element code is developed to numerically solve the resistive finite-pressure magnetohydrodynamic equilibria without the necessity of postulating nested magnetic flux surfaces in the non-axisymmetric toroidal systems. The adopted approach integrates a hyperbolic parallel damping equation for pressure updating, along with a dynamic resistive relaxation for magnetic field. To address the non-axisymmetry in toroidal geometry, a pseudo flux mapping is employed to relate the axisymmetric computational domain to the physical domain. On the computational mesh, an isoparametric C1-continuous triangular element is utilized to discretize the poloidal plane, which is complemented with a Fourier decomposition in the toroidal direction. The versatility of the code is demonstrated through its application to several different non-axisymmetric toroidal systems, including the inherently three-dimensional equilibria in stellarators, the helical-core equilibrium states in tokamak plasmas, and the quasi-single-helicity states in a reversed-field pinch.</div></div>","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 109961"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145571120","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-03-01Epub Date: 2025-12-26DOI: 10.1016/j.cpc.2025.110012
Dandan Liao , Lei Ye , Pengfei Zhao , Qilong Ren , Nong Xiang
A second-order semi-implicit time integration method has been developed for solving the linearized collision operator in multi-species plasmas. A key feature of this method is that it treats the collision operator as a single entity, avoiding the operator splitting between its test particle and field particle components. This scheme employs an implicit trapezoidal time integration scheme for the isothermal test particle part (including pitch-angle scattering and energy diffusion) with a finite volume discretization in (v∥, μ) velocity coordinates, while the non-isothermal model term and field particle part are treated explicitly. This approach avoids the species cross-coupling required by fully implicit schemes, enabling the collision term to be computed efficiently with a banded/sparse-matrix solver. The method has been implemented in the gyrokinetic semi-Lagrangian code NLT and verified through multi-species relaxation tests and neoclassical transport simulations. Numerical benchmarks against explicit methods confirm its robustness, achieving order-of-magnitude improvements in the allowable time-step size, particularly in simulations of electron-ion collisions. Furthermore, the numerical discretization rigorously preserves particle number, momentum, and energy conservation, maintains the self-adjoint property of the collision operator, and satisfies Boltzmann’s H-theorem.
{"title":"Semi-implicit scheme for multi-species collision operators in Tokamak Plasma Simulations","authors":"Dandan Liao , Lei Ye , Pengfei Zhao , Qilong Ren , Nong Xiang","doi":"10.1016/j.cpc.2025.110012","DOIUrl":"10.1016/j.cpc.2025.110012","url":null,"abstract":"<div><div>A second-order semi-implicit time integration method has been developed for solving the linearized collision operator in multi-species plasmas. A key feature of this method is that it treats the collision operator as a single entity, avoiding the operator splitting between its test particle and field particle components. This scheme employs an implicit trapezoidal time integration scheme for the isothermal test particle part (including pitch-angle scattering and energy diffusion) with a finite volume discretization in (<em>v</em><sub>∥</sub>, <em>μ</em>) velocity coordinates, while the non-isothermal model term and field particle part are treated explicitly. This approach avoids the species cross-coupling required by fully implicit schemes, enabling the collision term to be computed efficiently with a banded/sparse-matrix solver. The method has been implemented in the gyrokinetic semi-Lagrangian code NLT and verified through multi-species relaxation tests and neoclassical transport simulations. Numerical benchmarks against explicit methods confirm its robustness, achieving order-of-magnitude improvements in the allowable time-step size, particularly in simulations of electron-ion collisions. Furthermore, the numerical discretization rigorously preserves particle number, momentum, and energy conservation, maintains the self-adjoint property of the collision operator, and satisfies Boltzmann’s H-theorem.</div></div>","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 110012"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145880486","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-03-01Epub Date: 2025-11-28DOI: 10.1016/j.cpc.2025.109966
Seokhwan Min, Jonghwa Shin
Controlling the scattering of waves from multi-shell spherical systems and particles is a crucial aspect in many applications in photonics such as superdirective antennae and structural coloring. Nevertheless, the effective design of such systems is non-trivial due to the coexistence of topological (number of shells and their material composition) and shape (shell thicknesses) parameters. Thus far, general-purpose algorithms such as parameter sweeps, gradient descent, differential evolution, and deep neural networks have been used to optimize particle shape under one or a few fixed topologies, limiting the complexity and effectiveness of the resulting designs. To address this shortcoming, we present a topology nucleation algorithm that allows the concurrent design of particle topology and shape through the use of a topology derivative expression derived from the transfer matrix formulation of the analytical Mie scattering theory. The principle behind our algorithm can readily be applied to the design of multi-shell spherical systems in other fields such as acoustics and quantum transport.
{"title":"Eschallot: A topology nucleation algorithm for designing stratified, spherically symmetric systems that exhibit complex angular scattering of electromagnetic waves","authors":"Seokhwan Min, Jonghwa Shin","doi":"10.1016/j.cpc.2025.109966","DOIUrl":"10.1016/j.cpc.2025.109966","url":null,"abstract":"<div><div>Controlling the scattering of waves from multi-shell spherical systems and particles is a crucial aspect in many applications in photonics such as superdirective antennae and structural coloring. Nevertheless, the effective design of such systems is non-trivial due to the coexistence of topological (number of shells and their material composition) and shape (shell thicknesses) parameters. Thus far, general-purpose algorithms such as parameter sweeps, gradient descent, differential evolution, and deep neural networks have been used to optimize particle shape under one or a few fixed topologies, limiting the complexity and effectiveness of the resulting designs. To address this shortcoming, we present a topology nucleation algorithm that allows the concurrent design of particle topology and shape through the use of a topology derivative expression derived from the transfer matrix formulation of the analytical Mie scattering theory. The principle behind our algorithm can readily be applied to the design of multi-shell spherical systems in other fields such as acoustics and quantum transport.</div></div>","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 109966"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145733048","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-03-01Epub Date: 2025-12-07DOI: 10.1016/j.cpc.2025.109988
Jaemyung Kim , Yujiro Hayashi , Sung Soo Ha , Makina Yabashi
In X-ray diffraction-based orientation microscopy, reconstructed grain structures can exhibit unrealistic or erroneous features due to the broadening and overlapping of diffraction peaks. Accurate grain boundary determination based on physical models remains a critical challenge for reliable microstructural characterization. While Voronoi tessellation is widely used to represent microstructures, its accuracy is often limited by the lack of weighting factors, leading to biased results. To address this, we developed a grain extraction algorithm combining a variation of the label-equivalent connected components labeling method with the marching squares algorithm for precise grain boundary detection. Using the extracted grain shapes, additively weighted Voronoi tessellation (AWVT) was applied, with each grain’s center of mass (COM) and equivalent radius serving as weighting factors. The AWVT boundaries showed strong agreement with experimental data, outperforming conventional Voronoi and Laguerre tessellations. Furthermore, the relationship between AWVT and curvature-driven grain growth models is discussed, demonstrating the method’s potential for improved microstructure characterization and grain growth analysis.
{"title":"Tessellation-based grain boundary determination for X-ray orientation microscopies","authors":"Jaemyung Kim , Yujiro Hayashi , Sung Soo Ha , Makina Yabashi","doi":"10.1016/j.cpc.2025.109988","DOIUrl":"10.1016/j.cpc.2025.109988","url":null,"abstract":"<div><div>In X-ray diffraction-based orientation microscopy, reconstructed grain structures can exhibit unrealistic or erroneous features due to the broadening and overlapping of diffraction peaks. Accurate grain boundary determination based on physical models remains a critical challenge for reliable microstructural characterization. While Voronoi tessellation is widely used to represent microstructures, its accuracy is often limited by the lack of weighting factors, leading to biased results. To address this, we developed a grain extraction algorithm combining a variation of the label-equivalent connected components labeling method with the marching squares algorithm for precise grain boundary detection. Using the extracted grain shapes, additively weighted Voronoi tessellation (AWVT) was applied, with each grain’s center of mass (COM) and equivalent radius serving as weighting factors. The AWVT boundaries showed strong agreement with experimental data, outperforming conventional Voronoi and Laguerre tessellations. Furthermore, the relationship between AWVT and curvature-driven grain growth models is discussed, demonstrating the method’s potential for improved microstructure characterization and grain growth analysis.</div></div>","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 109988"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145732949","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
We introduce an out-of-equilibrium replica exchange Monte Carlo (REMC) scheme designed to accelerate the equilibration of hard body systems under high-pressure conditions. The method deliberately violates balance by significantly increasing the swap acceptance probability, thereby forcing replicas to explore the entire range of sampling pressures. As a result, replicas undergo abrupt pressure changes analogous to annealing processes that help reduce defects. We demonstrate the efficacy of the approach on systems of N hard disks confined within a circular cavity (with 400 ≤ N ≤ 450), where it consistently achieves higher maximal packing fractions than standard REMC within the same number of cycles. Furthermore, we establish new maximal packing records for cases where the best known maximal packings fall below Cantrell’s low-density limit and for N ≤ 1000. Similar improvements are observed for disks confined within a square cavity. Finally, we discuss the potential of alternating cycles of out-of-equilibrium and equilibrium REMC to further approach the equilibrium equation of state of hard body systems at very high densities.
我们介绍了一种非平衡复制交换蒙特卡罗(REMC)方案,旨在加速高压条件下硬体系统的平衡。该方法故意通过显著增加交换接受概率来破坏平衡,从而迫使副本探索采样压力的整个范围。因此,仿制品会经历类似于退火过程的突然压力变化,从而有助于减少缺陷。我们证明了该方法在圆形腔内(400 ≤ N ≤ 450)的N个硬盘系统上的有效性,在相同次数的循环内,它始终比标准REMC获得更高的最大填充分数。此外,我们建立了新的最大包装的情况下,最知名的最大包装低于Cantrell的低密度极限和N ≤ 1000。类似的改进被观察到圆盘限制在一个方形腔内。最后,我们讨论了非平衡和平衡REMC交替循环的潜力,以进一步接近非常高密度下硬体系统的状态平衡方程。
{"title":"Densest packings and accelerated equilibration of hard body systems via out-of-equilibrium replica exchange Monte Carlo! method","authors":"Eduardo Basurto , Peter Gurin , Szabolcs Varga , Gerardo Odriozola","doi":"10.1016/j.cpc.2025.109990","DOIUrl":"10.1016/j.cpc.2025.109990","url":null,"abstract":"<div><div>We introduce an out-of-equilibrium replica exchange Monte Carlo (REMC) scheme designed to accelerate the equilibration of hard body systems under high-pressure conditions. The method deliberately violates balance by significantly increasing the swap acceptance probability, thereby forcing replicas to explore the entire range of sampling pressures. As a result, replicas undergo abrupt pressure changes analogous to annealing processes that help reduce defects. We demonstrate the efficacy of the approach on systems of <em>N</em> hard disks confined within a circular cavity (with 400 ≤ <em>N</em> ≤ 450), where it consistently achieves higher maximal packing fractions than standard REMC within the same number of cycles. Furthermore, we establish new maximal packing records for cases where the best known maximal packings fall below Cantrell’s low-density limit and for <em>N</em> ≤ 1000. Similar improvements are observed for disks confined within a square cavity. Finally, we discuss the potential of alternating cycles of out-of-equilibrium and equilibrium REMC to further approach the equilibrium equation of state of hard body systems at very high densities.</div></div>","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 109990"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145836600","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-03-01Epub Date: 2025-11-27DOI: 10.1016/j.cpc.2025.109959
Zhixin Lu, Guo Meng, Roman Hatzky, Philipp Lauber, Matthias Hoelzl
The features of the TRIMEG-GKX code are described with emphasis on the exploration using novel/different schemes compared to other gyrokinetic codes, particularly the use of object-oriented programming, filter/buffer-free treatment, and a high-order piecewise field-aligned finite element method. The TRIMEG-GKX code solves the electromagnetic gyrokinetic equation using the particle-in-cell scheme, taking into account multi-species effects and shear Alfvén physics. The mixed-variable/pullback scheme has been implemented to enable electromagnetic studies. This code is parallelized using particle decomposition and domain cloning among computing nodes, replacing traditional domain decomposition techniques. The applications to study the micro- and macro-instabilities are demonstrated, including the energetic-particle-driven Alfvén eigenmode, ion temperature gradient mode, and kinetic ballooning mode. Good performance is achieved in both ad hoc and experimentally reconstructed equilibria, such as those of the ASDEX Upgrade (AUG), Tokamak á configuration variable (TCV), and the Joint European Torus (JET). Future studies of edge physics using the high-order C1 finite element method for triangular meshes in the TRIMEG-C1 code will be built upon the same numerical methods.
{"title":"TRIMEG-GKX: An electromagnetic gyrokinetic particle code with a piecewise field-aligned finite element method for micro- and macro-instability studies in tokamak core plasmas","authors":"Zhixin Lu, Guo Meng, Roman Hatzky, Philipp Lauber, Matthias Hoelzl","doi":"10.1016/j.cpc.2025.109959","DOIUrl":"10.1016/j.cpc.2025.109959","url":null,"abstract":"<div><div>The features of the TRIMEG-GKX code are described with emphasis on the exploration using novel/different schemes compared to other gyrokinetic codes, particularly the use of object-oriented programming, filter/buffer-free treatment, and a high-order piecewise field-aligned finite element method. The TRIMEG-GKX code solves the electromagnetic gyrokinetic equation using the particle-in-cell scheme, taking into account multi-species effects and shear Alfvén physics. The mixed-variable/pullback scheme has been implemented to enable electromagnetic studies. This code is parallelized using particle decomposition and domain cloning among computing nodes, replacing traditional domain decomposition techniques. The applications to study the micro- and macro-instabilities are demonstrated, including the energetic-particle-driven Alfvén eigenmode, ion temperature gradient mode, and kinetic ballooning mode. Good performance is achieved in both ad hoc and experimentally reconstructed equilibria, such as those of the ASDEX Upgrade (AUG), Tokamak á configuration variable (TCV), and the Joint European Torus (JET). Future studies of edge physics using the high-order <em>C</em><sup>1</sup> finite element method for triangular meshes in the TRIMEG-C1 code will be built upon the same numerical methods.</div></div>","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 109959"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145681673","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}
Pub Date : 2026-03-01Epub Date: 2025-11-19DOI: 10.1016/j.cpc.2025.109955
Laurent Hascoët , Matt Menickelly , Sri Hari Krishna Narayanan , Jared O’Neal , Nicolas Schunck , Stefan M. Wild
The HFBTHO code implements a nuclear energy density functional solver to model the structure of atomic nuclei. HFBTHO has previously been used to calibrate energy functionals and perform sensitivity analysis by using derivative-free methods. To enable derivative-based optimization and uncertainty quantification approaches, we must compute the derivatives of HFBTHO outputs with respect to the parameters of the energy functional, which are a subset of all input parameters of the code. We use the algorithmic/automatic differentiation (AD) tool Tapenade to differentiate HFBTHO. We compare the derivatives obtained using AD against finite-difference approximation and examine the performance of the derivative computation.
{"title":"HFBTHO-AD: Differentiation of a nuclear energy density functional code","authors":"Laurent Hascoët , Matt Menickelly , Sri Hari Krishna Narayanan , Jared O’Neal , Nicolas Schunck , Stefan M. Wild","doi":"10.1016/j.cpc.2025.109955","DOIUrl":"10.1016/j.cpc.2025.109955","url":null,"abstract":"<div><div>The HFBTHO code implements a nuclear energy density functional solver to model the structure of atomic nuclei. HFBTHO has previously been used to calibrate energy functionals and perform sensitivity analysis by using derivative-free methods. To enable derivative-based optimization and uncertainty quantification approaches, we must compute the derivatives of HFBTHO outputs with respect to the parameters of the energy functional, which are a subset of all input parameters of the code. We use the algorithmic/automatic differentiation (AD) tool Tapenade to differentiate HFBTHO. We compare the derivatives obtained using AD against finite-difference approximation and examine the performance of the derivative computation.</div></div>","PeriodicalId":285,"journal":{"name":"Computer Physics Communications","volume":"320 ","pages":"Article 109955"},"PeriodicalIF":3.4,"publicationDate":"2026-03-01","publicationTypes":"Journal Article","fieldsOfStudy":null,"isOpenAccess":false,"openAccessPdf":"","citationCount":null,"resultStr":null,"platform":"Semanticscholar","paperid":"145681674","PeriodicalName":null,"FirstCategoryId":null,"ListUrlMain":null,"RegionNum":2,"RegionCategory":"物理与天体物理","ArticlePicture":[],"TitleCN":null,"AbstractTextCN":null,"PMCID":"","EPubDate":null,"PubModel":null,"JCR":null,"JCRName":null,"Score":null,"Total":0}