Papers

We present a real-space formulation of DFT for extended systems in which localized Wannier-like functions are constructed directly from localized Gaussian basis functions without explicitly computing canonical Bloch-like states during the self-consistent cycle. Building on the formalism of Pederson and Lin [Phys. Rev. B 35, 2273 (1987)], a variational set of Wannier-like functions is generated self-consistently within a finite Wannier domain and used to construct the charge density, electrostatic potential, and per cell total energy of the extended system. The occupied space can be determined entirely from the Wannier-like functions. Electronic band structures can be recovered in a post-processing step by solving the full Hamiltonian in a Bloch-like basis constructed from Gaussian orbitals. A key feature of the method is that it can incorporate combined translation–rotation, or screw, symmetries, enabling efficient simulations of chiral and helical systems with finite twist angles at essentially the same computational cost as systems described by pure translational symmetry. The approach is validated through calculations on linear and twisted $\mathrm{{-}C{\equiv}C{-}}$ and $\mathrm{{-}Li{-}F{-}}$ chains, as well as graphene, where total energies and band structures show excellent agreement with reference periodic calculations. To illustrate the ability of the method to treat three-dimensional systems, it is further applied to AA graphite, in which carbon atoms in adjacent graphene layers are aligned directly above one another, as well as helically stacked AA graphite structures. The Wannier-first framework provides a practical route for treating extended systems with nontrivial translational, rotational, and screw symmetries, and provides a natural foundation for the implementation of orbital-dependent functionals such as the Perdew–Zunger self-interaction correction.

In many systems, the interaction between a deformable surface or interface and the surrounding bulk fluid is coupled with intrinsic spatiotemporal dynamics within the moving surface. Examples include tumor growth, biological tissue morphogenesis, cardiac mechanics, multi-phase surfactant chemistry, additive manufacturing, clothing wear-and-tear, and reactive combustion flows. Solving such problems requires both geometric computing algorithms to track and resolve the surface and numerical methods to solve the coupled governing equations in the surface and the surrounding bulk phase. Here, we present a fully meshfree numerical solver for such coupled bulk-surface problems with deformable interfaces. The presented solver tracks the surface implicitly, solving for the dynamic surface geometry based on stress balance coupled to surrounding fluid phases. We show convergence for a mass-conserving case on a growing sphere and solve bulk-surface problems with incompressible Navier-Stokes fluids coupled to in-surface nonlinear reaction-diffusion dynamics. Finally, we show a model of biological morphogenesis, solving simultaneously for the dynamic surface shape and the fields on the curved surface with two-way coupling.

Immersed-boundary lattice Boltzmann (IB-LBM) simulations of cylinder and particle flows are usually compared by ranking boundary-enforcement schemes, kernels, and collision models. Such comparisons can obscure the mechanism that controls each response. Here, canonical fixed-cylinder, oscillating-cylinder, and sedimenting-particle cases are used as controlled probes of Eulerian-Lagrangian interface momentum transfer. For prescribed bodies, differences among direct forcing (DF), multi-direct forcing (MDF), and distribution-function correction (DFC) are most clearly discriminated by local no-slip fidelity rather than by a universal drag ranking. For DFC, the kernel-dependent drag-coefficient ranking reversal between the hat and Peskin 4-point kernels is associated with the spatial redistribution of the marker-resolved correction and the resulting near-boundary slip and pressure deviation, rather than the total correction magnitude alone. For sedimenting particles, single-particle settling provides a moving-body baseline, whereas two-particle differential-density wake-interaction sedimentation shows the wake-exposed light particle to be comparatively more sensitive to the explicit internal-mass correction in the force evaluation, a configuration-dependent, finite-window moving-boundary closure response. Targeted two-relaxation-time (TRT) and central-moment multiple-relaxation-time (CM-MRT) collision controls remain secondary to the boundary and closure mechanisms within the tested regimes (prescribed-body comparisons at Re <= 200, two-dimensional, and the reported moving-particle cases). The resulting picture replaces a universal scheme ranking with a mechanism-resolved interpretation of IB-LBM interface momentum transfer.

