Discovery of Unstable Singularities (2509.14185v1)

Abstract: Whether singularities can form in fluids remains a foundational unanswered question in mathematics. This phenomenon occurs when solutions to governing equations, such as the 3D Euler equations, develop infinite gradients from smooth initial conditions. Historically, numerical approaches have primarily identified stable singularities. However, these are not expected to exist for key open problems, such as the boundary-free Euler and Navier-Stokes cases, where unstable singularities are hypothesized to play a crucial role. Here, we present the first systematic discovery of new families of unstable singularities. A stable singularity is a robust outcome, forming even if the initial state is slightly perturbed. In contrast, unstable singularities are exceptionally elusive; they require initial conditions tuned with infinite precision, being in a state of instability whereby infinitesimal perturbations immediately divert the solution from its blow-up trajectory. In particular, we present multiple new, unstable self-similar solutions for the incompressible porous media equation and the 3D Euler equation with boundary, revealing a simple empirical asymptotic formula relating the blow-up rate to the order of instability. Our approach combines curated machine learning architectures and training schemes with a high-precision Gauss-Newton optimizer, achieving accuracies that significantly surpass previous work across all discovered solutions. For specific solutions, we reach near double-float machine precision, attaining a level of accuracy constrained only by the round-off errors of the GPU hardware. This level of precision meets the requirements for rigorous mathematical validation via computer-assisted proofs. This work provides a new playbook for exploring the complex landscape of nonlinear partial differential equations (PDEs) and tackling long-standing challenges in mathematical physics.

• The paper presents a novel computational framework to discover unstable self-similar singularities in fluid PDEs using high-precision PINNs and full-matrix Gauss-Newton optimization.
• It demonstrates an empirical linear relationship between the inverse scaling rate and instability order in models like IPM, Boussinesq, and CCF.
• The study emphasizes rigorous error control and the integration of mathematical structure in neural representations to achieve validation at machine precision.

Discovery of Unstable Singularities in Nonlinear Fluid PDEs

Introduction and Context

The formation of singularities in solutions to nonlinear PDEs governing fluid dynamics, such as the 3D Euler and Navier-Stokes equations, remains a central open problem in mathematical physics. Singularities correspond to finite-time blow-up of gradients or velocities from smooth initial data, with implications for the predictive validity of these models. While stable singularities—those robust to perturbations—have been numerically identified in several settings, the existence and structure of unstable singularities, which require infinitely precise initial conditions and are destroyed by infinitesimal perturbations, are far less understood. The paper "Discovery of Unstable Singularities" (2509.14185) presents a systematic computational framework for discovering and validating unstable self-similar singularities in canonical fluid models, including the incompressible porous media (IPM) and Boussinesq equations, and the 1D Córdoba-Córdoba-Fontelos (CCF) model.

Self-Similar Singularities and Instability Hierarchies

The approach is based on reformulating the governing PDEs in self-similar coordinates, parameterized by a scaling exponent $\lambda$ that determines the spatiotemporal structure of the blow-up. The task reduces to identifying admissible $\lambda$ values for which smooth self-similar profiles exist. The paper demonstrates the existence of a hierarchy of singularities, indexed by instability order $n$, with each higher $n$ corresponding to an additional unstable direction in the linearized spectrum.

Figure 2: Self-similar singularities for IPM and Boussinesq equations, showing spatial vorticity profiles, cross-sections, and the empirical linear relationship between inverse scaling rate and instability order.

The results reveal a linear empirical relationship between the inverse scaling rate $1/\lambda$ and the instability order $n$ for both IPM and Boussinesq systems. This provides a predictive formula for the scaling exponents of higher-order unstable singularities, which is critical for guiding future numerical searches and for understanding the asymptotic structure of singularity hierarchies.

High-Precision PINN Framework and Solution Validation

