Parameter Estimation and Adaptive Solution of the Leray-Burgers Equation using Physics-Informed Neural Networks

We consider a problem of computing the solution $v:[0,T]\times\Omega\to\mathbb{R}$ of an evolution equation

	$\displaystyle v_{t}(t,x)+\mathcal{N}_{\alpha}\,[v](t,x)=0,$	$\displaystyle\quad\forall(t,x)\in[0,T]\times\Omega,$		(2.1)
	$\displaystyle v(0,x)=v_{0}(x),\quad$	$\displaystyle\quad\forall x\in\Omega,$

where $\mathcal{N}_{\alpha}$ is a nonlinear differential operator acting on $v$ with a small constant parameter $\alpha>0$,

$$\mathcal{N}_{\alpha}[v]=vv_{x}+\alpha^{2}v_{x}v_{xx}.$$

(2.2)

Here, the subscripts of $v$ mean partial derivatives in $t$ and $x$, $\Omega\subset\mathbb{R}$ is a bounded domain, $T$ denotes the final time and $v_{0}:\Omega\to\mathbb{R}$ is the prescribed initial data. Although the methodology allows for different types of boundary conditions, we restrict our discussion to Dirichlet or periodic cases and prescribe the boundary data as

$$v_{b}(t,x)=v(t,x),\quad\forall(t,x)\in[0,T]\times\partial\Omega,$$

where $\partial\Omega$ denotes the boundary of the domain $\Omega$.

Equation (2.1) is called the Leray-Burgers equation (LB). It is also known as Burgers-$\alpha$, connectively filtered Burgers equation, Leray regularized reduced order model, etc., in literature. Bhat and Fetecau [2] introduced (2.1) as a regularized approximation to the inviscid Burgers equation

$$v_{t}+vv_{x}=0.$$

(2.3)

They considered a special smoothing kernel associated with the Green function of the Helmholtz operator

$$u_{\alpha}=\mathcal{H}^{-1}_{\alpha}v=(I-\alpha^{2}\partial_{x}^{2})^{-1}v,\ \ (I=identity),$$

where $\alpha>0$ is interpreted as the characteristic wavelength scale below which the smaller physical phenomena are averaged out (see, for example, [8]) and it accelerates energy decay. Applying the smoothing kernel to the convective term in (2.3) yields

$$v_{t}+u_{\alpha}v_{x}=0,$$

(2.4)

where $v=v(t,x)$ is a vector field and $u_{\alpha}$ is the filtered vector field. The filtered vector $u_{\alpha}$ is smoother than $v$ and the equation (2.4) is a nonlinear Leray-type regularization [12] of the inviscid Burgers equation. Here and in the following, we abuse the notation of the filtered vector $u_{\alpha}$ with $u$. If we express the equation (2.4) in the filtered vector $u$, it becomes a quasilinear evolution equation that consists of the inviscid Burgers equation plus $\mathcal{O}(\alpha^{2})$ nonlinear terms [2, 3, 4]:

$$u_{t}+uu_{x}=\alpha^{2}(u_{txx}+uu_{xxx}).$$

(2.5)

In this paper, we follow Zhao and Mohseni [21] to expand the inverse Helmholz operator in $\alpha$ to higher orders of the Laplacian operator:

$$(1-\alpha^{2}\Delta)^{-1}=1+\alpha^{2}\Delta+\alpha^{4}\Delta^{2}+\cdots\ \mbox{if}\ \alpha\lambda_{\mathrm{max}}<1,$$

where $\lambda_{\mathrm{max}}$ is the highest eigenvalue of the discretized operator $\Delta$. Then we can write (2.4) in the unfiltered vector fields $v$ to obtain the equation (2.1)-(2.2) with $\mathcal{O}(\alpha^{4})$ truncation error.

