Conceptio › Archive › arXiv CS

arXiv CSopen access

Balance-Guided Sparse Identification of Multiscale Nonlinear PDEs with Small-coefficient Terms

2026 · arxiv_cs

arXiv CS · Papers · License: Open Access · 2026

Open Source ↗Direct PDF ↓

neural-networks

machine learning, deep learning, neural networks

Graphical Abstract Balance-Guided Sparse Identiﬁcation of Multiscale Nonlinear PDEs with Small-coeﬃcient Terms

(a) Data sampling & derivative estimation

(b) Linear system construction

Spatiotemporal data

(d) Progressive sparsification

…

l Numerical Differentiation l Automatic Differentiation

Neural Network

(e) Selection of the optimal model

...

arXiv:2604.18414v1 [cs.LG] 20 Apr 2026

Zhenhua Dang, Lei Zhang, Long Wang, Guowei He

Discovered equation

Remove the least active term

Importance for term j

Highlights Balance-Guided Sparse Identiﬁcation of Multiscale Nonlinear PDEs with Small-coeﬃcient Terms Zhenhua Dang, Lei Zhang, Long Wang, Guowei He • Balance-guided sparsiﬁcation enables identiﬁcation of physically important small-coeﬃcient terms. • Progressive pruning discards insigniﬁcant terms based on physical importance instead of predeﬁned coeﬃcient thresholding. • Residual-based model selection identiﬁes the optimal sparse governing equation.

Balance-Guided Sparse Identiﬁcation of Multiscale Nonlinear PDEs with Small-coeﬃcient Terms Zhenhua Danga,b,c , Lei Zhanga,b , Long Wanga,b,∗, Guowei Hea,b,∗∗ a

State Key Laboratory of Nonlinear Mechanics, Institute of Mechanics, Chinese Academy of Sciences, Beijing 100190, China b School of Engineering Sciences, University of Chinese Academy of Sciences, Beijing 100049, China c School of Mathematics and Statistics, Northwestern Polytechnical University, Xi’an 710072, China

Abstract Data-driven discovery of governing equations has advanced signiﬁcantly in recent years; however, existing methods often struggle in multiscale systems where dynamically signiﬁcant terms may have small coeﬃcients. Therefore, we propose Balance-Guided SINDy (BG-SINDy) inspired by the principle of dominant balance, which reformulates ℓ0 -constrained sparse regression as a term-level ℓ2,0 -regularized problem and solves it using a progressive pruning strategy. Terms are ranked according to their relative contributions to the governing equation balance rather than their absolute coeﬃcient magnitudes. Based on this criterion, BG-SINDy alternates between least-squares regression and elimination of negligible terms, thereby preserving dynamically signiﬁcant terms even when their coeﬃcients are small. Numerical experiments on the Korteweg–de Vries equation with a small dispersion coeﬃcient, a modiﬁed Burgers equation with vanishing hyperviscosity, a modiﬁed Kuramoto–Sivashinsky equation with multiple small-coeﬃcient terms, and a two-dimensional reaction–diﬀusion system demonstrate the validity of BG-SINDy in discovering small-coeﬃcient terms. The proposed method thus provides an eﬃcient approach for discovering governing equations that contain small-coeﬃcient terms. ∗

Corresponding author Corresponding author Email addresses: [email protected] (Long Wang), [email protected] (Guowei He) ∗∗

Keywords: SINDy, data-driven equation discovery, dominant balance principle, progressive pruning, multiscale systems 1. Introduction Data-driven discovery of governing equations has emerged as a powerful paradigm for uncovering the underlying mechanisms of complex dynamical systems, particularly when ﬁrst-principles derivations are intractable [1, 2, 3, 4, 5, 6, 7]. However, robustly identifying governing equations that contain small-coeﬃcient terms in multiscale systems remains a critical challenge [8]. Despite their small magnitudes, such terms cannot be neglected, as they may become locally dominant in speciﬁc spatiotemporal regions (e.g., viscous terms in high-Reynolds-number boundary layers [9, 10] and secular terms in the three-body problem [11]); omitting them disrupts the dominant balance and can lead to fundamentally incorrect models. Nevertheless, most existing equation-discovery methods rely on sparsity-promoting techniques, typically based on practical approximations to ℓ0 regularization of coeﬃcients, which select terms according to their coeﬃcient magnitudes [5, 6, 12]. Consequently, small-coeﬃcient terms are often inadvertently discarded during sparsiﬁcation, even if their physical contributions are essential in certain regimes [8]. To overcome this sparsity bias, we propose Balance-Guided SINDy (BGSINDy), inspired by the classical dominant balance principle [13, 14, 15, 16], to discover governing equations for multiscale systems with small-coeﬃcient terms. Over the past decade, substantial progresses have been achieved in discovering governing equations from observational data through sparse regression and symbolic regression approaches. The sparse identiﬁcation of nonlinear dynamics (SINDy) framework [5] pioneered sparse regression over libraries of candidate functions and has been extended to partial diﬀerential equations via PDE-FIND [6]. To improve robustness against noise and limited data, weak-form formulations [17] and ensemble learning [18] have been introduced. The widely used Python package PySINDy [19] has greatly facilitated the adoption of these methods. Multi-step-accumulation (MSA)-based sparse identiﬁcation has also been proposed to improve robustness to noisy and incomplete measurements by minimizing the accumulated prediction error over multiple time steps [20]. In addition, the integration of physics-informed neural networks (PINNs) [21, 22] has produced hybrid PINN-SINDy [23] 2

approaches that enable accurate discovery of governing equations even from scarce or highly noisy measurements. However, PINNs often struggle to accurately identify small perturbation parameters in singularly perturbed systems due to sharp boundary layers and severe loss imbalance; consequently, specialized variants such as the scale-decomposed PINN (SD-PINN) [24] have been developed to address these speciﬁc challenges. A powerful complement to sparse regression is symbolic regression, which discovers free-form analytical expressions by evolving symbolic structures from elementary operators and variables. A major breakthrough came with the AI Feynman algorithm, which integrates neural-network guidance with brute-force symbolic search [25]. This laid the foundation for advances in evolutionary computation, including symbolic genetic algorithms [26] and efﬁcient modern implementations such as PySR [27], which excel in chaotic and nonlinear dynamical systems. Kolmogorov–Arnold networks (KANs) [28, 29] have further bridged neural networks and symbolic regression by learning spline-based univariate functions that can be directly translated into explicit symbolic expressions. More recently, large-language-model-assisted approaches have emerged that dynamically update candidate libraries or generate novel expressions [30]. At the current frontier, generative transformer models such as EqGPT learn co-occurrence patterns of PDE terms directly from mathematical handbooks, enabling autonomous discovery of free-form equations [31]. Collectively, these methods constitute the primary avenues for interpretable and physics-consistent equation discovery from observational data. Despite these advances, accurately retaining small-coeﬃcient terms remains a signiﬁcant challenge in data-driven equation discovery. The SPIDER method [32] addresses this through residual-driven greedy pruning, which retains candidate terms by monitoring whether their removal causes a sharp increase in the residual. This approach has successfully recovered the Navier–Stokes equations with extremely low viscosity from high-Reynoldsnumber turbulent channel ﬂow data [32], as well as the governing equations of Rayleigh–Bénard convection [33], active nematics [34], and turbulent magnetohydrodynamic ﬂows [35]. The SPRINT algorithm [36] accelerates the pruning process using bisection with analytic bounds on rank-1 SVD updates while preserving sensitivity to minute coeﬃcients even at large library sizes. Similarly, the bi-level optimization framework BILLIE [37] decouples term selection from coeﬃcient estimation via reinforcement learning and independent validation splits, thereby preserving subtle yet dynamically essential 3