The computational framework employs physics-informed neural networks (PINNs) with carefully designed architectures that encode mathematical symmetries, boundary conditions, and asymptotic behaviors. The optimization is performed using a full-matrix Gauss-Newton method, which, in contrast to standard first-order or approximate second-order optimizers, enables convergence to solutions with residuals at or near double-precision machine epsilon.

Figure 1: Comparison of optimizer performance and convergence for the CCF first unstable solution, demonstrating the superiority of the Gauss-Newton method and the impact of multi-stage training.

Multi-stage training is used to further reduce high-frequency residuals by sequentially training correction networks, yielding solutions with maximum residuals as low as $O(10^{-13})$ for the CCF equation and $O(10^{-8})$ to $O(10^{-11})$ for IPM and Boussinesq. This level of accuracy is sufficient for subsequent computer-assisted proofs (CAPs) using interval arithmetic, which require rigorous control of numerical errors.

Figure 4: Spatial distribution of equation residuals for the first unstable IPM solution after multistage training, and summary of log-10 maximum residuals for all discovered solutions.

Mathematical Structure and Inductive Bias in Neural Representations

A key methodological advance is the explicit incorporation of mathematical structure into the neural network ansatz. This includes:

• Input coordinate transformations to handle infinite domains and enforce symmetries.
• Output envelopes to impose decay at infinity and regularity at the origin.
• Analytical relationships for $\lambda$ derived from smoothness constraints, reducing the search space and improving convergence.
• Iterative feedback between numerical experiments and mathematical analysis, allowing for the discovery and encoding of emergent structural properties (e.g., vanishing rates at the origin).

This approach ensures that the PINN does not merely interpolate data but is constrained to the physically and mathematically relevant solution manifold, which is essential for resolving highly unstable and non-generic singularities.

Stability Analysis and Completeness of Discovered Families

For each discovered singularity, the linearized spectrum of the PDE around the self-similar profile is computed. The number of unstable modes matches the instability order $n$, confirming the expected structure of the solution hierarchy. The empirical completeness of the family is supported by the absence of additional admissible $\lambda$ values within the explored range.

Figure 3: Maximum residuals as a function of $\lambda - \lambda_i$ for the Boussinesq equation, showing sharp minima at admissible $\lambda_i$ and rapid growth away from these values, indicating isolated admissible scaling exponents.

Implications and Future Directions

The discovery of unstable singularities in the IPM and Boussinesq equations, with validation at machine precision, provides a new paradigm for the paper of singularity formation in nonlinear PDEs. The empirical scaling laws for $\lambda_n$ offer a roadmap for constructing singularity hierarchies in more complex systems, including the 3D Euler and Navier-Stokes equations. The demonstrated precision and validation pipeline set a new standard for the use of machine learning in mathematical physics, emphasizing the necessity of domain-specific inductive bias and rigorous error control.

The approach is not intended as a general-purpose PDE solver but as a targeted discovery tool for elusive, non-generic solutions. The integration of PINNs with high-precision optimization and mathematical structure is likely to be broadly applicable to other problems in nonlinear analysis, including the search for unstable or non-unique solutions in other classes of PDEs.

The main open challenge remains the extension of these methods to boundary-free 3D Euler and Navier-Stokes equations, where the existence and structure of unstable singularities are directly tied to the Millennium Prize Problem. Further advances will require deeper understanding of the qualitative properties of singular solutions and continued development of mathematically informed neural architectures and training protocols.

Conclusion

This work establishes a computational and mathematical framework for the systematic discovery and validation of unstable self-similar singularities in nonlinear fluid PDEs. By combining PINNs with full-matrix Gauss-Newton optimization, multi-stage training, and mathematically structured neural representations, the authors achieve unprecedented accuracy and uncover new families of unstable singularities, along with empirical scaling laws for their blow-up rates. These results have significant implications for the analysis of singularity formation, the design of computer-assisted proofs, and the future integration of machine learning with rigorous mathematical analysis in PDE theory.