For smooth initial data $v(0,x)$ that decrease at least at one point (so there exists $y$ such that $v_{x}(0,y)<0$), the classical solution $v(t,x)$ of the inviscid Burgers equation (when $\alpha=0$) fails to exist beyond a specific finite break time $T_{s}>0$. It is because the characteristics of the inviscid equation intersect in finite time. The Leray-Burgers equation bends the characteristics to make them not intersect each other, avoiding any finite-time intersection and remedying the finite-time breakdown [2, 4]. So, the Leray-Burgers equation possesses a classical solution globally in time for smooth initial data for $\alpha>0$ [2]:

Furthermore, the Leray-Burgers solution $u_{\alpha}(t,x)$ with initial data $u_{\alpha}(0,x)=\mathcal{H}^{-1}_{\alpha}v_{0}(x)$ for $v_{0}\in W^{2,1}(\mathbb{R})$ converges strongly, as $\alpha\rightarrow 0^{+}$, to a global weak solution $v(t,x)$ of the following initial value problem for the inviscid Burgers equation (Theorem 2 in [2]):

$\displaystyle v_{t}+\frac{1}{2}\left(v^{2}\right)_{x}=0\ \ \ \mbox{with}\ \ v(0,x)=v_{0}(x).$

Bhat and Fetecau [2] found numerical evidence that the chosen weak solution in the zero-$\alpha$ limit satisfies the Oleinik entropy inequality, making the solution physically appropriate. The proof relies on uniform estimates of the unfiltered velocity $v$ rather than the filtered velocity $u$. It made possible the strong convergence of the Leray-Burgers solution to the correct entropy solution of the inviscid Burgers equation. In the context of the filtered velocity $u_{\alpha}$, they also showed that the Leray-Burgers equation captures the correct shock solution of the inviscid Burgers equation for Riemann data consisting of a single decreasing jump [4]. However, since $u_{\alpha}$ captures an unphysical solution for Riemann data comprised of a single increasing jump, it was necessary to control the behavior of the regularized equation by introducing an arbitrary mollification of the Riemann data to capture the correct rarefaction solution of the inviscid Burgers equation. With that modification, they extended the existence results to the case of discontinuous initial data $u_{\alpha}\in L^{\infty}$. However, it is still an open problem for the initial data $v_{0}\in L^{\infty}$. In [7], Guelmame, et all, derived a similar regularized equation to (2.5):

$$u_{t}+uu_{x}=\alpha^{2}(u_{txx}+2u_{x}u_{xx}+uu_{xxx}),$$

(2.6)

which has an additional term $2u_{x}u_{xx}$ on the right-handed side. Notice that $u$ in this equation is the filtered vector field in (2.4). When they were establishing the existence of the entropy solution, Guelmame, et all. resorted to altering the equation (2.6), as Bhat and Fetecau had to modify the initial data for their proof in [4]. Analysis in the context of the filtered vector field $u$ appears to induce an additional modification of the equations to achieve the desired results. Working with the actual vector field $v$ may avoid such arbitrary changes.

Equation (2.4) and related models have previously appeared in the literature. We refer [1, 2, 3, 4, 6, 7, 11, 13, 17, 20, 18] for more properties related to the Leray-Burgers equation. The paper [17] explores the role of $\alpha$ in regularizing Proper Orthogonal Decomposition (POD)-Galerkin models for the Kuramoto-Sivashinsky (KS) equation. The $\alpha$-regularization is introduced to enhance the stability and accuracy of these models by applying Helmholtz filtering to the eigenmodes of the quadratic terms. This filtering controls energy transfer between modes, specifically reducing the impact of high-wavenumber modes that contribute to instability, while preserving the system’s key dynamical features. The link between regularization procedures such as Helmholtz regularization and numerical schemes, for example, had been studied in [6, 13]. They argued that, in numerical computations, the parameter $\alpha^{2}$ cannot be interpreted solely as a length scale because it also depends on the numerical discretization scheme chosen. They observed that the choice of $\alpha$ depends on a relation between $\alpha$ and the mesh size that preserves stability and consistency with conservation conditions for the chosen numerical scheme [2, 6, 13]. Also, they found that, for a fixed number of grid points, there is a particular value of $\alpha\approx 0.02$ below which the solution becomes oscillatory (even with continuous initial profiles).