contributions even under noise or data scarcity. Despite their notable success, these methods typically incur substantial computational overhead due to exhaustive greedy searches, large-library SVD computations, or bi-level reinforcement learning. To provide a computationally eﬃcient and physically principled alternative, we draw inspiration from the dominant balance principle [13, 14, 15, 16], a cornerstone of asymptotic analysis in complex multiscale systems. By systematically comparing the relative magnitudes of terms in the governing equations, this principle identiﬁes the dominant physical contributions under speciﬁc limiting regimes. A classic example is the derivation of the Prandtl boundary-layer equations from the full Navier–Stokes equations in the highReynolds-number limit [9]. Complementary techniques, such as matched asymptotic expansions [38], multiple-scale analysis [39], and the WKB approximation [13], enable the systematic reduction of complex models into their essential dynamical structures. More recently, this paradigm has been extended to data-driven reduced-order models [40] when the complete governing equations are known a priori. Using unsupervised learning techniques such as Gaussian mixture models and sparse principal component analysis, local dominant balances can be identiﬁed as distinct clusters in an “equation space”, thereby automatically determining which terms may be safely neglected in diﬀerent spatiotemporal regions without compromising model ﬁdelity. In this paper, we propose Balance-Guided SINDy (BG-SINDy), a novel sparse regression framework for discovering governing equations of multiscale dynamical systems. Inspired by the principle of dominant balance, BGSINDy reformulates the conventional coeﬃcient-level ℓ0 regularization as a term-level ℓ2,0 problem and solves it via a progressive pruning strategy. By prioritizing the relative physical contributions of candidate terms to the overall equation balance rather than their magnitudes, the method iteratively alternates between constrained least-squares estimation and importance-based pruning. This enables accurate discovery of parsimonious models that retain essential small-coeﬃcient terms in multiscale systems. The main contributions of this work are summarized as follows: (i) A balance-guided sparse regression framework that reformulates coeﬃcientlevel ℓ0 regularization as term-level ℓ2,0 regularization. (ii) A progressive pruning strategy that sequentially eliminates candidate 4

terms according to their relative physical contributions, thereby avoiding heuristic coeﬃcient thresholds. (iii) A residual-based criterion that automatically selects the optimal sparse representation of the governing equation. The remainder of this paper is organized as follows. Section 2 presents the mathematical formulation and algorithmic implementation of the BG-SINDy framework. Section 3 validates the eﬀectiveness of the proposed method on challenging cases, including high-order derivative terms with small coeﬃcients, multiple small-coeﬃcient terms, and scalability to high-dimensional systems. Concluding remarks and future directions are provided in Section 4. 2. Method In this section, we ﬁrst review the SINDy method and then present the proposed BG-SINDy method. For simplicity, the following introduction of the two methods is illustrated using one-dimensional examples. 2.1. SINDy In the SINDy framework, the evolution of the state variable u(x, t) is assumed to satisfy ∂u = φ(u) · ξ, (1) ∂t where φ(u) is a row vector of candidate nonlinear functions constructed from polynomials of u and its spatial derivatives, and ξ is the coeﬃcient vector. The candidate library consists of polynomial terms up to degree P and their products with spatial derivatives up to order Q: q p p∂ u p = 0, . . . , P ; q = 1, . . . , Q . {u |p = 0, . . . , P } ∪ u ∂xq Accordingly, the total number of candidate terms is M = (P + 1)(Q + 1). Evaluating these terms at sampled spatiotemporal points (xi , tk ) yields the overdetermined regression system Ut = Φξ,

(2)

where Ut ∈ RN ×1 collects the time derivatives at N sampling points, and Φ ∈ RN ×M is the library matrix. Each row of Φ contains the evaluated candidate 5

terms at a given point, with Φij denoting the value of the j-th term at the i-th sample, and ξ ∈ RM ×1 . The required derivatives can be approximated either by classical numerical methods (ﬁnite diﬀerences, spectral methods, etc.) [6, 41, 42] or by automatic diﬀerentiation through neural-network surrogates [43, 44, 45]. The sparse vector ξ is then recovered by solving the ℓ0 -regularized leastsquares problem [5] ξ̂ = arg min ∥Ut − Φξ∥22 + λ∥ξ∥0 , ξ

(3)

where λ > 0 controls the level of sparsity, and ∥ξ∥0 counts the non-zero elements in the coeﬃcient vector ξ. In practice, this NP-hard problem [46] is approximately solved using sparsity-promoting algorithms within the SINDy framework. The original SINDy algorithm [5] employs the sequential thresholded least-squares (STLSQ) method to discover ordinary diﬀerential equations (ODEs), which iteratively alternates between ordinary least-squares regression and hard thresholding of small coeﬃcients. For the discovery of partial diﬀerential equations (PDEs), the PDE-FIND algorithm [6], a popular extension of the SINDy framework, adopts the TrainSTRidge method (train sequential thresholded ridge regression), which combines ridge regression with sequential thresholding. The threshold is selected via a trainvalidation split to improve robustness to noise. 2.2. BG-SINDy The BG-SINDy method promotes sparsity at the level of candidate terms according to their contributions to the equation balance evaluated at sampling points, rather than at the individual coeﬃcient level, in the spirit of the principle of dominant balance. The optimal sparse governing equation is obtained by solving ξ̂ = arg min ∥Ut − Φξ∥22 + λ diag(ξ)ΦT 2,0 , ξ

(4)

where Ut ∈ RN ×1 , Φ ∈ RN ×M , and ξ ∈ RM ×1 follow the deﬁnitions in Sec. 2.1. Here, diag(ξ) is a diagonal matrix formed from ξ, and ∥ · ∥2,0 counts the number of nonzero rows of a matrix [47, 48]. Since each nonzero row of diag(ξ)ΦT corresponds uniquely to one active candidate term, Eq. (4) reformulates the coeﬃcient-level ℓ0 regularization as an equivalent term-level ℓ2,0 regularization. In particular, diag(ξ)ΦT 2,0 6

counts the number of nonzero columns of Φ, i.e., the number of active candidate terms. This formulation enforces sparsity directly at the termlevel, enabling entire candidate functions to be selected or discarded based on their collective spatiotemporal contribution. This shifts the sparsity mechanism from coeﬃcient magnitude to term-level contribution, fundamentally altering the model selection criterion. In contrast to standard SINDy, which relies on coeﬃcient thresholding, BG-SINDy selects terms according to their dynamical importance, consistent with the dominant balance principle [13, 14, 15, 16]. To solve the l2,0 -norm constrained regression problem in Eq. (4) eﬀectively, we propose a progressive pruning workﬂow that realizes Occam’s razor [7, 49] through term-level sparsiﬁcation. As illustrated in Fig. 1, the procedure proceeds through the following iterative steps: (i) Spatio-temporal data sampling and derivative estimation (Fig. 1(a)): We ﬁrst collect spatio-temporal data and then compute the required derivatives accurately using either numerical diﬀerentiation or automatic diﬀerentiation [43, 45], for example via PyTorch-based neural networks [50]. The resulting dataset consists of N spatio-temporal sample points, forming the target time-derivative vector Ut ∈ RN . (ii) Linear system construction (Fig. 1(b))): Based on the state variables, an initial library of candidate terms is constructed. Progressive pruning assumes that no linear combination of candidate terms constitutes a valid equation on its own. This linear independence prevents term coupling and enables the individual elimination of redundant terms. If the initial library violates this, a maximal linearly independent subset is extracted a priori using singular value decomposition or QR factorization with column pivoting [51] to construct the candidate library matrix Φ ∈ RN ×M . We initialize the active set as A(0) = {1, 2, . . . , M } and set the iteration counter k = 0. The initial least-squares estimate is computed as: 2

ξ (0) = arg min Φ:,A(0) ξ − Ut 2 , ξ

(5)

with the initial residual of the ﬁtting model deﬁned as: Res0 =

1 2 Φ:,A(0) ξ (0) − Ut 2 . N 7