The Fokker-Planck partial differential equation (FP-PDE) governs uncertainty evolution in stochastic dynamical systems. In orbital dynamics, solving the FP-PDE is challenging because of nonlinear motion, high-dimensional states, and large space-time domains. We develop a physics-informed neural network (PINN) approach that approximates the FP-PDE solution as a single space-time probability density, while also quantifying its worst-case approximation error. This approach is, in principle, independent of the choice of state coordinates and neural network architecture. Specifically, to enforce probability density function (PDF) properties into the neural network, we design a Physics-informed Gaussian mixture model (PINN-GMM). Then a companion error PINN learns the dynamics of the approximation error and yields time-dependent bounds that define an ambiguity set of PDFs. This ambiguity set enables rigorous computation of upper and lower bounds on event probabilities through tractable linear programs. Numerical studies on illustrative 1D examples and several 4D–6D orbital test cases demonstrate accurate uncertainty propagation, correct and informative error bounds, and improved reliability over common uncertainty-propagation baseline methods (Gaussian approximation, unscented transform, and Gaussian mixture model). Constructing the PINN-GMM requires offline training, making it costlier than the baseline approximations; once trained, however, a single forward pass returns the density at any time in sub-millisecond time $(0.16~\mathrm{ms}$ in our implementation).

We develop a structure-preserving, data-driven collision operator for spatially homogeneous multi-species kinetic systems from molecular dynamics (MD). The operator consists of diagonal self-collision blocks and ordered off-diagonal cross-species blocks to describe intra- and inter-species momentum and energy exchange. Within a local and point-wise identifiable kernel class, we develop the necessary and sufficient condition for the admissible kernel class satisfying the conservation laws, the H-theorem, and the frame indifference. Unlike the classical Landau operator, the off-diagonal kernels are not restricted to be symmetric under permutation of the two velocity variables. This unique structural freedom captures the distinct responses of different species to unresolved correlations and many-body effects arising from micro-scale particle interactions. The equivalent parameterizable kernel formalization enables us to learn a generalized data-driven collision operator directly from MD, where the low-rank tensor representations and random sampling are used to achieve efficient kernel training and numerical simulation. Numerical experiments show that the learned operator accurately predicts transport coefficients and the non-equilibrium relaxation, while retaining discrete conservation and entropy production. In particular, it captures plasma kinetics in the moderately coupled regime, where the predictions of both the Landau and the data-driven model restricted to velocity-permutation symmetry show significant discrepancies.

We investigate the refractive index $n$ of LiF using DFT+G$_0$W$_0$+BSE, with ionic-temperature effects included through QMD. We calculate photon-energy dispersions $n(\omega)$ and $k(\omega)$ at ambient pressure and $n(\rho,T)$ at 532 and 1550 nm under shock compression, where $\rho$ and $T$ vary together. Quasiparticle band structures at ambient and compressed conditions and the ambient-pressure orbital-projected density of states connect the optical response to the electronic structure. At ambient pressure, G$_0$W$_0$ yields a quasiparticle gap of 14.25 eV, close to the experimental 14.2 eV, while BSE reproduces the main excitonic feature at 12.5 eV, close to the observed 12.6 eV. Quasiparticle and excitonic effects are thus essential for accurate optical dispersion; the hybrid HSE functional does not reproduce the optical spectra. Along the principal Hugoniot up to 140 GPa, $n(\rho)$ agrees closely with shock data at 1550 nm and acceptably at 532 nm, and agrees better overall than earlier first-principles calculations. Above approximately 110-120 GPa, QMD produces a downturn relative to cold $n(\rho)$ curves, suggesting that ionic dynamics may contribute to deviations from the linear $n$-$\rho$ Gladstone-Dale relation. QMD-sampling and finite-k-grid uncertainties are estimated, and transition-peak broadening is assessed using a recently proposed criterion. The gap increases under pressure, with a transition to a $\Gamma\to\mathrm L$ indirect gap near 50 GPa. An exploratory calculation at approximately 1400 GPa finds a gap of about 24 eV, disfavoring gap-closure metallization. Structural and elastic benchmarks of the underlying models are also provided.