The Leray-type regularization was first introduced by Jean Leray for the Navier-Stokes equations governing incompressible fluid flow [12]. This regularization scheme enhances stability and accuracy by applying a Helmholtz filter to the convective term, whose effect is evident in the Fourier domain:

$$\widehat{\mathcal{H}^{-1}v}=\frac{\widehat{v}}{1+\alpha^{2}k^{2}}.$$

This filtering mechanism regulates the energy transfer among different modes by attenuating high-wavenumber contributions, which are typically associated with instability. The parameter $\alpha$ thus serves as a regularization constant that sets the length scale of filtering, suppressing smaller spatial scales while preserving the dominant dynamical behavior of the flow.

In the realm of stochastic partial differential equations, Leray regularization has proven effective in enhancing the stability of reduced order models (ROMs), particularly for convection-dominated systems. Iliescu et al., [11], explored this in their study of a stochastic Burgers equation driven by linear multiplicative noise. They found that standard Galerkin ROMs (G-ROMs) produce spurious numerical oscillations in convection-dominated regimes, a problem exacerbated by increasing noise amplitude. To counter this, they applied an explicit spatial filter to the convective term creating the Leray ROM (L-ROM). This approach significantly mitigates oscillations, yielding more accurate and stable solutions compared to the G-ROM, especially under stochastic perturbations. The L-ROM’s robustness to noise variations suggests that Leray regularization may help preserve statistical properties or conserved quantities, such as energy or moments of the solution, in a stochastic context. This extends the utility of Leray regularization beyond deterministic settings, offering a practical tool for modeling complex stochastic dynamics while maintaining numerical fidelity.

Following the original work of Raissi et al. [15, 14], we use a Physics-Informed Neural Network (PINN) to determine physically meaningful $\alpha$-values closely approximating the entropy solutions to the inviscid Burgers equation. For the inverse problem, Alpha2net in Figure 1 is not used because we are looking for a fixed value of $\alpha$. The PINN enforces the physical constraint,

$$\mathcal{F}(t,x):=v_{t}+\lambda_{1}vv_{x}+\lambda_{2}v_{x}v_{xx}$$

on the MLP surrogate $\hat{v}(t,x)=v_{\theta}(t,x;\xi)$, where $\theta=\theta(W,b)$ denotes all parameters of the network (weights $W$ and biases $b$) and $\xi=(\lambda_{1},\lambda_{2})$ the physical parameters in (4.1), acting directly in the loss function

$$\mathcal{L}(\theta,\xi)=\mathcal{L}_{d}(\theta,\xi)+\mathcal{L}_{r}(\theta,\xi),$$

(4.4)

where $\mathcal{L}_{d}$ is the loss function on the available measurement data set that consists in the mean-squared-error (MSE) between the MLP’s predictions and training data and $\mathcal{L}_{r}$ is the additional residual term quantifying the discrepancy of the neural network surrogate $v_{\theta}$ with respect to the underlying differential operator in (4.1). Note that $\mathcal{L}_{d}=w_{0}\mathcal{L}_{0}+w_{b}\mathcal{L}_{b}$ with $w_{0}=w_{1}=1$ in Figure 1. We define the data residual at $(t_{i},x_{i},v_{i})$ in $D$:

$$\mathcal{R}_{d,\theta}(t_{i},x_{i};\xi):=v_{\theta}(t_{i},x_{i};\xi)-v_{i},$$

and the PDE residual at $(t_{i},x_{i})$ in $D$:

$$\mathcal{R}_{r,\theta}(t_{i},x_{i};\xi):=\partial_{t}v_{\theta}+\lambda_{1}v_{\theta}\partial_{x}v_{\theta}+\lambda_{2}\partial_{x}v_{\theta}\partial_{xx}v_{\theta},$$