(6)

(iii) Regression and term importance calculation (Fig. 1(c)): At iteration k, given the regression coeﬃcients ξ (k) on the current active set A(k) , we quantify the contribution of each term to the equation balance in the same spirit as the principle of dominant balance. For each sample point i = 1, . . . , N and each active term j ∈ A(k) , the local relative importance is deﬁned as (k)

wij =

Φij ξj

(k)

maxl∈A(k) |Φil ξl | + ϵ

(7)

where ϵ > 0 is a small constant for numerical stability. It follows that 0 ≤ wij ≤ 1. (iv) Progressive pruning (Fig. 1(d)): The local relative importance is averaged over all sample points to obtain the global importance: 1 X wij . N i=1 N

Wj =

(8)

Consequently, 0 ≤ Wj ≤ 1. The least important term is identiﬁed as j ∗ = arg min Wj j∈A(k)

and removed from the active set: A(k+1) ← A(k) \ {j ∗ }. The regression coeﬃcients are then updated on the reduced active set by solving the unconstrained least-squares problem 2

ξ (k+1) = arg min Φ:,A(k+1) ξ − Ut 2 . ξ

(9)

The residual of the ﬁtted model is given by Resk+1 =

1 2 Φ:,A(k+1) ξ (k+1) − Ut 2 . N

(10)

(v) Selection of the Optimal Model (Fig. 1(e)): As terms are progressively removed, the residual typically remains approximately constant until a 8

(a) Data sampling & derivative estimation

(b) Linear system construction

Spatiotemporal data

… ...

(d) Progressive sparsification

l Numerical Differentiation l Automatic Differentiation

Neural Network

(e) Selection of the optimal model

Discovered equation

Remove the least active term

Importance for term j

Figure 1: Schematic of the BG-SINDy method. (a) Spatio-temporal data sampling and derivative estimation via numerical diﬀerentiation or automatic diﬀerentiation. (b) Construction of the linear system Ut = Φξ. (c) Regression and term importance calculation. (d) Progressive pruning by iteratively discarding the least signiﬁcant terms. (e) Identiﬁcation of the optimal sparse model based on the evolution of residuals.

dynamically essential term is eliminated, at which point it exhibits a substantial increase. After computing Resk+1 , the termination condition is evaluated: Resk+1 > τ, (11) Resk where τ = 3 in this work. When this condition is satisﬁed, the current model is deemed to have reached an overly sparse representation, and the optimal model is taken from the previous iteration corresponding to the step immediately before the increase in the residual: ξ ∗ = ξ (k) ,

A∗ = A(k) .

Otherwise, the iteration counter is updated as k ← k + 1, and the algorithm returns to step (iii). By transforming the combinatorial ℓ2,0 -regularized regression problem into a deterministic iterative procedure guided by physically interpretable criteria of dynamical contribution, the proposed BG-SINDy method oﬀers several key advantages. It enables robust identiﬁcation of governing equations from observational data while successfully preserving dynamically essential terms with small coeﬃcients. Furthermore, unlike traditional sparsitypromoting techniques that require delicate manual tuning of absolute coeﬃcient thresholds, BG-SINDy systematically identiﬁes the optimal model 9

structure through a residual-change criterion. This fundamentally enhances the physical interpretability of the identiﬁed governing equations and facilitates more faithful recovery of underlying dynamics. 3. Results We demonstrate the performance of the BG-SINDy method on four benchmark problems for discovering nonlinear governing equations containing smallcoeﬃcient terms. The method is ﬁrst evaluated on the Korteweg–de Vries (KdV) equation with a small dispersion coeﬃcient to assess its robustness to measurement noise and varying numbers of sampling points. Next, a modiﬁed Burgers equation augmented with vanishing hyperviscosity is then considered to assess its ability to identify high-order derivative terms with small coeﬃcients. Then, the method’s capability to handle equations with multiple small-coeﬃcient terms is examined using a modiﬁed Kuramoto– Sivashinsky (KS) equation. Finally, the method’s scalability is investigated using a two-dimensional reaction–diﬀusion equation. 3.1. KdV equation with a small dispersion coeﬃcient In this example, the KdV equation is used to evaluate the performance of the proposed BG-SINDy method, focusing on its ability to identify terms associated with small coeﬃcients. The KdV equation [52] describes the propagation of nonlinear dispersive waves in one spatial dimension and has been widely used in equation discovery studies [6, 21, 53, 54]. Here we consider the KdV equation with a small dispersion coeﬃcient [55] ∂u ∂u ∂ 3u +u + ε 3 = 0, ∂t ∂x ∂x

(x, t) ∈ Ω × [0, T ],

(12)

where ε = 4.84 × 10−4 , Ω = [0, 2], and T = 3. The initial condition is u(x, 0) = 0.9 sech2 12.45(x − 0.5) + 0.3 sech2 7.1875(x − 0.85) , x ∈ Ω, and homogeneous Dirichlet boundary conditions are imposed. The reference dataset is generated by numerically solving the KdV equation (12) using a fourth-order Runge–Kutta scheme in time combined with a ﬁnite-diﬀerence discretization in space [56]. The spatial domain is discretized into NΩ = 260 uniform grid points, and the time step is ∆t = 0.001.

To discover the governing equation from sparse and noisy data, the required derivatives are estimated via automatic diﬀerentiation [43, 44, 45] on a neural-network surrogate (a fully connected network with ﬁve hidden layers of 50 neurons each and sin activation) trained on only 500 points (0.064% of the full dataset) obtained via Latin hypercube sampling [57]. The candidate library consists of polynomial terms in u up to degree 2, together with their products with spatial derivatives of u up to order 4, yielding a total of 15 candidate terms: q p p∂ u {u | p = 0, 1, 2} ∪ u p = 0, 1, 2; q = 1, . . . , 4 . ∂xq The performance of BG-SINDy on the KdV equation is summarized in Fig. 2. Following the workﬂow illustrated in Fig. 1, the candidate library is progressively pruned by iteratively eliminating the least signiﬁcant term according to the importance score deﬁned in Eqs. (7) and (8). As shown in Fig. 2(a), this procedure yields a sequence of increasingly sparse models accompanied by a monotonic increase in the residual. The selected model, indicated by the red circle, corresponds to the elbow point (pronounced change) in the residual. The evolution of term importance during the sparsiﬁcation process is visualized in the heatmap in Fig. 2(b). The validity of the discovered equation is further veriﬁed by numerically integrating the model and comparing the resulting solution with the reference solution shown in Fig. 2(c), yielding a relative L2 error of 1.19%. A comparison with SINDy and PDE-FIND under the same candidate library and derivative estimation is shown in Fig. 2(d). When the equation form is correctly identiﬁed, accuracy is quantiﬁed by the coeﬃcient error [31]: coeﬃcient error =

N term X

|ξi − ξitrue | , Nterm i=1 |ξitrue | 1

(13)

where Nterm is the number of terms that appear in both the discovered and reference equations, and ξi , ξitrue are the corresponding coeﬃcients in the discovered and reference models, respectively. BG-SINDy successfully recovers the KdV equation, with a coeﬃcient error of only 1.35%, whereas both SINDy and PDE-FIND fail to identify the equation structure. The robustness of BG-SINDy to noise and sampling data size is assessed by adding Gaussian white noise to the clean data following [6, 31]: ũ = u + γ · std(u) · N (0, 1), 11

(a)

Normalized residual

(b)

Residual history in sparsification

1 u u2 ux uux u 2 ux uxx uuxx u2uxx uxxx uuxxx u2uxxx uxxxx uuxxxx 2 u uxxxx

Residual history Discovered equation

10−1

ut = −0.9955uux − 4.7308 × 10−4uxxx 10−2

10−3

0.4 0.2

0.4

1.0

0.2

0.5 0.5

1.0 x

1.5

2.0

0.0