Nonadiabatic coupled electron-nuclear dynamics upon electronic excitation underpin the microscopic mechanism and rational modulation of diverse photoinduced functional phenomena in materials, yet their direct first-principles simulations remain computationally demanding. Here we develop a framework for nonadiabatic excited-state machine-learning molecular dynamics (EMLMD) simulations, where the nonequilibrium electronic information upon photoexcitation such as electron temperature is rigorously calibrated from high-precision real-time time-dependent density functional theory (rt-TDDFT) benchmark simulations, enabling accurate reconstruction of excited-state potential energy surfaces (PES). This framework natively incorporates the excited-state electron-phonon couplings and intrinsically captures photoinduced phonon anharmonicity, both of which are missing in standard machine learning molecular dynamics, thus delivering first-principles-level accuracy for excited-state atomic evolutions. Large-scale EMLMD simulations resolve time- and momentum-resolved phonon dynamics in photoexcited materials, directly uncovering the competition between photogenerated coherent phonons and thermal phonons during photoinduced phase transition of bismuth. It also simultaneously resolves elusive atomic-scale microscopic dynamics and global structural rearrangement for selenium photoamorphization. Balancing high accuracy and efficiency, EMLMD offers a versatile paradigm to tackle key challenges in the study of complex excited-state molecular dynamics.

Gausslets are among the few electronic-structure bases that permit the four-index electron–electron interaction to be replaced by an accurate two-index integral diagonal approximation (IDA). Near a many-electron nucleus, however, the nuclear attraction and core-electron Hartree field are individually large and substantially cancel. Treating the first as a full finite-basis matrix while treating the second with IDA leaves an avoidable imbalance. We introduce reference-density Hartree screening: the Hartree field of a chosen reference density is represented accurately, and IDA is applied only to density fluctuations about it. Tests on He, Ne, atomic F, F$_2$, and Cr$_2$ show large reductions in direct Hartree errors, including transfer of fitted neutral-atom fields to molecules. For Cr$_2$, atomic-core screening prevents the spurious HF collapse found with the unscreened $q=5$ and $q=7$ Hamiltonians, whereas finite-reference matching without screening does not. Screening leaves exchange and residual correlation unchanged. We therefore also introduce a low-rank one-particle correction that uses an accurate conventional Gaussian-basis Hartree–Fock calculation to match either occupied-space exchange information or the complete occupied Fock vectors. In F$_2$ and Cr$_2$, $X_{\rm HF}$ reproduces the finite-reference energy and occupied Fock vectors to numerical precision and the selected states return after orbital perturbations. For Cr$_2$, the corrected $q=5$ basis uses one quarter as many functions as the $q=7$ control while retaining sub-mHa mean-field accuracy. Screening provides the physical improvement to the direct field; the state-specific correction then restores the remaining accuracy of the Gaussian-basis mean-field reference.

The chaotic dynamics of a large one-dimensional lattice of nonlinear maps with a spatial coupling of increasing extent is explored using the covariant Lyapunov vectors (CLVs). A theoretical approach is developed that predicts the spectral variation of the CLVs and the spectrum of Lyapunov exponents. The theoretical approach is independent of the details of the mapping function and requires only a description of the lattice coupling. The general structure of the spatial power spectrum of the CLVs is determined using only the eigenvectors of the coupling matrix. This yields a decomposition of the CLVs into two regimes. The first regime, corresponding to the largest Lyapunov exponents, contains highly entangled CLVs with significant contributions from large scale spatial structures and aligns with the sorted eigenvectors, which occur in pairs with wavenumbers that increase linearly with their index. The second regime contains significant mixing of length scales with increasing index. The index separating these two regimes is approximately equal to the fractal dimension of the dynamics.

We investigate the Self-Healing Diffusion Monte Carlo (SHDMC) method using a one-dimensional model with periodic boundary conditions. An inversion symmetry is introduced to mimic the antisymmetry property of fermionic wave functions, with the bosonic and fermionic sectors being modeled by the even and odd eigenstates, respectively. As in realistic fermionic systems, the nodal structure is only partially constrained by symmetry, making this model a non-trivial testbed for nodal optimization algorithms such as SHDMC. We show that the nodal evolution under SHDMC iterations can be cast into a dynamical system exhibiting both attractive and repulsive fixed points. When applying SHDMC to this model in the absence of statistical noise, the fixed-node energy is found to increase during iterations and the node converges to an incorrect value. This occurs for any finite basis, even when the exact node is representable in it; the attractive fixed point approaches the exact node only as the basis becomes complete, while the convergence toward it becomes increasingly slow. We further show that this problem can be largely cured by modifying the nodal update criterion to give more importance to the nodal region. Our results suggest that, in the general case, an efficient nodal optimization may require both an improved nodal update criterion and a localized basis set allowing local distortions of the trial wave function, as is the case for this model.