where $v_{\theta}=v_{\theta}(t,x;\xi)$. Then, the data loss and residual loss functions in (4.4) can be written as

	$\displaystyle\mathcal{L}_{d}(\theta,\xi)$	$\displaystyle=\frac{1}{N_{d}}\sum_{i=1}^{N_{d}}\Big{|}\mathcal{R}_{d,\theta}(t_{i},x_{i};\xi)\Big{|}^{2},$
	$\displaystyle\mathcal{L}_{r}(\theta,\xi)$	$\displaystyle=\frac{1}{N_{r}}\sum_{i=1}^{N_{r}}\Big{|}\mathcal{R}_{r,\theta}(t_{i},x_{i};\xi)\Big{|}^{2}.$

The goal is to find the network and physical parameters $\theta$ and $\xi$ that minimize the loss function (4.4):

$$(\theta^{*},\xi^{*})=\underset{\theta\in\Theta,\xi\in\Xi}{\arg\min}\,\mathcal{L}(\theta,\xi)$$

over an admissible set $\Theta$ and $\Xi$ of training network parameters $\theta$ and $\xi$, respectively.

In practice, given the set of scattered data $v_{i}=v(t_{i},x_{i})$, the MLP takes the coordinate $(t_{i},x_{i})$ as input and produces output vectors $v_{\theta}(t_{i},x_{i};\xi)$ that have the same dimension as $v_{i}$. The PDE residual $\mathcal{R}_{r,\theta}(t,x;\xi)$ forces the output vector $v_{\theta}$ to comply with the physics imposed by the LB equation. The PDE residual network takes its derivatives with respect to the input variables $t$ and $x$ by applying the chain rule to differentiate the compositions of functions using the automatic differentiation integrated into TensorFlow. The residual of the underlying differential equation is evaluated using these gradients. The data loss and the residual loss are trained using input from across the entire domain of interest.

We consider the inviscid Burgers equation (2.3) with some standard Riemann initial data of the form