0.4

1.0

0.2

0.5 0.0 0.0

0.9

1.5

0.3

0.5

0.0 0.0 1.0 1.5 2.0 0.0 x −4 Reference equation: ut = −uux − 4.84 × 10 uxxx 0.5

0.6

1.0

0.5

1.0 x

1.5

Method

Discovered equation

Coefficient error

BG-SINDy

ut = −0.9955uux − 4.7308 × 10−4uxxx

1.35%

ut = −0.3003ux

—

SINDy PDE-FIND

ut = −0.5279uux

×10−2 1.2

2.0

0.6

1.5

0.0

2.5

0.8

2.0

0.6

1.5

0.0 0.0

2.5

0.8

2.0

(d)

0.6

151413121110 9 8 7 6 5 4 3 2 1 0 15 14 13 12 11 10 9 8 7 6 5 4 3 2 1 Number of candidate terms Number of candidate terms Solution of reference equation Solution of discovered equation Absolute error 3.0 1.0 3.0 1.0 3.0 2.5

1.0 0.8

(c)

Term importance in sparsification

2.0

0.0

—

Figure 2: Results for the KdV equation. (a) History of residual in the sparsiﬁcation process. The red circle denotes the discovered equation. (b) Heatmap of term importance in the sparsiﬁcation process. (c) Comparison of the solutions obtained from the reference equation and the discovered equation. The relative L2 error is 1.19%. (d) Discovered equation by BG-SINDy and comparison with SINDy and PDE-FIND. BG-SINDy accurately recovers the reference KdV equation, while SINDy and PDE-FIND fail.

where γ denotes the noise amplitude. The coeﬃcient error, computed according to Eq. (13) for noise levels ranging from 0% to 25% and sampling sample sizes from 1,000 to 100,000, is presented in the heatmap of Fig. 3. Even at 25% noise, BG-SINDy maintains low coeﬃcient errors when suﬃcient data are available.

2.36

2.20

2.25

50000 (6.408%)

2.25

2.60

4.67

25000 (3.204%)

2.59

2.31

3.55

10000 (1.282%)

2.37

3.95

5000 (0.641%)

2.75

6.01

1000 (0.128%)

3.23

5.04

6% 5% 4% 3% 2% 1%

5% 10% Noise levels

25%

Coefficient error (%)

Samples

100000 (12.816%)

Figure 3: Heatmap of the coeﬃcient error in the BG-SINDy discovery of the KdV equation across diﬀerent noise levels and sampling data sizes, computed according to Eq. (13). The colormap ranges from light red (low error) to dark red (high error). Gray cells indicate failed discovery of the correct equation structure.

3.2. Modiﬁed Burgers equation with vanishing hyperviscosity This example aims to highlight BG-SINDy’s ability to discover governing equations that contain small-coeﬃcient terms involving high-order derivatives. We consider the modiﬁed Burgers equation augmented with a vanishing hyperviscosity term [58, 59] (i.e., a fourth-order derivative with a small coeﬃcient): ∂u ∂u ∂ 2u ∂ 4u +u − 0.5 2 + ε 4 = 0, ∂t ∂x ∂x ∂x

(x, t) ∈ Ω × [0, T ],

(14)

subject to the initial condition u(x, 0) = cos

x 16

x ∈ Ω,

and periodic boundary conditions, where ε = 1 × 10−3 , Ω = [0, 32π], and T = 100. The dataset is generated by numerically solving Eq. (14) using a fourthorder Runge–Kutta integrator with time step ∆t = 0.1 and Fourier spectral discretization with 4048 modes [17, 56]. Temporal derivatives are computed via ﬁnite diﬀerences, while spatial derivatives are evaluated spectrally. The candidate library comprises the polynomial terms up (p = 0, 1, 2) and all q products of these terms with spatial derivatives ∂∂xuq (q = 1, . . . , 4), yielding a 13

total of 15 terms:

{u | p = 0, 1, 2} ∪ p

∂qu up q ∂x

p = 0, 1, 2;

q = 1, . . . , 4 .

The performance of BG-SINDy on the modiﬁed Burgers equation with vanishing hyperviscosity is summarized in Fig. 2. As shown in Fig. 4(a), the residual increases monotonically as the candidate library is pruned. The selected model, marked by the red circle, corresponds to the elbow point, indicating a pronounced change in the residual. The evolution of term importance during the sparsiﬁcation process is visualized in the heatmap in Fig. 4(b). The validity of the discovered equation is further veriﬁed by numerically integrating the model and comparing the resulting solution with the reference solution shown in Fig. 4(c), yielding a relative L2 error below 0.01%. A comparison with SINDy and PDE-FIND under the same candidate library and derivative estimation is shown in Fig. 4(d). BG-SINDy accurately recovers the modiﬁed Burgers equation with a coeﬃcient error below 0.01%, whereas both SINDy and PDE-FIND identify only the nonlinear and viscous terms and miss the small-coeﬃcient hyperviscosity term. This example demonstrates the robustness of BG-SINDy to small-coeﬃcient high-order derivative terms. 3.3. Modiﬁed KS equation with multiple small-coeﬃcient nonlinear terms This example aims to highlight the ability of BG-SINDy to discover governing equations containing multiple small-coeﬃcient terms from data. We consider the modiﬁed KS equation augmented with multiple small-coeﬃcient nonlinear terms [36] 6 X ∂u ∂u ∂ 2 u ∂ 4 u ∂ +u + 2 + 4 +ε uk = 0, ∂t ∂x ∂x ∂x ∂x k=3

(x, t) ∈ Ω × [0, T ], (15)

with the initial condition u(x, 0) = cos(3x) − 21 sin(x),

x ∈ Ω.

and the periodic boundary condition, where ε = 1 × 10−6 , Ω = [0, 22] and T = 200. The candidate library comprises polynomial terms in u up to degree 10, together with their products with spatial derivatives of u up to 10th order, 14

(a)

(b)

Residual history in sparsification

1 u u2 ux uux u 2 ux uxx uuxx u2uxx uxxx uuxxx u2uxxx uxxxx uuxxxx 2 u uxxxx

Residual history Discovered equation

10−1

ut = − 1.0000uux + 0.5000uxx − 1.0000 × 10−3uxxxx

10−5 10−7 10−9 10−11

151413121110 9 8 7 6 5 4 3 2 1 0 Number of candidate terms

0.4 0.2

15 14 13 12 11 10 9 8 7 6 5 4 3 2 1 Number of candidate terms

Solution of discovered equation 100 1.0

100

0.5

1.5

0.0

1.0

−0.5

0.5

50 x

−1.0 0 50 75 100 0 25 x Reference equation: ut = −uux + 0.5uxx − 1 × 10−3uxxxx

100

−1.0

Absolute error

0.0

1.0

(d)

0.6

100

1.0 0.8

Normalized residual

10−3

Term importance in sparsification

50 x

Method

Discovered equation

Coefficient error

BG-SINDy

ut = − 1.0000uux + 0.5000uxx − 1.0000 × 10−3uxxxx

< 0.01%

SINDy

ut = − 0.9998uux + 0.5019uxx

—

PDE-FIND

ut = − 0.9998uux + 0.5019uxx

—

×10−9 2.0

100

0.0

Figure 4: Results for the modiﬁed Burgers equation. (a) History of residual in the sparsiﬁcation process. The red circle denotes the discovered equation. (b) Heatmap of term importance in the sparsiﬁcation process. (c) Comparison of the solutions obtained from the reference equation and the discovered equation. The relative L2 error is less than 0.01%. (d) Discovered equation by BG-SINDy and comparison with SINDy and PDE-FIND. BGSINDy accurately recovers the reference equation, while SINDy and PDE-FIND fail to identify the small-coeﬃcient term.

resulting in 121 candidate terms: ∂qu p {u | p = 0, 1, . . . , 10} ∪ up q p = 0, 1, . . . , 10; ∂x 15