Clayey slip surfaces control the reactivation and long-term deformation of slow-moving landslides and may experience thermal fluctuations from climate, seasonal ground-temperature changes, or subsurface heat sources. Experiments show that residual shear strength depends on temperature and shearing rate, yet most numerical approaches use temperature-independent strength parameters. We present a non-isothermal viscoplastic constitutive model for clayey slip surfaces implemented with zero-thickness interface elements. It includes temperature-dependent normal and tangential stiffness, progressive degradation of cohesion and friction angle, and rate-dependent viscoplastic slip governed by a non-associated flow rule. Coupling with hydraulic and thermal balance equations allows the interface response to evolve with stress state, temperature, aperture, and accumulated irreversible displacement. Validation against temperature-controlled drained ring-shear tests on bentonite and smectite-rich soils covers heating–cooling, cooling–heating, and combined thermal paths. The simulations reproduce thermal strengthening at slow shearing rates and thermal weakening or limited sensitivity at higher rates. Application to the Congress Street cut benchmark shows that zero-thickness elements improve the representation of strain localization and progressive failure. Increasing temperature progressively degrades interface strength, increases displacement, joint aperture, and shear strain, and accelerates sliding. Temperature-dependent interface degradation can therefore reduce the apparent stability margin of clayey slopes and should be included in slope-stability assessments involving thermal fluctuations.

Crystalline materials discovery is essential for energy, electronics, and catalysis, but the vast chemical and structural space makes exhaustive screening infeasible. Existing voxel-based methods are limited by the local receptive fields of three-dimensional convolutional neural networks and the posterior collapse of high-dimensional variational autoencoders. Here, we develop a Fourier Neural Operator (FNO)-based crystal-field solver that maps a prescribed chemical formula and lattice parameters to periodic number-density and electron-density fields. By operating on global Fourier modes, the solver captures long-range correlations in periodic crystal fields beyond conventional local convolutions. Building on this solver, we construct a coupled generation-solving framework in which a conditional variational autoencoder generates diverse candidate lattice parameters in a low-dimensional basis-coefficient space, followed by density-field prediction and atomic reconstruction through peak detection, position optimization, and weight optimization. The reconstructed structures are further screened using voxel-level filtering, machine-learning interatomic-potential relaxation, and first-principle calculations. The framework generates novel structures across 104 chemical formulas with competitive reconstruction accuracy, demonstrating high generative diversity and structural validity. By extending Fourier neural operators to periodic crystal fields and coupling them with composition-conditioned lattice generation, our approach provides a scalable route to crystal structure discovery from prescribed chemical compositions.

Tensor networks, originally developed for quantum many-body physics, are promising models for machine learning. We derive stochastic Riemannian optimizers for tree tensor networks (TTNs) on both their parameter and quotient manifolds, including adaptive and learning-rate-free schemes suitable for minibatch training. Using a hybrid CNN-TTN architecture, we evaluate the methods on Fashion-MNIST, CIFAR10, and Imagenette. The proposed optimizers achieve predictive performance comparable to unconstrained optimization while enabling numerically stable downstream compression.

Optical waveguides are key elements for high fidelity, long distance optical communications. Coupled waveguide arrays allow for higher information density, steerting the propagation direction, and for encoding information. However, due to the mixing of relative phases for short pulses containing multiple waveguide-mode frequencies, a process for retrieving an encoded input state once these signals undergo coherent propagation remains elusive. A concept is presented to extract with high fidelity the phase-encrypted input signal from spatio-temporal propagated states. As a realization, an array of coupled waveguides is suggested with the retrieval mechanism being realized by local phase shifts that comply with the identified retrieval concept. Three dimensional full-wave electromagnetic simulations for broadband optical signals in coupled dielectric waveguides confirm the validity of the scheme and the high fidelity of information retrieval pointing to potential applications, for instance in ultrafast coherent coding and decoding of information imprinted on pulse sequences.

We assess computed lift coefficient data for flow past a circular cylinder to evaluate their suitability for practical applications. Specifically, we consider lift coefficient data for a circular cylinder over Reynolds numbers from 120 to 8000. The results are obtained from two-dimensional finite element simulations of the incompressible Navier-Stokes equations using pressure robust discretizations. We compare the computed lift coefficients with published experimental and numerical results, finding good agreement in some cases but significant disagreement in others. Because lift fluctuations are central to vortex-induced vibration concepts, these data therefore provide input for the analysis and preliminary design of bladeless turbines.

