2026-11-10 | | Total: 12
Isogeometric beam formulations provide an efficient reduced-dimensional representation of slender structures, combining accurate geometry and high-order continuity with few structural degrees of freedom. Their application to fluid--structure interaction (FSI) requires transferring motion and loads between the one-dimensional beam centerline and the surrounding fluid interface. While mapping strategies exist for conventional finite element beams, geometrically nonlinear isogeometric formulations require operators consistent with their spline-based kinematics and distinct descriptions of cross-sectional rotation. This work develops consistent beam-to-surface mappings for geometrically nonlinear isogeometric Bernoulli and Timoshenko beams in partitioned FSI. Fluid-interface points are attached to the beam in the reference configuration and reconstructed using the nonlinear kinematics of each formulation. Consistent linearization yields a configuration-dependent mapping operator whose transpose ensures work-conjugate load transfer. The formulations are verified through analytical deformation tests, convergence studies, and virtual-work checks, and evaluated in beam-based adaptations of the Mok and Turek--Hron benchmarks and the Goland wing flutter benchmark. The results confirm accurate interface reconstruction and work-conjugate load transfer under large deformations, while demonstrating the applicability of both beam formulations to coupled FSI problems.
Numerical simulations of combustion processes rely on chemical kinetic models, whose predictive capability is typically assessed by comparing measurements of quantities of interest (QoIs) in idealized experimental configurations with the corresponding simulated values. Sensitivity analysis quantifies the influence of model parameters on QoIs, making it a key tool for model development and optimization. Given the large number of parameters involved in detailed kinetic models compared to the number of QoIs, adjoint methods are a promising approach for efficient sensitivity analysis. Ignition delay times (IDTs) and extinction strain rates (ESRs) are among the QoIs of primary importance. However, applying adjoint methods to these quantities is not straightforward and is addressed in this work. The system of governing equations is reformulated via a coordinate transformation that makes the QoI an explicit component of the solution vector. This reformulation enables the application of adjoint sensitivity analysis for both IDT and ESR. A computational complexity analysis demonstrates that the proposed formulation significantly improves algorithmic efficiency compared to conventional approaches. The methodology is verified by comparing adjoint sensitivity coefficients with those obtained from other methods. Execution time measurements show computational cost reductions of up to two orders of magnitude relative to standard direct sensitivity analysis. A benchmark study is conducted across kinetic models of increasing size to assess how far theoretical efficiency gains are realized in practice. While the results show that implementation details can limit performance for smaller problem sizes, the expected asymptotic scaling is observed for larger models. The benchmark thus provides a reference framework for identifying key performance bottlenecks and guiding subsequent optimization efforts.
We integrate GEMMul8, a library that emulates FP64 matrix multiplication on INT8 Tensor Cores, into the GPU version of the numerical-atomic-orbital first-principles code OpenMX 4.0. The GEMMul8 configuration is compared against the cuBLAS FP64 configuration within each of two environments, a consumer GeForce RTX 5080 and an HPC-oriented H100 PCIe, for the four computational paths formed by collinear (CL) or noncollinear (NCL) spin and band or $Γ$-point-cluster solvers. In fixed-work measurements that completed 25 SCF iterations, the diagonalization section of a 384-atom CL band calculation on the RTX 5080 is accelerated by a factor of 2.1 and the total run time by 1.8; a 512-atom CL cluster calculation, which uses only real matrix products in contrast to the CL band calculation, still gains a factor of 1.5 in the diagonalization section. On the H100, GEMMul8 provides little additional acceleration of the diagonalization section. In an accuracy test on a 216-atom Si system in which one atom was displaced by 0.05 angstrom, all 24 configurations were converged using an SCF convergence criterion of $10^{-8}$; every configuration converged at the same iteration within each path, and the differences between GEMMul8 and cuBLAS are at most $2\times10^{-12}$ Ha in the total energy and $7\times10^{-12}$ Ha/Bohr in force components.
Unsteady compressible simulations with high-order methods on unstructured meshes require boundaries that let acoustic, entropy and vorticity waves leave the domain without spurious reflection. For the spectral difference (SD) method, the Navier-Stokes characteristic boundary conditions (NSCBC) have only been formulated on quadrilaterals and hexahedra. Their strong-form formalism relies on a line of flux points normal to the boundary face, which does not exist on triangles and tetrahedra. This paper extends the NSCBC to the SD method on simplex elements, including multi-species reacting gases. The characteristic analysis is performed in the physical frame of each flux point of the boundary face, and the normal derivative of the normal flux is evaluated using the Raviart-Thomas divergence operator. The boundary fluxes are recovered from a linear system local to each boundary element, whose matrix depends on the geometry only. The treatment is compared with the tensor-product formulation at a comparable number of degrees of freedom on a series of test cases covering acoustic, vortical, viscous and reacting flows on triangles and tetrahedra. In all cases, the results are in close agreement with those obtained on quadrilaterals and hexahedra, so that the NSCBC become available for the SD method on simplex meshes, which are well suited to complex geometries.
At finite temperature and pressure, the equilibrium structure of a crystal minimizes the Gibbs free energy rather than the static Born-Oppenheimer energy. We present WARP (Wyckoff-constrained Anharmonic Relaxations at Pressure and temperature), a symmetry-constrained, variable-cell relaxation based on the stochastic temperature-dependent effective potential method, including anharmonicity and zero-point motion. Every position update is projected onto Wyckoff positions and every strain onto symmetric strains, so a relaxation can raise the symmetry but never lower it. For aluminum, the thermal expansion follows measurements at cryogenic temperatures, including the quantum plateau below 50 K. For lanthanum hydrides at megabar pressures, equations of state agree with phonon-based and path-integral calculations, and LaH$_{10}$ relaxes from its rhombohedral static minimum to the observed cubic phase, also with a machine-learned potential fine-tuned on 150 supercells. For brucite, Mg(OD)$_2$, lattice parameters, bond distances and thermal expansion agree with neutron diffraction, and the O-D isotope effect is captured.
The method of projection onto localized states intended for describing reactions and clustering is examined. In particular, the method offers a simple extension of the no-core shell model to continuum- and clusterized bound states. A critical analysis of a previous realization of the method is provided. A way to minimize the number of basis functions required for convergence is given. An easy way to obtain exactly unitary and symmetric $S$ matrix from calculations is pointed out. The method is tested in one- and two-channel $np$ scattering problems with the Daejeon $NN$ interaction. It is shown that excellent convergence can be reached even when only oscillator states with low number of excitation quanta are retained.
Asynchronous parallel simulations of stochastic processes with absorbing states can exhibit systematic sampling bias when trajectory durations depend on their eventual absorbing outcomes. Workers that continuously restart upon absorption preferentially record rapidly terminating trajectories within a fixed computational budget, distorting the observed ensemble statistics. We investigate this effect using a minimal renewal model with two absorbing states and state-dependent exponential durations. Renewal theory yields exact expressions for the expected number and composition of completed trajectories, showing that the resulting bias decreases algebraically with execution time but persists as the number of parallel workers increases at fixed budget. We compare continuous sampling with First-Trajectory Filtering, which retains only the initial trajectory attempted by each worker. Although finite execution deadlines introduce residual censoring bias, this bias decreases exponentially with the computational budget. An analytical mean-squared error calculation reveals budget-dependent crossover worker allocations beyond which First-Trajectory Filtering becomes more accurate despite collecting fewer trajectories. We also examine the limitations of statistical inference under severe censoring when the absorbing-state probabilities and characteristic durations are unknown. Finally, simulations of the Axelrod model near its nonequilibrium transition demonstrate the preferential sampling of rapidly absorbing fragmented configurations, illustrating the relevance of the renewal mechanism in a spatially interacting system. Our results show that increasing computational throughput does not necessarily improve statistical accuracy and highlight the importance of accounting for trajectory-duration heterogeneity in parallel simulations of absorbing processes.
The defining features of a physical system are its symmetries and its invariants. For a network, the symmetries are the automorphisms, and the invariants are the vertex orbits they define. Two vertices of a graph may appear structurally indistinguishable without being interchangeable under any graph automorphism. Distinguishing apparent similarity from genuine symmetry lies at the heart of the vertex-orbit partition problem. In this work, we develop a hierarchy of increasingly expressive, automorphism-invariant signatures that draw on the spectral properties of the graph Laplacian and the Dirac operator to reveal subtle structural differences. Guided by a sequence of counterexamples, we progressively strengthen these signatures, ultimately incorporating probing techniques, which are temporary structural perturbations, that expose the hidden distinctions in the original unperturbed graph. We evaluate our approach across diverse benchmark families, recovering exact orbit partitions in challenging instances, including graphs that resist low-dimensional Weisfeiler--Leman refinement. We also examine the computational and numerical limitations of our constructions, distinguishing empirical success from theoretical completeness. Our findings highlight the potential of spectral and probing-based methods for uncovering the structural distinctions that define vertex orbit.
Causal discovery from observational data is central to scientific reasoning, yet established constraint-based and score-based methods often recover only a Markov equivalence class and provide no signal in bivariate or strongly autocorrelated regimes typical in some areas of science. Here, we recast causal discovery as a Bayesian symbolic regression problem. Given a set of candidate parents for each target variable, the Bayesian symbolic regression samples closed-form expressions from a posterior over symbolic models. The fraction of sampled expressions in which each candidate variable appears defines a posterior inclusion probability that we use as the causal score. This replaces pointwise model selection by an ensemble of plausible mechanisms. We evaluate the method on coupled nonlinear bivariate dynamics, on the CauseMe benchmark suite spanning linear, chaotic, nonlinear, higher-dimensional and high-noise settings, and on a four-variable biogeoscience structural causal model under Gaussian, heavy-tailed and mixed-noise regimes. Our approach is competitive with, and often outperforms, specialized causal discovery algorithms. Additionally, the method returns closed-form mechanisms leading to explicit causality and interpretability.
We present \textit{Hydrograd} v2.0, a significantly improved differentiable 2D shallow-water model written in Python with JAX. The model solves the shallow-water equations on the unstructured meshes of production models such as SRH-2D and HEC-RAS 2D with a well-balanced, positivity-preserving finite-volume scheme that handles wetting and drying. Gradients are tailored to the type of problem. Steady states are found by Newton's method with an exact sparse Jacobian obtained by graph colouring, and their gradients follow from the implicit function theorem at the cost of one linear solve. Unsteady flows are differentiated in reverse mode with checkpointing, or with an implicit time integrator whose gradient costs about two forward runs. All gradients agree with finite differences, to about $10^{-8}$ for steady flows. On the Savannah River case, the steady water surface differs from that of SRH-2D by 4.7\,mm (root mean square). Steady closure learning was 96 times faster than differentiating through time marching and reached a lower loss. An inflow hydrograph was inverted from stage records with a peak error of 2\,\%. In an unsteady roughness calibration, a singular value decomposition analysis of the sensitivity matrix showed that the stage records at the gauges are 70 to 160 times more sensitive to the total resistance of the reach than to its distribution among the roughness zones, which explains why only a gauge layout with gauges in every zone recovered all five roughness values. On one GPU, a mesh of 512,000 cells takes 7.6\,s for a forward run and 51\,s for a gradient.
We prove universal approximation for the bilinear-readout equivariant message-passing network (BREMP), a reference architecture for equivariant graph neural network interatomic potentials. Its layers combine equivariant convolutions with affine maps, and its readout forms invariants from products of the equivariant features. We construct a one-layer feature that determines the local environment up to permutations of atoms of the same type. The construction uses spherical harmonic filters of degree no larger than the largest number of neighbors at a common distance from the center. Iterating the construction recovers the environment reachable in as many hops as there are layers, using the same bound on filter degrees at the layer cutoff. The readout approximates every continuous invariant function of the features uniformly on compact sets. Together these give universal approximation of atomic energies at a prescribed physical range using a shorter layer cutoff. We assume that, throughout the admissible class, message passing reaches all neighbors in each physical environment within a common number of layers. Finally, DPA4, an interatomic potential architecture in production use, inherits this result by approximating the BREMP layers constructed in the proof and realizing the bilinear readout.
Single precision (FP32) halves the memory and data movement of double precision (FP64), and on accelerators is often several times faster. Its rounding errors are harmless except in the stored state. A modified equation for Runge-Kutta methods in Shu-Osher form and linear multistep methods separates the rounding of the stored state, a forcing proportional to $1/Δt$, from truncation and from right-hand-side rounding, which sets a $Δt$-independent floor. We bound this forcing for strong-stability-preserving schemes and identify two mechanisms that keep its sign fixed, so it accumulates as a bias: past a turnaround step, refining the step increases the error, in FP64 as in FP32. For a model oscillator the analysis predicts the floor and the FP32-to-FP64 ratios of turnaround step and minimum error within 30%. Our remedy needs no FP64 arithmetic: each state variable is stored as the unevaluated sum of two FP32 numbers, updated with error-free transformations, while every right-hand-side evaluation stays in FP32. Unlike compensated summation, it needs no rewriting of the scheme. In model problems from geophysical flow, mechanical dynamics and celestial mechanics, with five Runge-Kutta and multistep schemes, it restores design-order convergence to the floor. In large-eddy simulations of a high-lift airplane wing on three refined grids, up to 544 million cells, with time steps $Δt$, $Δt/2$ and $Δt/4$, FP32 overpredicts the lift by 2.1, 4.5 and 7.7%, while the two-float lift agrees with FP64 within 0.3% on every grid, and on graphics processing units the production solver runs 1.6 to 2.7 times faster than FP64.