q = 1, . . . , 10 .

These terms are constructed from the solution u obtained by numerically solving the modiﬁed KS equation (15) using a sixth-order Runge–Kutta scheme with time step ∆t = 0.004 and Fourier spectral discretization with 128 modes [36, 60]. Temporal derivatives are computed using the ﬁnite diﬀerence method, and spatial derivatives are computed using the spectral method. The performance of BG-SINDy on the modiﬁed KS equation with multiple small-coeﬃcient nonlinear terms is summarized in Fig. 5. As shown in Fig. 5(a), the residual increases monotonically as the candidate library is pruned. The selected model, marked by the red circle, corresponds to the elbow point, indicating a pronounced change in the residual. The validity of the discovered equation is further veriﬁed by numerically integrating the model and comparing the resulting solution with the reference solution shown in Fig. 5(b), yielding a relative L2 error below 0.01%. A comparison with SINDy and PDE-FIND under the same candidate library and derivative estimation is shown in Fig. 5(c). BG-SINDy accurately recovers the modiﬁed KS equation with a coeﬃcient error of 0.02%. In contrast, SINDy only recovers the two terms with small coeﬃcients, and PDE-FIND fails to identify any of the small-coeﬃcient terms. This example demonstrates the robustness of BG-SINDy to multiple small-coeﬃcient terms. 3.4. Two-dimensional reaction–diﬀusion system This example aims to highlight the ability of BG-SINDy to discover coupled governing equations in two spatial dimensions. The two-dimensional reaction–diﬀusion system [61, 62] is a widely studied model that exhibits complex pattern formation through the interplay of nonlinear reaction kinetics and diﬀusion. Here we consider the reaction–diﬀusion equation with small diﬀusion coeﬃcients 2 ∂u ∂ u ∂ 2u 3 2 2 3 = u + 0.5v − uv + 0.5u v − u + ε + , ∂t ∂x2 ∂y 2 2 ∂v ∂ v ∂ 2v (16) 3 2 2 3 = v − v − 0.5uv − u v − 0.5u + ε + 2 , 2 ∂t ∂x ∂y (x, y, t) ∈ Ω × [0, T ],

(a)

Residual history in sparsification Residual history Discovered equation

Normalized residual

10−1 10−3 10−5

ut = − 1.0000uux − 1.0000uxx − 1.0000uxxxx − 3.0005 × 10−6u2ux − 3.9960 × 10−6u3ux − 4.9998 × 10−6u4ux − 6.0009 × 10−6u5ux

10−7 10−9 10−11 10−13 120

110

100

150

1.5

150

1.5

150

100

0.0

100

0.0

100

−1.5

−3.0

10 15 x

Reference equation:

3.0

200

(b) Solution of reference equation

(c)

Number of candidate terms Solution of discovered equation 200 3.0

200

Absolute error

×10−3 1.0 0.8 0.6 0.4 0.2

10 15 20 0 5 10 15 20 x x ut = − uux − uxx − uxxxx − 3 × 10−6u2ux − 4 × 10−6u3ux −5 × 10−6u4ux − 6 × 10−6u5ux

Method

Discovered equation

Coefficient error

BG-SINDy

ut = − 1.0000uux − 1.0000uxx − 1.0000uxxxx − 3.0005 × 10−6u2ux − 3.9960 × 10−6u3ux − 4.9998 × 10−6u4ux − 6.0009 × 10−6u5ux

0.02%

SINDy

ut = − 1.0001uux − 1.0000uxx − 1.0000uxxxx − 5.6881 × 10−6u4ux − 6.5543 × 10−6u5ux

—

PDE-FIND

ut = −0.0595uux

—

0.0

Figure 5: Results for the modiﬁed KS equation. (a) History of residual in the sparsiﬁcation process. The red circle denotes the discovered equation. (b) Comparison of the solutions obtained from the reference equation and the discovered equation. The relative L2 error below 0.01%. (c) Discovered equation by BG-SINDy and comparison with SINDy and PDE-FIND. BG-SINDy accurately recovers the reference modiﬁed KS equation, whereas SINDy only recovers the two terms with small coeﬃcients, and PDE-FIND fails to identify any of the small-coeﬃcient terms.

where Ω = [−1.5, 1.5] × [−1.5, 1.5], ε = 1 × 10−3 , and T = 5. The initial conditions are given by p p 2 2 2 2 u(x, y, 0) = tanh x + y cos 2θ(x + iy) − x + y , (17) p p v(x, y, 0) = tanh x2 + y 2 sin 2θ(x + iy) − x2 + y 2 , where θ(z) is the principal angle of z ∈ C, and periodic boundary conditions are imposed [19]. The reference dataset is generated by numerically solving the 2D reaction– diﬀusion system (16) using the adaptive RK45 solver [63, 64] in time, combined with a Fourier spectral discretization in space. The spatial domain is discretized on a 256 × 256 periodic grid, and the output time interval is ∆t = 0.05. The candidate library consists of all monomials in u and v up to total degree 3, together with the ﬁrst- and second-order spatial derivatives of both variables, yielding a total of 20 candidate terms [19]: i j u v i, j ≥ 0, i + j ≤ 3 ∪ {ux , uy , uxx , uxy , uyy , vx , vy , vxx , vxy , vyy } . The required derivatives for constructing this library are estimated via automatic diﬀerentiation [43, 44, 45] on a neural-network surrogate (a fully connected network with ﬁve hidden layers of 50 neurons each using sinusoidal activation functions) trained on 7,500 points (0.114% of the full dataset) obtained via Latin hypercube sampling [43, 44, 45]. The results for the u and v components are summarized in Fig. 6 and Fig 7, respectively. For each component, Fig. 6(a) and Fig. 7(a) shows the evolution of the residual during the sparsiﬁcation process, with the red circles indicating the selected models. Fig. 6(b) and Fig. 7(b) compares solution snapshots from the reference simulation and the discovered model at t = 1.0 and t = 3.0. The relative L2 errors for the u and v components, evaluated over the entire time domain, are 0.22% and 0.14%, respectively. As shown in Fig. 6(c) and Fig. 7(c), BG-SINDy accurately recovers the reference equations for both components, whereas standard SINDy and PDE-FIND fail to identify the correct equation structures. This example demonstrates the scalability of BG-SINDy to two-dimensional coupled systems.

(a)

Residual history in sparsification

Normalized residual

100

Residual history Discovered equation

10−1

ut = 0.9950u + 0.5001v 3 − 0.9950uv 2 + 0.5023u2v − 0.9956u3 + 1.013 × 10−3uxx + 9.7439 × 10−4uyy

10−2

10−3 17

Solution of reference equation 1.0 1

t = 1.0 y

0.0

−1

0 x

−0.5

−1

0.0

−1

0 x

−0.5

−1 −1

(c)

0 x

−1.0

t = 3.0 y

0.0

1 0

3.0 1.5

−1

−1.0

−1

0 x

0.0

−0.5

−1 −1

0 x

0.0 ×10−2 2.4

0.5

×10−3 6.0 4.5

1.0

0.5

Absolute error

−0.5

−1

−1.0

0.5

1.0 1

Number of candidate terms Solution of discovered equation 1.0

0.5

(b)

−1.0

1.8

1.2 0.6

−1 −1

0 x

0.0

Reference equation: ut = u + 0.5v 3 − uv 2 + 0.5u2v − u3 + 1 × 10−3uxx + 1 × 10−3uyy Method

Discovered equation

Coefficient error

BG-SINDy

ut = 0.9950u + 0.5001v 3 − 0.9950uv 2 + 0.5023u2v − 0.9956u3 + 1.013 × 10−3uxx + 9.7439 × 10−4uyy

0.83%

SINDy

ut = 0.8342u + 0.5017v 3 − 0.8392uv 2 + 0.5025u2v − 0.8389u3