Machine-learning interatomic potentials (MLIPs) can reproduce potential energies and forces accurately, but their internal energy allocation is often difficult to interpret. E3D-IQA is introduced as a diagnostic framework connecting the latent edge-energy representation of an Allegro-type MLIP with Interacting Quantum Atoms (IQA) energy decomposition. The Allegro edge-energy path is retained as a latent pair contribution, while a node-energy path is trained against the IQA intra-atomic energy. IQA interatomic energies are not direct training targets; instead, the learned edge energies are evaluated after training against the IQA pair terms. Tests on H/C/N/O organic reaction structures show that intra-atomic supervision is essential: energy and force training alone does not recover an IQA-like one-body/two-body allocation. With intra-atomic supervision, node energies reproduce IQA intra-atomic terms, and latent edge energies show meaningful correspondence with IQA interatomic terms. Residual errors are concentrated in positive or weak pair interactions, exposing internal allocation failures that remain hidden in total-energy and force metrics. Adding structures labeled only with energies and forces improves transfer to larger molecules and reduces decomposition errors. E3D-IQA therefore provides a route for diagnosing and guiding latent scientific representations using partially labeled quantum-chemical datasets.

The electron self-energy is central to quasiparticle theory, yet how an optical cavity enters it remains unclear. We address this question for a molecule in a single-mode cavity using the dipole-gauge Pauli-Fierz Hamiltonian and a coherent-state QED Hartree-Fock reference. The cavity enters through three channels: the static dipole self-energy (DSE) shift of reference orbital energies, direct DSE augmentation of the screened interaction, and the polariton pole carrying the bilinear electron-photon coupling. We benchmark QED-$GW$ ionization potentials (IPs) and electron affinities (EAs) against a cavity $\Delta$-method ladder from QED-HF to correlated wave-function methods, whose cavity-induced shifts agree within 1 meV where directly comparable. For closed-shell molecules with unbound anions, $GW$ systematically overestimates cavity-induced IP redshifts, whereas EA shifts are reproduced nearly quantitatively, although this does not imply comparable accuracy for absolute EAs. For ionic molecules with bound anions, this ordering reverses, consistent with published QED coupled-cluster results. Coupling and detuning scans show that the error is predominantly quadratic in $\lambda$ and DSE-driven rather than resonant. The spectral function develops a polariton-replica photoemission sideband with weight scaling as $\lambda^2$. In the static screened interaction used in the Bethe-Salpeter equation, bare-photon exchange cancels the matching DSE contribution to the direct interaction, while exchange and polariton-screening corrections remain. Their net effect on the lowest excitation is appreciable only for ammonia in the molecules studied. Exciton-binding energies involving unbound anions are strongly basis-dependent and should therefore be viewed as diagnostics of electron-hole interactions rather than basis-converged molecular quantities.

We present a spectral element lattice Boltzmann method (LBM) for partial wetting on curved geometries. A non-ideal gas phase-field model is incorporated into the LBM framework to enable phase separation with a constant interface thickness and the potential form of surface tension force is used. We adopt the force-splitting approach, yielding significantly improved stability and accuracy. Complex boundaries are naturally handled using a flux bounce-back scheme, which resolves inconsistencies in normal vectors across adjacent elements. Additionally, a general wetting boundary condition is implemented to capture static contact line in a thermodynamically consistent manner. The method is validated through simulations of droplets on flat surfaces, 2/3-dimensional curved surfaces, and equilibrium droplets without boundaries. Results demonstrate that parasitic currents are significantly reduced on unstructured meshes with complex geometries, reaching residual kinetic energy levels on the order of $10^{-24}$ for wetting configurations and $10^{-30}$ for isolated droplets.

We prove that the Hypergraph Neural Network, an invariant architecture with 3-body message passing, is a universal approximator for potential energy surfaces. Our main contribution is a multi-layer completeness theory. We show that $L$ layers of message passing on sparse, cutoff-based graphs achieve the same representational power as having access to the full $L$-hop neighborhood, provided the configurations are generic, satisfy an overlap condition and a connectivity condition. This provides the first rigorous justification for the common practice of using multi-layer message passing with a per-layer cutoff smaller than the physical interaction range, the setting used by virtually all practical graph neural network based machine-learned interatomic potentials. As immediate consequences, we show that both DPA3 and CHGNet architectures inherit universal approximation.