$$v_{0}(x)=\left\{\begin{array}[]{cc}v_{L},&x\leq 0\\ v_{M},&0<x\leq 1\\ v_{R},&x>1\end{array}\right.$$

We used the conservative upwind difference scheme to generate training data. For each initial profile, we computed $256\times 101=25856$ data points throughout the entire spatiotemporal domain. We modified the code in [14] and, for each case, performed ten computational simulations with 2000 training data randomly sampled for each computation. We adopted the Limited-Memory BFGS (L-BFGS) optimizer with a learning rate of 0.01 to minimize MSE (4.4). When the L-BFGS optimizer diverged, we preprocessed with the ADAM optimizer and finalized the optimization with the L-BFGS. We manually checked with random sets of hyperparameters by training the algorithm and selected the best set of parameters that fit our objective, 8 hidden layers and 20 units per layer, $10000$ epochs. We trained the other models with the same parameters, which might not be the best but reasonable fit for them. One remark is that our problem is identifying the model parameter $\alpha$ rather than inferencing solutions, and it is unnecessary to consider physical causality in our loss function (4.4) as pointed out in [19].

Upon training, the network is calibrated to predict the entire solution $v(t,x)$, as well as the unknown parameters $\theta$ and $\xi$. Along with the relative $L^{2}$-norm of the difference between the exact solution and the corresponding trial solution

$$E_{r}(=E_{r}(\hat{v})):=\frac{||v-\hat{v}||_{2}}{||v||_{2}},$$

we used the absolute error of $\lambda_{1}$,

$$\epsilon(\lambda_{1})=|1-\lambda_{1}|$$

in determining the validity of each computational result.

The practical range of $\alpha$ was determined by ensuring the relative $L^{2}$ error $E_{r}$ remained below $10^{-2}$ while $\epsilon(\lambda_{1})<0.01$, aligning the LB solution with the inviscid Burgers entropy solution. The results show that the $\alpha$ value depends on the initial data, with the effective range of $\alpha$ being between 0.01 and 0.05 for continuous initial profiles and between 0.01 and 0.03 for discontinuous initial profiles.

When appropriate, we will also measure the averaged relative $L^{2}$ error in time,

$$\bar{E}_{r}=\frac{1}{T}\int_{0}^{T}\frac{||v-\hat{v}||_{2}}{||v||_{2}}\,dt.$$

In this section, we consider the following viscous Burgers equation for a training data set:

$$v_{t}+vv_{x}=\nu v_{xx},\ \ \ \forall(t,x)\in(0,T]\times\Omega$$

(4.8)

with

		$\displaystyle\mathrm{(A)}\ \ \left\{\begin{array}[]{cc}\nu=\frac{0.01}{\pi},\ T=1,\Omega=[-1,1]&\\ v(0,x)=-\sin(\pi x),\ \forall x\in\Omega&\\ v(t,-1)=v(t,1)=0,\ \forall t\in[0,1]\end{array}\right.\hskip 21.68121pt$
	$\displaystyle\mathrm{or}\,\,$	$\displaystyle\mathrm{(B)}\ \ \left\{\begin{array}[]{cc}\nu=0.07,\ T\approx 0.4327,\ \Omega=[0,2\pi]&\\ v(0,x)=-2\nu\frac{\phi^{\prime}(x)}{\phi(x)}+4,\ \forall x\in\Omega&\\ \phi(x)=\exp\left(\frac{-x^{2}}{4\nu}\right)+\exp\left(\frac{-(x-2\pi)^{2}}{4\nu}\right).\end{array}\right.$

The corresponding LB equation is

$$v_{t}+vv_{x}=-\alpha^{2}v_{x}v_{xx}.$$

For the initial and boundary data (A), Rudy, et al., [16] proposed the data set that can correctly identify the viscous Burgers equation solely from time series data. It contains 101-time snapshots of a solution to the Burgers equation with a Gaussian initial condition propagating into a traveling wave. Each snapshot has 256 uniform spatial grids. For our experiment, we adopt the data set prepared by Raissi, et all, in [15, 14] based on [16], $101\times 256=25856$ data points, generated from the exact solution to (4.8). For training, $N_{d}=2000$ collocation points are randomly sampled and we use the L-BFGS optimizer with a learning rate of 0.8. The average of ten experiments is $\alpha=0.0158$ with $E_{r}\approx 3.8\times 10^{-2}$. The computational simulation shows that the equation develops a shock properly (Figure 5). Note that $\nu=0.01/\pi\approx 12.7\alpha^{2}$.

For the initial and periodic boundary condition (B), we generate $256\times 500=128000$ training data from the exact solution formula for the whole dynamics over time. With $N_{d}=2000$ training data, the PINN diverges frequently. We experiment with the model with 4000 or more data points to determine an appropriate number of training data. $L^{2}$ error remains around $10^{-2}$ for all cases, which does not provide a clear cut. So, we use the absolute error of $\lambda_{1}$ to determine the appropriate number of training data. For each case of $N_{d}$, we perform the computation 5 to 10 times (Table 2).

As $N_{d}$ increases, the results get better and need to do until it reaches the upper limit. Errors between 14000 and 18000 look better than other ranges. More than 18000 does not seem to improve the results. $N_{d}=16000$ (12.5% of total data) is chosen. More than this does not seem to be better. More likely almost the same. The average of ten computations is $\alpha\approx 0.0894$ with $E_{r}=1.64\times 10^{-2}$. Observe that $\nu=0.07\approx 8.8\alpha^{2}$.

Every part of the solution for (B) moves to the right at the same speed, which differs from (A) (Fig. 6). In (A), the left side of a peak moves faster than the right side, developing a steeper middle. It resulted in a higher value of $\alpha$ with (B) than with (A).

In summary, we observe that $\nu=0.01/\pi\approx 12.7\alpha^{2}$ with the profile (A) and $\nu=0.07\approx 8.8\alpha^{2}$ with the profile (B). These results demonstrate that the LB equation can capture nonlinear interactions at significantly smaller length scales compared to the viscous Burgers equation. Notably, numerical schemes for the viscous Burgers equation become unstable at lower $\nu$ values, whereas the LB equation maintains stability and convergence under these conditions. This observation will be clearer when we compare the forward inferred solutions of two equations in Section 5 (Part C).

We write Equation (2.1) in the filtered vector $u_{\alpha}=u$, which is a quasilinear evolution equation that consists of the inviscid Burgers equation plus $\mathcal{O}(\alpha^{2})$ nonlinear terms [2, 3, 4]:

$$u_{t}+uu_{x}=\alpha^{2}(u_{txx}+uu_{xxx}).$$

(4.9)

We compute the equation with the same conditions as in the previous corresponding experiments. The results show that the filtered equation (4.9) also tends to depend on the continuity of the initial profile as shown in Table 3.

When initial profiles contain discontinuities, the $\alpha$ values are much smaller than those with continuous initial profiles. Compared to the unfiltered equation (2.1), the $\alpha$ values for the filtered equation (4.9) are smaller, which may cause more oscillation in forward inference.

With the initial profile (II), the parameter $\lambda_{1}$ for the filtered velocity is not close to 1 with $\epsilon(\lambda_{1})\approx 0.1076$ on average. By increasing the number of epochs from 10000 to 50000 we get a better result. $\lambda_{1}$ gets closer to 1 with $\epsilon(\lambda_{1})\approx 0.0531$, slightly better relative error and loss, which makes the solution better at later time. The oscillation near the discontinuity gets reduced. This verifies that $u$ needs very small $\alpha$ values to approximate the inviscid Burgers solution.

Having established $\alpha$’s practical range, we next explore its application in forward inference.

Figure 8 demonstrates that the Leray-Burgers equation effectively captures the shock formation with the continuous initial profile (I) within the range of $0<\alpha<0.05$. As $\alpha\rightarrow 0^{+}$, the LB solution converges to the inviscid Burgers solution (the last graph in Figure 8).

With the discontinuous initial file (II), the Leray-Burgers equation still accurately captures the shock formation within the range of $0.01<\alpha<0.03$ (Figure 9). However, the MLP-based PINN generates spurious oscillations near the discontinuity at the beginning. Although the network quickly recovers and fits the oscillations as time progresses, the oscillations worsen, and nonlinear instability arises as the $\alpha$ scale becomes smaller than 0.01, which leads to the deviation of the network solution from the actual inviscid Burgers solution (the last graph in Figure 9).

In this section, we employ the MLP-based Physics-Informed Neural Network (MLP-PINN) to effectively learn the nonlinear operator $\mathcal{N}_{\alpha}[v]$, wherein $v$ represents the primary variable and $\alpha$ denotes a parameter. Coutinho, et all., [5], introduced the idea of adaptive artificial viscosity that can be learned during the training procedure and does not depend on the a priori choice of artificial viscosity coefficient. Instead of incorporating the parameter $\alpha$ in place of the artificial viscosity as in [5], we set up a dedicated subnetwork, Alpha2Net depicted in Figure 1, to find the optimal $\alpha(t)$ value. The integration of the subnetwork into the main PINN architecture makes PINN train both $v$ and $\alpha$ to achieve a robust fit with the LB equation. Two examples highlight the ability of the LB equation to capture shock and rarefaction waves as well as the corresponding optimal values of $\alpha$, which are presented in Figure 10.

For the computations, we generated $100\times 1000=100000$ training data in the domain $[0,2]\times[-2,4]$ from the corresponding analytical solution for each case. With $N_{0}=N_{b}=1000$, $N_{r}=10000$, and epochs =20000. The first graph presents computational snapshots of the system’s evolution with the initial profile (II) over the time interval [0, 2]. The computational outputs are $\mathcal{L}\approx 4.1\times 10^{-4}$ and the averaged relative $L^{2}$-error in time is around $3.5\times 10^{-2}$ with $\alpha\approx 0.0169$. Note that $\alpha$ is the average of trained values of $\alpha(t)$ over time. The second graph illustrates snapshots of the evolution of a rarefaction wave with the initial profile (IV). The computational outputs are $\mathcal{L}\approx 1.9\times 10^{-4}$ and the averaged relative $L^{2}$ error in time is around $7.9\times 10^{-3}$ with the averaged $\alpha=0.0032$ in time.

When comparing the Leray-Burgers equation (2.1) with the viscous Burgers equation (4.8), the term $\alpha^{2}v_{x}v_{xx}$ in (2.1) serves as a nonlinear regularization mechanism, acting as a substitute for the linear diffusion term in the viscous Burgers equation. Unlike linear diffusion, the $\alpha$ term in Equation (2.1) depends on both the first derivative $v_{x}$ and the second derivative $v_{xx}$, suggesting that its smoothing effect is more pronounced in regions with high gradients, modulated by the parameter $\alpha$. Thus, it is valuable to assess the performance of these two equations in relation to the inviscid Burgers equation.

Both equations are solved using PINNs with consistent training configurations: 20,000 epochs, fixed weights, and identical network architectures (8 hidden layers, 20 neurons per layer). The key metric for comparison is the $L^{2}$ error, which quantifies the difference between the predicted and exact solutions, with lower values indicating better accuracy. Computations provide $L^{2}$ errors for both equations across different values of $\alpha$, with $\nu$ set equal to $\alpha^{2}$ in the viscous Burgers equation (4.8). The averaged $L^{2}$ errors over time are summarized in Table 4.

For $\alpha$ values ranging from 0.025 to 0.033 ($\nu$ from 0.000625 to 0.001089), the LB equation consistently outperforms the viscous Burgers equation in the averaged $L^{2}$ error. The averaged $L^{2}$ error for LB equation remains relatively stable, ranging from $2.8025\times 10^{-2}$ to $3.2280\times 10^{-2}$. In contrast, the Burgers equation exhibits higher errors, ranging from $7.8998\times 10^{-2}$ to $1.2258\times 10^{-1}$. There is no clear monotonic trend, indicating variability in the neural network’s ability to approximate the solution.

These results indicate that for small values of $\nu$, the viscous Burgers equation is prone to developing shocks due to its hyperbolic nature. PINNs may struggle to accurately capture these discontinuities. In contrast, the LB equation, through its nonlinear regularization effect (dependent on $\alpha$), likely smooths these discontinuities, leading to improved accuracy.

The data also suggest a tipping point between $\alpha=0.032$ ($\nu=0.001024$), where the performance of the two models being to shift. A more definitive transition appears to occur between $\alpha=0.033$ and $\alpha=0.035$ ($\nu=0.001225$), at which point the viscous Burgers equation begins to outperform the LB equation. This transition is illustrated in Figure 11.

To further examine the differences in solution behavior, Figure 12 presents heatmaps of the difference between the solutions of the LB equation and the inviscid Burgers equation, as well as the difference between the viscous Burgers equation and the inviscid Burgers equation, for $\alpha=0.025$ ($\nu=0.000625$). The LB equation exhibits a more gradual transition in error distribution, while the viscous Burgers equation shows sharper localized discrepancies along the shock region. This suggests that the nonlinear regularization in LB equation helps mitigate sharp discontinuities, leading to improved prediction accuracy.

In summary, the parameter $\alpha$ (through $\nu=\alpha^{2}$) controls the regularization strength. Smaller values of $\alpha$ correspond to finer scales where regularization enhances accuracy, while larger values increase $\nu$, potentially leading to over-smoothing compared to the standard viscous Burgers equation. In practice, the LB equation may be preferable for smaller length scales (low $\alpha$), while the viscous Burgers equation may be more suitable for larger scales (higher $\alpha$). The transition appears to occur near $\alpha=0.035$. These findings emphasize the interplay between physical regularization, viscosity, and the numerical approximation capabilities of PINNs.