—

PDE-FIND

ut = 0.8349u + 0.5017v 3 − 0.8401uv 2 + 0.5022u2v − 0.8396u3

—

Figure 6: Results for the 2D reaction–diﬀusion equation (u component). (a) History of residual in the sparsiﬁcation process. The red circle denotes the discovered equation. (b) Comparison of solution snapshots from the reference and discovered equations at t = 1.0 and t = 3.0. The relative L2 error, computed over the entire time domain, is 0.22%. 19 (c) Discovered equation by BG-SINDy and comparison with SINDy and PDE-FIND. BGSINDy accurately recovers the reference equation, whereas SINDy and PDE-FIND fail to identify the equation structure.

(a)

Residual history in sparsification

Normalized residual

100

Residual history Discovered equation

10−1

vt = 0.9981v − 0.9984v 3 − 0.4992uv 2 − 0.9979u2v − 0.5001u3 + 9.9491 × 10−4vxx + 9.8530 × 10−4vyy

10−2

10−3 17

Solution of reference equation 1.0 1

t = 1.0 y

0.0

−1

0 x

−0.5

−1

0.0

−1

0 x

−0.5

−1 −1

(c)

0 x

−1.0

t = 3.0 y

0.0

1 0

3.0 1.5

−1

−1.0

−1

0 x

0.0 −0.5

−1 −1

0 x

−1.0

0.0 ×10−3 8.0

0.5

×10−3 6.0 4.5

1.0

0.5

Absolute error

−0.5

−1

−1.0

0.5

1.0 1

Number of candidate terms Solution of discovered equation 1.0

0.5

(b)

6.0

4.0 2.0

−1 −1

0 x

0.0

Reference equation: vt = v − v 3 − 0.5uv 2 − u2v − 0.5u3 + 1 × 10−3vxx + 1 × 10−3vyy Method

Discovered equation

Coefficient error

BG-SINDy

vt = 0.9981v − 0.9984v 3 − 0.4992uv 2 − 0.9979u2v − 0.5001u3 + 9.9491 × 10−4vxx + 9.8530 × 10−4vyy

0.39%

SINDy

vt = 0.8372v − 0.8424v 3 − 0.5016uv 2 − 0.8414u2v − 0.5018u3

—

PDE-FIND

vt = 0.8380v − 0.8432v 3 − 0.5018uv 2 − 0.8424u2v − 0.5018u3

—

Figure 7: Results for the 2D reaction–diﬀusion equation (v component). (a) History of residual in the sparsiﬁcation process. The red circle denotes the discovered equation. (b) Comparison of solution snapshots from the reference and discovered equations at t = 1.0 and t = 3.0. The relative L2 error, computed over the entire time domain, is 0.14%. 20 (c) Discovered equation by BG-SINDy and comparison with SINDy and PDE-FIND. BGSINDy accurately recovers the reference equation, whereas SINDy and PDE-FIND fail to identify the equation structure.

4. Conclusions In this paper, we propose BG-SINDy, a data-driven framework for discovering nonlinear governing equations in multiscale systems with smallcoeﬃcient terms. Motivated by the principle of dominant balance, the method reformulates the ℓ0 -constrained sparse regression (Eq. (3)) as a term-level ℓ2,0 regularized formulation (Eq. (4)) and solves it using a progressive pruning strategy. Candidate terms are ranked according to their relative contributions to the equation balance and progressively pruned, with coeﬃcients reestimated via iterative least-squares regression, enabling the identiﬁcation of dynamically essential terms even when their coeﬃcients are small. The performance of BG-SINDy is illustrated on a set of benchmark problems, including the KdV equation with a small dispersion coeﬃcient, a modiﬁed Burgers equation with vanishing hyperviscosity, a modiﬁed KS equation with multiple small-coeﬃcient terms, and a two-dimensional reaction– diﬀusion system. In these examples, the method is able to recover the underlying partial diﬀerential equations, including terms that are often diﬃcult to identify using conventional sparse-regression approaches. This behavior is associated with a balance-aware term selection strategy, which diﬀers from standard approaches that rely on coeﬃcient magnitude and instead evaluates the contributions of candidate terms to the governing dynamics. More broadly, incorporating dominant-balance concepts into equationdiscovery algorithms may provide a useful perspective for constructing interpretable models of complex physical systems. Future work will explore extensions to other equation-discovery frameworks, including conventional sparse regression and symbolic regression, as well as integration with deep learning approaches for learning governing equations and surrogate models from sparse or noisy data. These directions are particularly relevant for data-driven modeling of challenging multiscale phenomena, such as nonNewtonian constitutive relations [65] and subgrid-scale closures in turbulence [66, 67, 68, 69]. Declaration of competing interest The authors declare that they have no known competing ﬁnancial interests or personal relationships that could have appeared to inﬂuence the work reported in this paper.

Acknowledgements This work was supported by the National Natural Science Foundation of China (NSFC) under the Excellence Research Group Program “Multiscale Problems in Nonlinear Mechanics” (Grant No. 12588201). Lei Zhang was also supported by the NSFC under Grant No. 12202451. Zhenhua Dang gratefully acknowledges Zihao Yang and Junfeng Zhao at Northwestern Polytechnical University for helpful discussions and support. Long Wang also acknowledges valuable discussions with Fengshun Zhang at the Institute of Mechanics, Chinese Academy of Sciences.

References [1] H. Voss, M. J. Bünner, M. Abel, Identiﬁcation of continuous, spatiotemporal systems, Phys. Rev. E 57 (1998) 2820–2823. https://doi.org/ 10.1103/PhysRevE.57.2820. [2] M. Bär, R. Hegger, H. Kantz, Fitting partial diﬀerential equations to space-time dynamics, Phys. Rev. E 59 (1999) 337–342. https://doi. org/10.1103/PhysRevE.59.337. [3] M. Schmidt, H. Lipson, Distilling Free-Form Natural Laws from Experimental Data, Science 324 (5923) (2009) 81–85. https://doi.org/10. 1126/science.1165893. [4] J. Bongard, H. Lipson, Automated reverse engineering of nonlinear dynamical systems, Proc. Natl. Acad. Sci. U.S.A. 104 (24) (2007) 9943– 9948. https://doi.org/10.1073/pnas.0609476104. [5] S. L. Brunton, J. L. Proctor, J. N. Kutz, Discovering governing equations from data by sparse identiﬁcation of nonlinear dynamical systems, Proc. Natl. Acad. Sci. U.S.A. 113 (15) (2016) 3932–3937. https://doi.org/ 10.1073/pnas.1517384113. [6] S. H. Rudy, S. L. Brunton, J. L. Proctor, J. N. Kutz, Data-driven discovery of partial diﬀerential equations, Sci. Adv. 3 (4) (2017) e1602614. https://doi.org/10.1126/sciadv.1602614. [7] K. Champion, B. Lusch, J. N. Kutz, S. L. Brunton, Data-driven discovery of coordinates and governing equations, Proc. Natl. Acad. Sci. U.S.A. 116 (45) (2019) 22445–22451. https://doi.org/10.1073/pnas. 1906995116. [8] L. Fung, U. Fasel, M. Juniper, Rapid Bayesian identiﬁcation of sparse nonlinear dynamics from scarce and noisy data, Proc. R. Soc. A, Math. Phys. Eng. Sci. 481 (2307) (2025) 20240200. https://doi.org/10. 1098/rspa.2024.0200. [9] L. Prandtl, über ﬂüssigkeitsbewegung bei sehr kleiner reibung, in: Verhandlungen des III. Internationalen Mathematiker-Kongresses, Heidelberg 1904, B. G. Teubner, Leipzig, 1905, pp. 484–491.