Hyperelastic constitutive models enable modeling large deformations in elastic solids. In common practice, a strain energy density function is prescribed in advance and model-specific parameters are calibrated from experiments. However, many applications require constitutive models for a family of related materials whose mechanical behavior varies with composition. A fixed constitutive model-form may not capture the full range of behavior across the family, while fitting separate forms does not provide a direct way to predict the response of new compositions. Recent work has developed data-driven constitutive models that learn flexible strain energy functions while incorporating key physical constraints. In this work, we propose using mixtures of convex neural potentials based on input convex neural networks as a modular and data efficient approach to modeling material families. Each potential is convex and monotonic with respect to polyconvex strain invariants, while a conditioning network maps material descriptors to mixture weights. We compare the approach with a monolithic partially input-convex neural network using experimental data from PolyJet 3D-printed materials and a synthetic Gent-type benchmark. Across both benchmarks, the mixture architecture generalized better to material descriptors not seen during training. In the PolyJet experimental benchmark, we showed that the mixture architecture is less sensitive to model hyperparameters, while in the Gent-type benchmark it generalized more reliably with sparse data in the material-descriptor space. These results suggest that representing a material family through a small set of shared convex neural potentials provides a useful structural prior for learning descriptor-dependent constitutive behavior from limited data.

We present a meshfree numerical solver for the incompressible Navier-Stokes equations on oriented curved surfaces that are represented by surface point clouds. On curved surfaces, numerical challenges pertaining to stiffness and pressure-velocity coupling are exacerbated. Moreover, vector calculus on curved surfaces differs from its Euclidean counterpart. The presented method operates on surface point clouds in an Eulerian frame of reference without requiring a computational grid or mesh. It achieves consistent approximation in space and time with high order of accuracy; we demonstrate up to order six. The incompressibility constraint is locally imposed as a weak artificial compressibility approximation, avoiding global matrix inversion. We show that the method provides consistent and convergent approximations of surface vector fields and differential operators. We study the relationship between error, spatial resolution, and artificial Mach number and characterize the frequency spectrum of the artificial oscillations. We provide numerical solutions of the incompressible Navier-Stokes equations on symmetric surfaces, such as the sphere and torus, and on parametric and non-parametric asymmetric surfaces. Since the proposed method works directly on unstructured surface point clouds, it provides a promising approach for simulations on image-derived geometries, such as in biological morphogenesis from microscopy videos.

We describe the staged modernization of a local multilayer TRIM.SP implementation, denoted TRIM.SP-NL. The work included behavior-preserving refactoring, correction of localized bookkeeping and initialization defects, support for more than five elements per layer, reproducible build modes, and a runtime-selectable random-number-generator (RNG) interface. The RANLUX RNG remains the default reference backend, while PCG32 and xoshiro256** provide faster alternatives without changing the input-deck format. The code validation combined exact regression testing with statistical comparisons over 13 representative target multilayer configurations, comprising 143 energy-configuration points per backend and 429 simulations in total. Across the test suite, PCG32 and xoshiro256** reduced runtime by factors of 1.98 and 2.01, respectively. The implanted, backscattered, and transmitted fractions, as well as the mean implantation depths, range straggling, and representative depth profiles showed no systematic RNG dependence. Together with a Node.js-based graphical interface for setup, scans, and result collection, the Fortran engine forms the open-source TRIM.SP-NL Workbench.

We study how colloidal inclusions modify the formation energy of dislocation pairs in lamellar block copolymer systems. Using a hybrid particle/Ginzburg–Landau model, we calculate defect formation energies by comparing defect-free and defect-containing states with and without embedded colloids. Finite-size scaling is used to obtain formation energies in the thermodynamic limit. The effect of colloid insertion depends strongly on particle sizes and surface patterning. Homogeneous particles increasingly stabilize dislocation pairs with increasing particle sizes. Particles larger than one lamellar domain preferentially occupy the dislocation cores, where they replace strained polymer rather than deforming defect-free lamellae. The magnitude of this stabilization depends on surface affinity. Balanced Janus particles instead increase the formation energy, because their competing surface preferences cannot be satisfied simultaneously near the curved core. Varying the patch ratio interpolates between these limits. Surface patterning has little effect for particles smaller than one lamellar domain but changes the formation energy by several tens of $k_BT$ for larger particles. These results provide quantitative guidelines for controlling topological defect stability in lamellar block copolymer systems.