[10] H. Schlichting, K. Gersten, Boundary-Layer Theory, 9th Edition, Springer, Berlin, 2017. https://doi.org/10.1007/ 978-3-662-52919-5. [11] C. D. Murray, S. F. Dermott, Solar System Dynamics, Cambridge University Press, Cambridge, 2000. https://doi.org/10.1017/ CBO9781139174817. [12] G.-J. Both, S. Choudhury, P. Sens, R. Kusters, DeepMoD: Deep learning for model discovery in noisy data, J. Comput. Phys. 428 (2021) 109985. https://doi.org/10.1016/j.jcp.2020.109985. [13] C. M. Bender, S. A. Orszag, Advanced mathematical methods for scientists and engineers: Asymptotic methods and perturbation theory, Springer, New York, 1999. https://doi.org/10.1007/ 978-1-4757-3069-2. [14] A. H. Nayfeh, Perturbation methods, John Wiley & Sons, Hoboken, 2024. [15] M. H. Holmes, Introduction to perturbation methods, Springer, New York, 2012. [16] J. K. Kevorkian, J. D. Cole, Multiple scale and singular perturbation methods, Vol. 114, Springer, New York, 2012. [17] D. A. Messenger, D. M. Bortz, Weak SINDy for partial diﬀerential equations, J. Comput. Phys. 443 (2021) 110525. https://doi.org/ 10.1016/j.jcp.2021.110525. [18] U. Fasel, J. N. Kutz, B. W. Brunton, S. L. Brunton, Ensemble-SINDy: Robust sparse model discovery in the low-data, high-noise limit, with active learning and control, Proc. R. Soc. A, Math. Phys. Eng. Sci. 478 (2260) (2022) 20210904. https://doi.org/10.1098/rspa.2021. 0904. [19] B. M. de Silva, K. Champion, M. Quade, J.-C. Loiseau, J. N. Kutz, S. L. Brunton, PySINDy: A Python package for the sparse identiﬁcation of nonlinear dynamical systems from data, J. Open Source Softw. 5 (49) (2020) 2104. https://doi.org/10.21105/joss.02104.

[20] L. Zhang, S. Tang, G. He, Learning chaotic systems from noisy data via multi-step optimization and adaptive training, Chaos 32 (12) (2022) 123134. https://doi.org/10.1063/5.0114542. [21] M. Raissi, P. Perdikaris, G. Karniadakis, Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial diﬀerential equations, J. Comput. Phys. 378 (2019) 686–707. https://doi.org/10.1016/j.jcp.2018.10.045. [22] G. E. Karniadakis, I. G. Kevrekidis, L. Lu, P. Perdikaris, S. Wang, L. Yang, Physics-informed machine learning, Nat. Rev. Phys. 3 (6) (2021) 422–440. https://doi.org/10.1038/s42254-021-00314-5. [23] Z. Chen, Y. Liu, H. Sun, Physics-informed learning of governing equations from scarce data, Nat. Commun. 12 (1) (2021) 6136. https: //doi.org/10.1038/s41467-021-26434-1. [24] L. Zhang, G. He, Scale-Decomposed Physics-Informed Neural Networks for singular perturbation problems, Commun. Comput. Phys. 39 (4) (2026) 1267–1298. https://doi.org/10.4208/cicp.OA-2024-0149. [25] S.-M. Udrescu, M. Tegmark, AI Feynman: A physics-inspired method for symbolic regression, Sci. Adv. 6 (16) (2020) eaay2631. https:// doi.org/10.1126/sciadv.aay2631. [26] Y. Chen, Y. Luo, Q. Liu, H. Xu, D. Zhang, Symbolic genetic algorithm for discovering open-form partial diﬀerential equations (SGAPDE), Phys. Rev. Res. 4 (2) (2022) 023174. https://doi.org/10. 1103/PhysRevResearch.4.023174. [27] M. Cranmer, Interpretable Machine Learning for Science with PySR and SymbolicRegression.jl (2023). arXiv:2305.01582. [28] Z. Liu, Y. Wang, S. Vaidya, F. Ruehle, J. Halverson, M. Soljacic, T. Y. Hou, M. Tegmark, KAN: Kolmogorov–Arnold networks, in: The Thirteenth International Conference on Learning Representations, 2025. [29] Z. Zhou, Z. Xu, Y. Liu, S. Wang, asKAN: Active subspace embedded Kolmogorov–Arnold network, Neural Netw. 195 (2026) 108280. https: //doi.org/10.1016/j.neunet.2025.108280.

[30] M. Du, Y. Chen, Z. Wang, L. Nie, D. Zhang, Large language models for automatic equation discovery of nonlinear dynamics, Phys. Fluids 36 (9) (2024) 097121. https://doi.org/10.1063/5.0224297. [31] H. Xu, Y. Chen, R. Cao, T. Tang, M. Du, J. Li, A. H. Callaghan, D. Zhang, Generative discovery of partial diﬀerential equations by learning from math handbooks, Nat. Commun. 16 (1) (2025) 10255. https://doi.org/10.1038/s41467-025-65114-2. [32] D. R. Gurevich, M. R. Golden, P. A. Reinbold, R. O. Grigoriev, Learning ﬂuid physics from highly turbulent data using sparse physics-informed discovery of empirical relations (SPIDER), J. Fluid Mech. 996 (2024) A25. https://doi.org/10.1017/jfm.2024.813. [33] C. J. Wareing, A. T. Roy, M. Golden, R. O. Grigoriev, S. M. Tobias, Data-driven discovery of the equations of turbulent convection, Geophys. Astrophys. Fluid Dyn. 119 (3–4) (2025) 550–582. https://doi.org/10. 1080/03091929.2025.2509469. [34] M. Golden, R. O. Grigoriev, J. Nambisan, A. Fernandez-Nieves, Physically informed data-driven modeling of active nematics, Sci. Adv. 9 (27) (2023) eabq6120. https://doi.org/10.1126/sciadv.abq6120. [35] M. Golden, K. Satapathy, D. Psaltis, Scalable Discovery of Fundamental Physical Laws: Learning Magnetohydrodynamics from 3D Turbulence Data (2025). arXiv:2501.04094. [36] M. Golden, Scalable Sparse Regression for Model Discovery: The Fast Lane to Insight (2024). arXiv:2405.09579. [37] Z. Li, H. Yuan, W. Han, Y. Hou, H. Li, H. Ding, Z. Jiang, L. Yang, Bilevel identiﬁcation of governing equations for nonlinear physical systems, Nat. Comput. Sci. 5 (6) (2025) 456–466. https://doi.org/10.1038/ s43588-025-00804-x. [38] E. J. Hinch, Perturbation Methods, Cambridge University Press, Cambridge, 1991. https://doi.org/10.1017/CBO9781139172189. [39] J. Kevorkian, J. D. Cole, Multiple Scale and Singular Perturbation Methods, Springer, New York, 1996. https://doi.org/10.1007/ 978-1-4612-3968-0. 26