The Poisson–Nernst–Planck (PNP) equations strongly couple ionic transport and electrostatic interactions through the Poisson equation, posing substantial numerical challenges under small permittivity and complex potential boundary conditions. Underlying these equations is a natural Wasserstein gradient-flow structure, in which the Poisson equation serves as a local realization of the nonlocal electrostatic interaction energy. Exploiting this structure, we formulate each time step as a constrained convex minimization problem where the ionic continuity equations and the Poisson equation are incorporated as linear constraints, allowing the concentrations, fluxes, and electrostatic potential to be updated simultaneously. The variational structure of the scheme intrinsically guarantees the dissipation of the original free energy, mass conservation, and nonnegativity of ionic concentrations under general electrostatic boundary conditions. Moreover, the framework is structurally modular: extending from classical to modified PNP models with steric interactions and concentration-gradient corrections requires only modifying the energy functional, while all structure-preserving properties are automatically retained. To efficiently solve the resulting large-scale constrained problems, we develop preconditioned and transformed primal–dual algorithms equipped with tailored fast dual solvers, namely DCT-based direct and Schur-complement iterative methods, that exploit the coupled block structure of the PDE constraints. Numerical experiments on classical and modified PNP systems demonstrate the accuracy and structure-preserving properties of the scheme, and show that the proposed algorithms converge reliably in strongly coupled small-permittivity regimes without significant growth in computational cost.

Mesoscopic physical models, including the Hydrodynamic Drude Model (HDM), the Generalized Nonlocal Optical Response (GNOR) Model, and the Surface Response Model (SRM), have been proposed to investigate nonlocal effects in nanometric structures. The combination of classical electromagnetism with these mesoscopic material models calls for new computational electromagnetic (CEM) algorithms, or update of conventional ones, in what is termed computational mesoscopic electromagnetics (CMEM). In this work, we present a MATLAB toolbox for the simulation of multiple spherical interfaces with arbitrary relative positions, with the incorporation of the aforementioned mesoscopic models. The method exploits vector spherical wave functions to properly express electric and magnetic fields, an S matrix formulation for the efficient treatment of incident and scattered fields at spherical interfaces, and a translation matrix to deal with propagating waves with different expansion centers. Excitation sources can be chosen among arbitrarily polarized plane waves, dipoles and electron beams. For the post-processing part, the calculation of cross sections and far/near-field mapping; fluorescence enhancement, Purcell factor and quantum yield; and cathodoluminescence and electron energy-loss probability is implemented. The toolbox is built in a modular manner, and each part (routine) has its own important functionality. This paper provides a full explanation of the proposed highly efficient and general toolbox, and a detailed guideline for researchers in the nanoplasmonics community.

Dark excitons play a central role in the nonequilibrium dynamics of two-dimensional semiconductors, but remain difficult to characterize. Transient-absorption experiments, with probes tuned in exciton-exciton transitions energy range (exc-tr-abs), can detect excitations from any populated dark excitons, including symmetry-forbidden, spin-forbidden, and finite-momentum ones. In this work, we develop a $GW$+BSE scheme for computing exc-tr-abs spectra from arbitrary populated exciton distributions. These dark excitons are included equally by evaluating exciton-exciton dipoles in a locally smooth gauge, including intra- and inter-band contributions. For monolayer MoS$_2$, the exc-tr-abs signal arises from the $\Gamma$, $K$, $M$, and $Q$ valleys, differing substantially from the $\Gamma$-only interpretation. Exciton-exciton dipoles show similar intensities across these valleys, while their weights are dictated by initial excitonic populations. State- and spin-resolved analyses assign the peaks to $1s \rightarrow 2p$ and $1s \rightarrow 3p$ transitions from both spin-flip and spin-conserving A and B excitons across valleys and momenta.