[40] J. L. Callaham, J. V. Koch, B. W. Brunton, J. N. Kutz, S. L. Brunton, Learning dominant physical processes with data-driven balance models, Nat. Commun. 12 (1) (2021) 1016. https://doi.org/10.1038/ s41467-021-21331-z. [41] K. Fukami, T. Murata, K. Zhang, K. Fukagata, Sparse identiﬁcation of nonlinear dynamics with low-dimensionalized ﬂow representations, J. Fluid Mech. 926 (2021) A10. https://doi.org/10.1017/jfm.2021. 697. [42] K. Kaheman, J. N. Kutz, S. L. Brunton, SINDy-PI: a robust algorithm for parallel implicit sparse identiﬁcation of nonlinear dynamics, Proc. R. Soc. A, Math. Phys. Eng. Sci. 476 (2234) (2020) 20200279. https: //doi.org/10.1098/rspa.2020.0279. [43] A. G. Baydin, B. A. Pearlmutter, A. A. Radul, J. M. Siskind, Automatic Diﬀerentiation in Machine Learning: a Survey, J. Mach. Learn. Res. 18 (153) (2018) 1–43. [44] H. Xu, H. Chang, D. Zhang, DL-PDE: Deep-learning based data-driven discovery of partial diﬀerential equations from discrete and noisy data, Commun. Comput. Phys. 29 (3) (2021) 698–728. https://doi.org/ 10.4208/cicp.OA-2020-0142. [45] K. Kaheman, S. L. Brunton, J. Nathan Kutz, Automatic diﬀerentiation to simultaneously identify nonlinear dynamics and extract noise probability distributions from data, Mach. Learn.: Sci. Technol. 3 (1) (2022) 015031. https://doi.org/10.1088/2632-2153/ac567a. [46] B. K. Natarajan, Sparse Approximate Solutions to Linear Systems, SIAM J. Comput. 24 (2) (1995) 227–234. https://doi.org/10.1137/ S0097539792240406. [47] X. Cai, F. Nie, H. Huang, Exact top-k feature selection via ℓ2,0 -norm constraint, in: Proceedings of the Twenty-Third International Joint Conference on Artiﬁcial Intelligence, Vol. 13, 2013, pp. 1240–1246. [48] T. Pang, F. Nie, J. Han, X. Li, Eﬃcient Feature Selection via ℓ2,0 -norm Constrained Sparse Regression, IEEE Trans. Knowl. Data Eng. 31 (5) (2019) 880–893. https://doi.org/10.1109/TKDE.2018.2847685. 27

[49] P. Lemos, N. Jeﬀrey, M. Cranmer, S. Ho, P. Battaglia, Rediscovering orbital mechanics with machine learning, Mach. Learn.: Sci. Technol. 4 (4) (2022) 045021. https://doi.org/10.1088/2632-2153/acfa63. [50] A. Paszke, S. Gross, F. Massa, A. Lerer, J. Bradbury, G. Chanan, T. Killeen, Z. Lin, N. Gimelshein, L. Antiga, A. Desmaison, A. Kopf, E. Yang, Z. DeVito, M. Raison, A. Tejani, S. Chilamkurthy, B. Steiner, L. Fang, J. Bai, S. Chintala, PyTorch: An Imperative Style, HighPerformance Deep Learning Library, in: Advances in Neural Information Processing Systems, Vol. 32, 2019, pp. 8024–8035. [51] G. H. Golub, C. F. Van Loan, Matrix Computations, 4th Edition, Johns Hopkins University Press, Baltimore, MD, 2013. [52] D. J. Korteweg, G. de Vries, XLI. On the change of form of long waves advancing in a rectangular canal, and on a new type of long stationary waves, Phil. Mag. 39 (240) (1895) 422–443. https://doi.org/10.1080/ 14786449508620739. [53] H. Xu, D. Zhang, Robust discovery of partial diﬀerential equations in complex situations, Phys. Rev. Res. 3 (3) (2021) 033270. https://doi. org/10.1103/PhysRevResearch.3.033270. [54] R. Stephany, C. Earls, Weak-PDE-LEARN: A weak form-based approach to discovering PDEs from noisy, limited data, J. Comput. Phys. 506 (2024) 112950. https://doi.org/10.1016/j.jcp.2024.112950. [55] İ. Dağ, Y. Dereli, Numerical solutions of KdV equation using radial basis functions, Appl. Math. Modell. 32 (4) (2008) 535–546. https: //doi.org/10.1016/j.apm.2007.02.001. [56] A.-K. Kassam, L. N. Trefethen, Fourth-Order Time-Stepping for Stiﬀ PDEs, SIAM J. Sci. Comput. 26 (4) (2005) 1214–1233. https://doi. org/10.1137/S1064827502410633. [57] M. D. McKay, R. J. Beckman, W. J. Conover, Comparison of Three Methods for Selecting Values of Input Variables in the Analysis of Output from a Computer Code, Technometrics 21 (2) (1979) 239–245. https://doi.org/10.1080/00401706.1979.10489755.

[58] E. Tadmor, Burgers’ equation with vanishing hyper-viscosity, Commun. Math. Sci. 2 (2) (2004) 317–324. https://doi.org/10.4310/ CMS.2004.V2.N2.A9. [59] U. Frisch, S. Kurien, R. Pandit, W. Pauls, S. S. Ray, A. Wirth, J.-Z. Zhu, Hyperviscosity, galerkin truncation, and bottlenecks in turbulence, Phys. Rev. Lett. 101 (4) (2008) 144501. https://doi.org/10.1103/ PhysRevLett.101.144501. [60] J. C. Butcher, Implicit Runge-Kutta Processes, Math. Comp. 18 (85) (1964) 50–64. https://doi.org/10.1090/ S0025-5718-1964-0159424-9. [61] A. M. Turing, The chemical basis of morphogenesis, Philos. Trans. R. Soc. Lond. B Biol. Sci. 237 (641) (1952) 37–72. https://doi.org/10. 1098/rstb.1952.0012. [62] J. D. Murray, Mathematical Biology II: Spatial Models and Biomedical Applications, 3rd Edition, Springer, New York, 2003. [63] J. R. Dormand, P. J. Prince, A Family of Embedded Runge-Kutta Formulae, J. Comput. Appl. Math. 6 (1) (1980) 19–26. https://doi.org/ 10.1016/0771-050X(80)90013-3. [64] P. Virtanen, R. Gommers, T. E. Oliphant, M. Haberland, T. Reddy, D. Cournapeau, E. Burovski, P. Peterson, W. Weckesser, J. Bright, S. J. van der Walt, M. Brett, J. Wilson, K. J. Millman, N. Mayorov, A. R. J. Nelson, E. Jones, R. Kern, E. Larson, C. J. Carey, İ. Polat, Y. Feng, E. W. Moore, J. VanderPlas, D. Laxalde, J. Perktold, R. Cimrman, I. Henriksen, E. A. Quintero, C. R. Harris, A. M. Archibald, A. H. Ribeiro, F. Pedregosa, P. van Mulbregt, SciPy 1.0 Contributors, SciPy 1.0: Fundamental Algorithms for Scientiﬁc Computing in Python, Nat. Methods 17 (3) (2020) 261–272. https://doi.org/10. 1038/s41592-019-0686-2. [65] L. Zhao, Z. Li, B. Caswell, J. Ouyang, G. E. Karniadakis, Active learning of constitutive relation from mesoscopic dynamics for macroscopic modeling of non-Newtonian ﬂows, J. Comput. Phys. 363 (2018) 116–127. https://doi.org/10.1016/j.jcp.2018.02.039.

[66] L. Zanna, T. Bolton, Data-Driven Equation Discovery of Ocean Mesoscale Closures, Geophys. Res. Lett. 47 (17) (2020) e2020GL088376. https://doi.org/10.1029/2020GL088376. [67] B. Sanderse, P. Stinis, R. Maulik, S. E. Ahmed, Scientiﬁc machine learning for closure models in multiscale problems: A review, Found. Data Sci. 7 (1) (2025) 298–337. https://doi.org/10.3934/fods.2024043. [68] B. Choi, M. Ugliotti, M. Reynoso, D. R. Gurevich, R. O. Grigoriev, Data-driven modeling of multiscale phenomena with applications to ﬂuid turbulence, Phys. Rev. E (2026). https://doi.org/10.1103/ rdj9-cjm9. [69] K. Jakhar, Y. Guan, P. Hassanzadeh, Analytical and AI-Discovered Stable, Accurate, and Generalizable Subgrid-Scale Closure for Geophysical Turbulence, Phys. Rev. Lett. 136 (2026) 064201. https: //doi.org/10.1103/v28b-5qmp.

Record · ID 120530 · SHA-256 097e4918b29fa68e
Retrieved via Conceptio — every document is proof-bundled with source, license, and retrieval metadata.