Magnetic materials exhibit an intricate coupling between atomic structure and spin degrees of freedom, posing a fundamental challenge for atomistic simulations across experimentally relevant length and time scales. Here we introduce HotPP-Spin, a spin-dependent extension of HotPP for magnetic machine learning interatomic potentials, built on Cartesian tensor equivariant message passing. Atomic magnetic moments are treated as explicit axial-vector degrees of freedom, while spatial-inversion and time-reversal parities are propagated through the tensor couplings. This construction provides a unified representation of exchange-dominated and spin-orbit-induced interactions without imposing predefined analytical interaction forms. A scalar spin-dependent potential energy surface yields energy-conserving atomic forces and magnetic effective fields through differentiation. Benchmarks spanning collinear magnetism, noncollinear magnetism, and spin-orbit-coupling-induced magnetic anisotropy show that HotPP-Spin accurately describes magnetic energy landscapes, magnetic forces, and magnetic-order-dependent energy-volume relations within the same general framework. For H-phase monolayer VSe\(_2\), stochastic spin-dynamics simulations using the learned magnetic effective fields locate the finite-size magnetic ordering crossover at 415–435~K, in close numerical agreement with the reported experimental value of \(418.5\pm7.8\)~K. These results establish Cartesian tensor message passing as a general route for connecting first-principles magnetic energetics with large-scale atomistic simulations of coupled structural and spin phenomena.

The chemical potentials of charged multi-component mixtures are central to electrolyte thermodynamics, but remain difficult to compute from atomistic simulations. The S0 method enables computing chemical potentials of mixtures from equilibrium molecular dynamics simulations. Here we extend the S0 method to charged mixtures by combining the composition-space framework developed in Part I: Neutral Multi-component Mixtures with a Coulombic treatment of the small-wavenumber limits of static structure factors. This approach separates thermodynamically relevant neutral composition fluctuations from forbidden macroscopic charge fluctuations, and further accounts for charge-neutrality constraints. We use the method to compute the chemical potentials of multiple-halide aqueous salt solutions, elucidating the ion-specific thermodynamic effects. We also calculate the mixing free energies of molten salt mixtures, demonstrating the importance of correctly describing long-wavelength electrostatic correlations.

This study presents a theoretical investigation of generalized solute dispersion in magnetohydrodynamic multiphase tube flow with porous layers. A two-fluid analytical model is developed for applications in biofluid and environmental fluid dynamics. The model comprises a micropolar (non-Newtonian) fluid core representing the rotational behaviour of red blood cells and a Newtonian plasma periphery embedded with Brinkman and Darcy porous structures, corresponding to the glycocalyx and endothelial layers with distinct permeability characteristics. A transverse magnetic field is incorporated to investigate how magnetic-field-induced modifications of the carrier flow influence solute localisation, with potential relevance to magnetic nanoparticle-mediated drug delivery. Using the generalised dispersion framework of Sankarasubramanian & Gill, analytical solutions are derived to investigate how the coupled axial velocity field and associated microrotational dynamics influence solute transport. The analytical predictions are independently validated through Brownian dynamics simulations, demonstrating excellent agreement for the temporal evolution of the zeroth and first transport moments. The results reveal the previously unexplored influence of microrotational dynamics on solute concentration, convection coefficients and effective dispersion, providing new insights into the coupled roles of translational and rotational fluid motion in biofluid transport. This work bridges an important gap in the literature and establishes a generalized theoretical framework linking magnetic fields, micropolar fluids and porous arterial structures for biofluid transport, targeted drug delivery and clinical engineering applications.

Future wireless systems are expected to transform the surrounding space from a passive propagation medium into a smart electromagnetic environment, where engineered surfaces control wave propagation, support wireless sensing, and create programmable electromagnetic fingerprints. A key challenge in realizing this vision is the inverse design of metasurfaces for tailored electromagnetic propagation. While forward analysis evaluates the response of a known geometry, the inverse task starts from a prescribed scattering signature and seeks a physically realizable structure that produces it. This inverse task is inherently nonlinear and often high-dimensional, while candidate solutions may be non-unique and provide no direct indication of practical realizability. Here, we introduce a conditional diffusion framework for inverse design of dielectric resonator metasurfaces from target angular scattering patterns. Trained on T-matrix simulated geometry-response pairs, the model learns a conditional distribution of geometries instead of a deterministic mapping, enabling multiple candidate designs for the ill-posed inverse problem. The best generated metasurface achieves a mean percentage error of 1.39%, outperforming CMA-ES optimization (4.1% after 10 h) while requiring only about one minute for after-training inference. The model also produces lower error distributions than deterministic neural baselines for out-of-distribution spectra, highlighting the potential of diffusion models for efficient metasurface design.