Two approximate symmetry frameworks for nonlinear partial differential equations with a small parameter: Comparisons, relations, approximate solutions

Abstract

The frameworks of Baikov–Gazizov–Ibragimov (BGI) and Fushchich–Shtelen (FS) approximate symmetries are used to study symmetry properties of partial differential equations with a small parameter. In general, it is shown that unlike the case of ordinary differential equations (ODEs), unstable BGI point symmetries of unperturbed partial differential equations (PDEs) do not necessarily yield local approximate symmetries for the perturbed model. While some relations between the BGI and FS approaches can be established, the two methods yield different approximate symmetry classifications. Detailed classifications are presented for two nonlinear PDE families. The second family includes a one-dimensional wave equation describing the wave motion in a hyperelastic material with a single family of fibers. For this model, approximate symmetries can be used to compute approximate closed-form solutions. Wave breaking times are found numerically and using the approximate solutions, which yield comparable results.

1. Introduction

Perturbed differential equations (DE) are differential equations of the form

\begin{equation*} F[u;\,\epsilon]=F_0[u]+\epsilon F_1[u]+\ldots=0,\end{equation*}

involving a small parameter $\epsilon$ and may be considered a perturbation of some given differential equation $F_0[u]=0$ that corresponds to $\epsilon=0$ . Such small perturbations can disturb the Lie symmetry properties of the unperturbed model in the sense that some or all exact point or local symmetries of the unperturbed equation disappear from the symmetry classification of the perturbed PDE. In particular, the perturbed model cannot have more symmetries holding for all values of $\epsilon$ than for a given value of $\epsilon$ , including $\epsilon=0$ . In order to find extra symmetry-like structures for models with a small parameter, several approximate symmetry methods have been developed. Baikov, Gazizov and Ibragimov [4–6] introduced the approximate symmetry transformation method that is based on expanding approximate symmetry generators with respect to the small parameter (the BGI method). A different approach to approximate symmetries, developed by Fushchich and Shtelen (the FS method) [21], combines a perturbation technique with the symmetry group method by expanding the dependent variables in a Taylor series in the small parameter and approximately replacing the original equations by a system of equations that are the coefficients at different powers of the small parameter; the classical Lie symmetry method is consequently applied to obtain symmetries of the new system. The BGI and FS approaches are not equivalent. They have been compared and used to obtain approximate symmetries and approximate solutions for several PDE models [23, 46]. Burde [11] developed a new approach for approximate symmetries by constructing equations that could be reduced by exact transformations to an unperturbed equation and at the same time would coincide approximately with the perturbed equation.

In this paper, we consider the BGI and FS approximate symmetry frameworks for perturbed partial differential equations (PDEs), with the main goal to investigate and compare the BGI and FS approximate symmetry structures of perturbed PDEs in more detail, and study their relations, and relate them to the symmetry structure of the unperturbed equations.

The majority of differential equations involve arbitrary parameters or arbitrary functions. These parameters have physical meanings and assume values in some ranges or belong to certain classes. Thus, to study the symmetry properties of system of differential equations involving arbitrary elements, one needs to investigate what happens to symmetries as these parameters assume special values. Namely, one classifies admitted symmetries depending on the forms of the arbitrary elements. At the first step towards the classification of symmetries of a differential equation-based model involving arbitrary elements, it is essential to find equivalence transformations for this model. Equivalence transformations map the given DE system to another differential equation from the same general class [26, 30]. Lie groups of equivalence transformations can be sought systematically using techniques similar to those used to determine admitted Lie symmetry groups. The first group classification was given by Sophus Lie for a class of linear second-order partial differential equations in two independent variables [29]. Later, Ovsiannikov [35] found a complete classification for the nonlinear diffusion equation. More recently, classification problems for PDEs and ODEs were widely considered by many researchers (see, e.g., [1, 19, 20, 24, 27, 31, 45]).

The main set of partial differential equations used in this paper as testbeds and physical examples of interest for approximate symmetry computations and comparisons are one-dimensional nonlinear wave equations. Extremely broad literature exists dedicated to various analytical methods and studies of such models, including the studies of symmetry-related properties. In [2], group properties of the nonlinear wave equation

(1.1) \begin{equation} u_{tt}=(f(u)u_x)_x,\qquad u=u(x,t)\end{equation}

were studied (here and below, subscripts denote partial derivatives). Bluman and Cheviakov [8] extended the group classification of (1.1) through a systematic construction of nonlocal symmetries. Point symmetry classifications for the generalized PDE classes

(1.2) \begin{equation} u_{tt}=(f(x,u)u_x)_x,\qquad u_{tt}=\left(f(u)u_x+g(x,u)\right)_x\end{equation}

were considered in [42, 43]. The complete group classification for the PDE family

(1.3) \begin{equation} u_{tt}=f(x,u_x)u_{xx}+g(x,u_x)\end{equation}

was implemented in [7]. Further classifications of different classes of one-dimensional wave equation can be found, for example, in [22, 28, 44, 48].

The study of wave propagation in nonlinear elastic materials has numerous applications in the study of complex materials [38], medical imaging [37], and other areas [10]. Of particular interest are hyperelastic solids, a class of materials that act as ideal elastic solids, and more general classes of models of this type. In particular, the stress within a hyperelastic solid is nonlinearly related to the deformation through a strain energy density function. In many settings, finite displacements in hyperelastic materials in one dimension reduce to nonlinear wave equations of the form

(1.4) \begin{equation}u_{tt}=R(u_x)u_{xx},\end{equation}

where R is the stored energy function [14, 18]. Among others, our goal here is to study and compare exact and approximate symmetry properties and to construct approximate solutions for an important special form of the PDE (1.4) describing fully nonlinear shear waves in a fiber-reinforced material:

(1.5) \begin{equation} u_{tt}=(c^2+\epsilon\, T(u_x))u_{xx},\quad u=u(x,t), \end{equation}

where $c>0$ is the linear wave speed and $\epsilon$ is a real parameter that controls the nonlinear term. The above perturbed model describes waves in anisotropic fiber-reinforced elastic solids, such as biological membranes, in general, and including the cases when fiber effects are relatively small (see, e.g., [18, 36]). While BGI and FS approximate symmetry classifications have been found for some models (see, e.g., [27, 50] and references therein), the exact and approximate (BGI and FS) symmetry classification for PDEs of the form (1.5) has not been considered and compared in detail in the literature.

As hyperbolic systems, wave equations possess characteristic curves (or surfaces) along which the solutions to the equation are simplified. If the characteristic curves intersect, the solution may become multi-valued. This is referred to as a shock or break in the wave and can have the physical meaning of a discontinuous solution [25]. In elementary models, breaks often form when the leading edge of the wave moves slower than the trailing edge, leading to a steeper and eventually vertical wave front. There are multiple approaches to determine the breaking time of a wave. If one can explicitly determine the characteristic curves, then their intersections and therefore the breaking time can be analytically calculated. When characteristic curves are not known explicitly, a breaking criterion can be developed based on the shape or velocity of the wave (see [47] and references therein). Tissier et al. [41] used local energy dissipation in addition to a maximal slope criterion to determine the breaking time. In this paper, we use approximate solution and finite-difference numerical simulations to find the approximate and numerical wave breaking times for the perturbed one-dimensional wave equation (1.5) with $T(u_x)=u_x^2$ .

The paper is organized as follows.

In Section 2, we review the framework of Lie point and local symmetries of perturbed and unperturbed PDEs in detail necessary for the further analysis. We also provide a minimal introduction to BGI [4–6] and FS [21] approximate transformations and approximate symmetries for PDE models with a small parameter.

In Section 3, we study and compare the BGI and FS approximate symmetry frameworks for a perturbed PDE. In particular, we discuss types of BGI and FS approximate symmetries that can arise, including trivial BGI and FS approximate symmetries that have trivial action on O(1) solution components. We classify stable and unstable point symmetries of a family of wave-type equation with respect to the forms of the perturbation defined by an arbitrary function of two variables. For this example, it is observed that FS symmetries are stable in more classification cases than BGI symmetries. It is also shown (Section 3.4) that it is possible to have an unstable point symmetry of an unperturbed PDE that yields a point or a higher-order approximate symmetry neither in BGI nor in FS framework. This makes the PDE case radically different from the ODE case [40]. In Section 3.5, we derive a relation between BGI and FS approximate symmetries, holding in the case of specific forms of symmetry components.

In Section 4, we classify exact and approximate (BGI and FS) point symmetries of the perturbed wave equation (1.5) and use the results to construct approximate solutions of (1.5) in the case of power T,

(1.6) \begin{equation} u_{tt}=\left(c^2+\epsilon u_x^s\right)u_{xx}. \end{equation}

The approximate solutions are compared with numerical solutions for the most physically relevant case $s=2$ and show good agreement. Even for smooth initial conditions, the wave dynamics leads to finite-time wave breaking. We use the numerically computed characteristic curves to approximate the wave breaking time. It is shown that the behavior of the approximate solution arising from approximate symmetry approach, near the time of breaking, can be used to provide an alternative estimate of the wave breaking time. There appears to be an inverse relationship between the wave breaking time and the small parameter $\epsilon$ .

An overview of the results of the current paper and related open problems that would require extensions of approximate symmetry methods are discussed in Section 6.

The current work uses GeM module [15–17] for the Maple symbolic suite to generate and simplify determining equations for equivalence transformations and approximate symmetries.

2. Lie groups of exact and approximate point and local symmetries

We now briefly review the framework of Lie point and local symmetries in comparison with the BGI [4–6] and Fushchich–Shtelen [21] approximate symmetry frameworks for PDE models involving a small parameter.

For simplicity of presentation and analysis, in this work, we consider single PDEs rather than PDE systems.

Let $x = (x^1 ,x^2 ,...,x^n)$ , $n>1$ , and $u=u(x)$ denote respectively the independent variables and the dependent variable of a given problem. We also denote partial derivatives by subscripts: $\partial u/\partial x_j \equiv u_j$ , etc., and the set of all partial derivatives of u of order k by $\partial^k u$ . A general $k^{\rm th}$ -order scalar PDE on u has the form

(2.1) \begin{equation} F_0[u]\equiv F_0(x,u,\partial u,\ldots, \partial^k u)=0,\quad k\geq 1.\end{equation}

We assume that the PDE (2.1) as it stands, or after a point transformation, is of generalized Kovalevskaya type [33], that is can be written in a solved form with respect to the highest pure derivative of u by one of the independent variables. In (2.1) and below, $F_0[u]$ and similar notation denotes differential functions (i.e., functions depending on x, u, and derivatives of u up to some prescribed order s), defined in a domain of the jet space $J^s(x|u)$ . (The latter is viewed as a multi-dimensional space with coordinates $x,u,\partial u,\ldots, \partial^s u$ .) We also note that repeated indices, where appropriate, assume summation.

The solution set $\mathcal{S}$ of the PDE (2.1) in $J^s(x|u)$ , $s\geq k$ , is a hypersurface defined by the relations $F_0[u]=0$ and its differential consequences $\partial F_0[u]=0$ , $\ldots$ , solved for the corresponding differential consequences of the leading derivative, up to the highest order s. Any differential function f[u] can be evaluated on the solution set of (2.1) by substituting the expressions of the leading derivative and its differential consequences into f[u], and the result is denoted by $f[u]|_{\mathcal{S}}$ or $f[u]|_{F_0[u]=0}$ .

2.1. Local symmetries of a PDE

A one-parameter (a) Lie group of transformations

(2.2) \begin{equation}\begin{array}{rcl} (x^*)^i & = & f^i (x,u;\,a) =x^i+a\xi_{0}^i(x,u)+ O\left(a^2\right) ,\quad i = 1,2, \ldots ,n, \\[1ex] u^* &= & g (x,u;\,a)= u +a \eta_0(x,u)+ O\left(a^2\right), \end{array}\end{equation}

in the space of the problem variables (x,u), with the corresponding infinitesimal generator

(2.3) \begin{equation}X^0 = \xi_{0}^i(x,u)\dfrac{\partial }{{\partial x^i }}+ \eta_0(x,u) \dfrac{\partial }{{\partial u }},\end{equation}

defines a point symmetry of the PDE (2.1) when the solution set of (2.1) is invariant under the action of the mapping (2.2); that is, each solution u(x) is mapped to a solution $u^*(x^*)$ of the same PDE. The local condition of this invariance is given by the determining equation

(2.4) \begin{equation}X^{0(k)} F_0[u]\big|_{F_0[u]=0}=0\end{equation}

in terms of the prolonged generator $X^{0(k)}$ . The determining equation splits into an overdetermined PDE system on the unknown symmetry components $\xi_{0}^i, \eta_0$ (see, e.g., [9, 33]).

For a point symmetry (2.2) of the PDE (2.1), the evolutionary (characteristic) form providing the same mapping between solutions has the form

(2.5) \begin{equation}\begin{array}{rcl}(x^*)^i &= &x^i, \quad i = 1, 2, \ldots ,n,\\[1ex] u^* &= &u +a \zeta_0[u]+ O\left(a^2\right),\end{array}\end{equation}

with the evolutionary component $\zeta_0[u]={\eta_0}(x,u)-u_i \xi_{0}^i(x,u)$ , and the generator is given by

(2.6) \begin{equation}\hat{X}^0=\zeta_0[u]\dfrac{\partial}{\partial u}\,.\end{equation}

Local (point, contact and higher-order) transformations and the related point, contact and higher-order local symmetries of the PDE (2.1) generalize (2.5), (2.6) by allowing the evolutionary infinitesimal components $\zeta_0=\zeta_0[u]$ to be general differential functions of u, depending on first- and/or higher-order derivatives of u (see, e.g., [9, 33] and references therein). The invariance condition (2.4) is replaced by

(2.7) \begin{equation}\hat{X}^{0(k)} F_0[u]\big|_{F_0[u]=0}=0.\end{equation}

The prolongation of $\hat{X}^0$ is defined by [33]

(2.8) \begin{equation} \hat{X}^{0(k)}=\zeta_0 \dfrac{\partial}{\partial u}+{\zeta_0}_{i}^{(1)} \dfrac{\partial }{{\partial u_i }}+\ldots + {\zeta_0}_{i_1 i_2 \ldots i_k }^{(k)} \dfrac{\partial}{{\partial u_{i_1 i_2 \ldots i_k } }}, \end{equation}

where the higher-order components are computed using

\begin{equation*} {\zeta_0}_{i}^{(1)} = D_i\zeta_0, \qquad {\zeta_0}_{i_1 i_2 \ldots i_p }^{(p)} =D_{i_p}{\zeta_0}_{i_1 i_2 \ldots i_{p-1} }^{(p-1)}, \end{equation*}

for $i,i_j=1,\ldots,n$ , $p=2,3,\ldots, k$ .

A point or local symmetry of a PDE is trivial if its components vanish on the solution set of the PDE (2.1). Trivial symmetries provide identity transformations (2.2), (2.5). Symmetry generators are commonly simplified modulo trivial symmetries.

The following elementary example will serve as a basis of further examples involving PDEs with a small parameter, and their exact and approximate symmetries.

Example 2.1. Consider a nonlinear wave-type equation [34]

(2.9) \begin{equation}u_{tt}=u_xu_{xx},\quad u=u(x,t).\end{equation}

The exact symmetry generator for the PDE (2.9) is given by

(2.10) \begin{equation}X^0 = \xi_{0}^1 (x,t,u)\dfrac{\partial }{{\partial x }} +\xi_{0}^2 (x,t,u)\dfrac{\partial }{{\partial t }}+ \eta_0 (x,t,u) \dfrac{\partial }{{\partial u }}.\end{equation}

The determining equations (2.4) yield the solution

(2.11) \begin{equation} \xi_{0}^1 = C_4+C_6 x,\quad \xi_{0}^2 = C_3+\left(C_6-\dfrac{C_5}{2}\right) t,\quad \eta_0 = C_1+C_2t+(C_5+C_6) u.\end{equation}

Consequently, the PDE (2.9) admits a six-dimensional Lie algebra of point symmetries spanned by

(2.12) \begin{equation}\begin{array}{rcl} X_1^0 = \dfrac{\partial}{\partial u},\quad X_2^0 = t\dfrac{\partial}{\partial u},\quad X_3^0 = \dfrac{\partial}{\partial t},\quad X_4^0 = \dfrac{\partial}{\partial x},\\[2ex] X_5^0 = u\dfrac{\partial}{\partial u}-\dfrac{t}{2}\dfrac{\partial}{\partial t},\quad X_6^0 =t\dfrac{\partial}{\partial t}+x\dfrac{\partial}{\partial x}+u\dfrac{\partial}{\partial u}.\end{array}\end{equation}

2.2. Exact and approximate symmetries of a PDE with a small parameter

A general first-order perturbation of a PDE (2.1) is a partial differential equation

(2.13) \begin{equation} F[u;\,\epsilon]=F_0[u]+\epsilon F_1[u]=o(\epsilon)\end{equation}

involving a small parameter $\epsilon$ . We assume that the perturbation $F_1[u]$ is regular, in the sense that the Kovalevskaya forms of the unperturbed PDE (2.1) and its perturbation (2.13) have the same leading derivatives.

2.2.1. Exact local symmetries

Exact point and local symmetry generators of (2.13) have the forms

(2.14) \begin{equation} Y=\alpha^{i}(x,u;\,\epsilon)\dfrac{\partial}{\partial x^i}+ \beta(x,u;\,\epsilon)\dfrac{\partial}{\partial u},\qquad \hat{Y}=\zeta[u;\,\epsilon]\dfrac{\partial}{\partial u}.\end{equation}

To find the exact symmetries of (2.13) holding for an arbitrary $\epsilon$ , one solves the determining equations (2.4) or (2.7).

Since (2.13) is a PDE family that includes an arbitrary element $\epsilon$ , the dimension of Lie algebra of point or local symmetries holding for a general $\epsilon$ cannot exceed that for some fixed $\epsilon$ , including $\epsilon=0$ . Therefore, the family (2.13) of perturbed PDEs will admit the same or smaller number of local symmetries than its unperturbed version (2.1).

Example 2.2. We compute exact point symmetries of a perturbed version of the PDE (2.9),

(2.15) \begin{equation} u_{tt}+\epsilon u u_t=u_xu_{xx}\end{equation}

holding for an arbitrary $\epsilon$ . The leading derivative $u_{tt}$ can be chosen for both (2.9) and (2.15). We obtain a Lie algebra of point symmetries spanned by

(2.16) \begin{equation} Y_1=X_3^0= \dfrac{\partial}{\partial t},\quad Y_2=X_4^0=\dfrac{\partial}{\partial x},\quad Y_3=\dfrac{4}{3}X_5^0-\dfrac{1}{3}X_6^0 = -t\dfrac{\partial}{\partial t}-\dfrac{x}{3}\dfrac{\partial}{\partial x} + u \dfrac{\partial}{\partial u},\end{equation}

a three-dimensional subalgebra of the six-dimensional Lie algebra of point symmetries (2.12) of the unperturbed wave equation (2.9).

2.2.2. Baikov–Gazizov–Ibragimov approximate point and local symmetries

Approximate symmetries can be useful for finding additional symmetry-like structures for the perturbed equation (2.13) that are not its exact symmetries but rather preserve (2.13) approximately, up to $o(\epsilon)$ [27]. A one-parameter Lie group of Baikov–Gazizov–Ibragimov (BGI) approximate point transformations with the group parameter a, acting on the (x,u)-space, is given by

(2.17) \begin{equation}\begin{array}{rcll} {(x^*)}^{i} & = &f^i(x,u;\,a,\epsilon) = f_0^i(x,u;\,a)+\epsilon f_1^i(x,u;\,a) + o(\epsilon),& i=1,...,n, \\[1ex] {(u^*)} & = &g(x,u;\,a,\epsilon) = g_0(x,u;\,a)+\epsilon g_1(x,u;\,a) + o(\epsilon),&\end{array}\end{equation}

where, as in (2.2), $f_j^i$ , $g_j$ are sufficiently smooth functions. The generator of the approximate group (2.17), up to $o(\epsilon)$ , is

(2.18) \begin{equation}X = X^0+\epsilon X^1 = \left(\xi_0^i(x,u)+\epsilon \xi_1^i(x,u)\right)\dfrac{\partial}{\partial x^i}+\left(\eta_0(x,u)+\epsilon \eta_1(x,u)\right)\dfrac{\partial}{\partial u}.\end{equation}

Definition 2.1. The approximate group (2.17) of point transformations defines a BGI approximate point symmetry of the PDE (2.13) if it satisfies the approximate invariance condition of (2.13) under the action of (2.18):

(2.19) \begin{equation} (X^{0\,{(k)}}+\epsilon X^{1\,{(k)}})(F_0[u] +\epsilon F_1[u])\bigg |_{F_0[u] +\epsilon F_1[u]=o(\epsilon)}=o(\epsilon).\end{equation}

In (2.19), O(1) and $O(\epsilon)$ terms must vanish independently. It is easy to see that the O(1) term yields the determining equation (2.4) for the invariance of the unperturbed equation (2.1) under a point transformation $X^0$ (2.3). Hence, the following result holds.

Theorem 2.1. If the PDE (2.13) is approximately invariant under an approximate group of BGI point transformations with the generator (2.18) such that its O(1) part $X^0\not\equiv 0$ , then the infinitesimal operator $X^0$ (2.3) is a generator of an exact point symmetry group of the unperturbed PDE (2.1).

The converse of the above result does not always hold. Indeed, as it will be seen in examples below, if $X^0$ (2.3) generates an exact point symmetry group of the unperturbed PDE (2.1), there may be no corresponding BGI transformation (2.17) that approximately preserves the perturbed PDE (2.13). The following definition is important.

Definition 2.2. Suppose that the vector field $X^0$ given by (2.3) is a generator of an exact point symmetry group of the unperturbed PDE (2.1). If the perturbed PDE (2.13) admits an approximate generator X (2.18) with its O(1) part given by $X^0$ , then $X^0$ corresponds to a stable point symmetry of the unperturbed PDE (2.1) (in the BGI sense). Otherwise, it corresponds to an unstable point symmetry of (2.1).

Below in this paper, Definition 2.2 will be used not only for BGI approximate point symmetries, but more generally, for BGI and FS approximate point and local symmetries.

It is straightforward to show that $O(\epsilon)$ term of the determining equation (2.19) leads to the PDEs

(2.20) \begin{equation} X^{1\,{(k)}}F_0[u]\bigg|_{F_0[u]=0}=H[u],\end{equation}

where H is obtained from the coefficients of $\epsilon$ in the expression

(2.21) \begin{equation} -X^{0\,{(k)}}(F_0[u]+\epsilon F_1[u])\bigg|_{F_0[u]+\epsilon F_1[u]=o(\epsilon)}.\end{equation}

It follows that in order to calculate all approximate BGI symmetries (2.17), (2.18) of the perturbed PDE (2.13), one can take the following steps.

1. 1. Compute all generators $X^0$ (2.3) of exact point symmetry groups of the unperturbed PDE (2.1).

2. 2. Use each $X^0$ and the perturbation term $F_1[u]$ to compute H using (2.21).

3. 3. Solve (2.20) to find the components of $X^1$ .

An alternative procedure for the calculation of BGI symmetries involves writing down exact symmetry determining equations for $F[u;\,\epsilon]=0$ (cf. (2.13)), substituting $\zeta[u]=\zeta_0[u]+\epsilon \zeta_1[u]$ , and collecting O(1) and $O(\epsilon)$ coefficients of each split determining equation.

Remark 2.1. The first-order condition (2.20) may (or may not) contain additional conditions on the components $\xi_{0}^i, \eta_0$ of the unperturbed symmetry generator $X^0$ (2.3). This leads to the symmetry generated by $X^0$ being unstable (or respectively, stable). If all symmetries of the equations (2.1) are stable, the perturbed equations (2.13) are said to inherit the symmetry structure of the unperturbed equations [27].

Remark 2.2. Similarly to exact local transformations with generators of the form (2.3), one can define more general local approximate BGI transformations with generators in evolutionary form given by

(2.22) \begin{equation}\hat{X} = \hat{X}^0+\epsilon \hat{X}^1 = \left(\zeta_0[u]+\epsilon \zeta_1[u]\right)\dfrac{\partial}{\partial u}\, .\end{equation}

Approximate local (including point, contact and higher-order) BGI symmetries of the perturbed PDE (2.13) can be found using the same procedure as described above for point approximate BGI symmetries. In particular, the analog of the first-order condition (2.20) takes the form

(2.23) \begin{equation} \left(\zeta_1 \dfrac{\partial}{\partial u}+{\zeta_1}_{i}^{(1)} \dfrac{\partial }{{\partial u_i }}+\ldots + {\zeta_1}_{i_1 i_2 \ldots i_k }^{(k)} \dfrac{\partial}{{\partial u_{i_1 i_2 \ldots i_k } }}\right)F_0[u]\bigg|_{F_0[u]=0}=H[u], \end{equation}

Theorem 2.1, the stability Definition 2.2 and Remark 2.1 concerning stability conditions of approximate symmetries directly carry over to the case of general local BGI symmetries.

2.2.3. Fushchich–Shtelen approximate point and local symmetries

Unlike the BGI approach where the symmetry generator is expanded in a power series in terms of the small parameter, the Fushchich–Shtelen method [21] applies the perturbation technique to the solution u(x) and the given PDE. In particular, the solution is written as

(2.24) \begin{equation} u(x) = v(x)+\epsilon w(x) + o(\epsilon).\end{equation}

Substituting (6.3) into the PDE (2.13) with a small parameter, expanding the result and setting to zero the O(1) and $O(\epsilon)$ terms independently, one obtains a coupled system of equations on v(x) and w(x) without the small parameter, given by

(2.25a) \begin{align} &G_1[v,w]\equiv F_0[v] = 0,\\[-5pt]\nonumber \end{align}

(2.25b) \begin{align} &G_2[v,w]\equiv {F_0}_v w+{F_0}_{v_i}w_i+{F_0}_{v_{ij}}w_{ij}+...+{F_0}_{v_{i_1i_2...i_k}}w_{i_1i_2...i_k}+F_1[v] = 0.\end{align}

It is clear that the first PDE (2.25a) is independent of w, and the second PDE (2.25b) is linear in w, with the linear operator being the Frèchet derivative of the $F_0[u]$ . We refer to equations (2.25) as the Fushchich–Shtelen system for the PDE (2.13). The PDE system (2.25) approximates the given PDE (2.13), in the sense that each exact solution pair (v(x), w(x)) of (2.25) yields an approximate solution (2.24) of the given PDE (2.13) up to the order $o(\epsilon)$ . Note that the Fushchich–Shtelen system (2.25) is in the extended Kovalevskaya form with respect to the similar leading derivatives of v and w as the leading derivatives of u in the original PDE (2.13).

Definition 2.3. A Lie group of point transformations with the group parameter a

(2.26) \begin{equation}\begin{array}{rcll} {(x^*)}^{i} & = &f^i(x,v,w;\,a) = x^i+a\, \lambda^i(x,v,w) + O\left(a^2\right),& i=1,...,n, \\[1ex] {(v^*)} & = &g(x,v,w;\,a) = v+ a\, \phi_1(x,v,w) + O\left(a^2\right),\\[1ex] {(w^*)} & = &h(x,v,w;\,a) = w+ a\, \phi_2(x,v,w) + O\left(a^2\right)\end{array}\end{equation}

with the generator

(2.27) \begin{equation} Z= \lambda^i(x,v,w)\dfrac{\partial}{\partial x^i}+ \phi_1(x,v,w) \dfrac{\partial}{\partial v}+ \phi_2(x,v,w) \dfrac{\partial}{\partial w}\end{equation}

defines a FS approximate point symmetry of the PDE (2.13) if it is an exact Lie point symmetry group of the Fushchich–Shtelen system (2.25).

In a similar manner, a generalized local (point or higher-order) transformation group in the evolutionary form

(2.28) \begin{equation}\begin{array}{ll} {(x^*)}^{i} &= x^i, \qquad i=1,...,n, \\[1ex] {(v^*)} &= v+ a \psi_1[v,w] + O\left(a^2\right),\\[1ex] {(w^*)} &= w+ a \psi_2[v,w] + O\left(a^2\right)\end{array}\end{equation}

with the generator

(2.29) \begin{equation} \hat{Z}=\psi_1[v,w]\dfrac{\partial}{\partial v}+\psi_2[v,w] \dfrac{\partial}{\partial w}\end{equation}

defines a local (point or higher-order) FS approximate symmetry of the PDE (2.13) if it is a local symmetry of the Fushchich–Shtelen system (2.25).

It is important to know whether the FS approximate symmetry structure of a PDE (2.13) with a small parameter is in some sense inherited from exact local symmetries of the unperturbed PDE (2.1). Similarly to the BGI case, one can define stable and unstable symmetries in the Fushchich–Shtelen framework.

Definition 2.4. Suppose $\hat{X}_0=\zeta_0[u]\,{\partial}/{\partial u}$ (2.22) is a generator of an exact local symmetry group of the unperturbed PDE (2.1). If the perturbed PDE (2.13) admits an approximate FS symmetry with generator (2.29) where the v-component $\psi_1[v,w]\equiv\zeta_0[v]$ , then $\hat{X}_0$ corresponds to a stable point symmetry of the unperturbed PDE (2.1) (in the FS sense). Otherwise, it corresponds to an unstable point symmetry of (2.1).

Similarly to the case for BGI approximate symmetries, a FS approximate symmetry of a PDE (2.13) given by (2.29) may be unstable because the second symmetry determining equation for the Fushchich–Shtelen system (2.25)

\begin{equation*} \hat{Z}^{(k)}G_2[v,w]\bigg|_{G_1[v,w]=G_2[v,w]=0} = 0\end{equation*}

could contain additional conditions on the v-component $\psi_1$ .

Remark 2.3. An important feature of the Fushchich–Shtelen approximate symmetry framework is the possibility of existence of approximate FS symmetries where the O(1) component $\psi_1[v,w]$ of the generator (2.29) depends on $O(\epsilon)$ solution component w. An example is provided by the linear PDE $u_{tt}+\epsilon u_t = u_{xx}$ which admits a point FS symmetry (2.29) with

(2.30) \begin{equation}\begin{array}{ll}\psi_1[v,w] = tv+2w +2xt v_x +(x^2+t^2)v_t,\\[1ex]\psi_2[v,w] = \dfrac{1}{2}x^2 v -tw +2xt w_x +(x^2+t^2)w_t\,.\end{array}\end{equation}

Such FS symmetries do not originate from stable local symmetries of the unperturbed PDE $u_{tt} = u_{xx}$ and cannot arise in the BGI framework.

3. BGI and FS approximate symmetries: Properties, connections and examples

3.1. Trivial approximate symmetries

3.1.1. Trivial BGI approximate symmetries

Consider a local BGI approximate transformation with the evolutionary generator (2.22):

\begin{equation*}\hat{X} = \hat{X}^0+\epsilon \hat{X}^1 = \left(\zeta_0[u]+\epsilon \zeta_1[u]\right)\dfrac{\partial}{\partial u}\, .\end{equation*}

The determining equations (2.19) for the generator (2.22) to define a local symmetry of the PDE (2.13) with a small parameter split into the O(1) part (2.4) and $O(\epsilon)$ part (2.20) with H defined by (2.21). Suppose that the O(1) part of the generator vanishes: $\hat{X}_0=0$ . In that case, the O(1) part (2.4) of the approximate symmetry determining equations is satisfied identically, and (2.21) yields $H=0$ . Consequently, the $O(\epsilon)$ part (2.20) of the determining equations (2.19) becomes

(3.1) \begin{equation} \hat{X}^{1\,{(k)}}F_0[u]\bigg|_{F_0[u]=0}=0,\end{equation}

which means that such $\hat{X}^{1}$ must be a local symmetry generator of the unperturbed equation (2.1). The opposite is also true: if $\hat{X}^{0}$ is a local symmetry generator of the unperturbed equation (2.1), then

(3.2) \begin{equation}\hat{X} = \epsilon \hat{X}^0\end{equation}

is a BGI approximate symmetry generator of the perturbed PDE (2.13). In the light of the above, we call a BGI approximate symmetry that has a generator with vanishing O(1) part

(3.3) \begin{equation}\hat{X} = \epsilon \hat{X}^1 = \epsilon \zeta_1[u]\dfrac{\partial}{\partial u}\end{equation}

a trivial BGI approximate symmetry. This triviality relates not to the trivial action of such symmetries but rather to the fact that every local symmetry $\hat{X}^0$ of the unperturbed equation (2.1) is guaranteed to yield a BGI approximate symmetry of the perturbed equation (2.13) having the form (3.2). The local action of a trivial BGI approximate symmetry in the evolutionary form defined by (3.3) is given by

(3.4) \begin{equation}\begin{array}{rcl}(x^*)^i &= &x^i, \quad i = 1, 2, \ldots ,n,\\[1ex] u^* &= &u + a\epsilon\zeta_1[u]+ O\left(a^2\right),\end{array}\end{equation}

with the first Taylor term of the transformation having the order $\sim a\epsilon = o(a,\epsilon)$ .

3.1.2. Trivial FS approximate symmetries

In a parallel fashion, one can define a trivial FS approximate symmetry of the perturbed PDE (2.13) as one for which the local generator (2.29) has a special form with the vanishing transformation component of the O(1) part of the solution $\psi_1=0$ , and $\psi_2[v,w]=\psi_2[v]$ :

(3.5) \begin{equation} \hat{Z}=0+\psi_2[v] \dfrac{\partial}{\partial w}\,.\end{equation}

(The respective trivial point FS symmetry generator (2.27) has $\lambda^i=\lambda^i(x,v)$ , $\phi_1=0$ , $\phi_2=\phi_2(x,v)$ .) For FS local symmetries with the generator of the form (3.5), it is straightforward to show that $\psi_2[v]$ is an evolutionary component of the local symmetry of the unperturbed equation (2.1) generated by

(3.6) \begin{equation} \hat{X}^0= \psi_2[u] \dfrac{\partial}{\partial u}.\end{equation}

Indeed, the action of (3.5) on the first PDE (2.25a) of Fushchich–Shtelen system is trivial, and the action on the linear equation (2.25b) is equivalent to the local symmetry determining equation (2.7) of the unperturbed PDE (2.1). The converse also holds: every local symmetry generator (3.6) of the unperturbed equation yields a trivial FS approximate symmetry generator (3.5).

3.2. Types of approximate symmetries arising in classifications

In the computation of BGI approximate symmetries of a PDE (2.13) with a small parameter, the following three types of symmetries can arise.

1. 1. BGI approximate symmetries with generators (2.22) having $\hat{X}^0 \ne 0$ , $\hat{X}^1=0$ correspond to exact local symmetries of the perturbed equation (2.13) (see Section 2.2.1).

2. 2. BGI approximate symmetries with generators having $\hat{X}^0 = 0$ , $\hat{X}^1 \ne 0$ correspond to trivial BGI approximate symmetries (Section 3.1.1).

3. 3. Genuine BGI approximate symmetries have generators with both $\hat{X}^0\ne 0$ and $\hat{X}^1 \ne 0$ .

For FS approximate symmetries, the following three types can arise.

1. 1. Symmetries with the same action on O(1) solution part v and $O(\epsilon)$ solution part w correspond to exact local symmetries of the perturbed equation (2.13). For example, an exact scaling symmetry with the generator $u\,\partial/\partial u$ admitted by the perturbed equation (2.13) is equivalent to a FS scaling symmetry with the generator $v\,\partial/\partial v + w\,\partial/\partial w$ .

2. 2. Trivial FS approximate symmetries, as defined in Section 3.1.2.

3. 3. Genuine FS approximate symmetries.

Genuine BGI and FS approximate symmetries are the main focus of the approximate symmetry study.

3.3. A computational example: Exact and approximate point symmetry classification for a second-order nonlinear PDE with a small parameter

In this section, we compare exact point symmetries (2.12) of the (1+1)-dimensional wave-type equation (2.9)

\begin{equation*}u_{tt}=u_xu_{xx}\end{equation*}

with BGI and FS approximate point symmetry classifications for the family of perturbed equations

(3.7) \begin{equation} u_{tt}+\epsilon F_1(u, u_t)=u_xu_{xx},\end{equation}

where $F_1(u, u_t)$ is an arbitrary function. The computation is limited to point symmetries for simplicity of presentation. The cases for the arbitrary function $F_1(u, u_t)$ are given modulo equivalence transformations

(3.8) \begin{equation} t=C_1^{-1}\,\tilde{t},\quad x=C_1^{-2/3}\,\tilde{x},\quad u=\tilde{u}+\dfrac{1}{2}\,\epsilon \,\left(1-C_1^{-2}\right)\, C_2 \,\tilde{t}^2,\quad F_1 = C_1^2\tilde{F}_1-\left(C_1^2-1\right)C_2,\end{equation}

where $C_1$ and $C_2$ are arbitrary constants.

Table 1. Stability of point symmetries (2.12) of the wave-type equation in terms of BGI and FS approximate point symmetries of the perturbed PDE (3.7), depending on the form of the arbitrary function $F_1$

The PDE (2.9) admits six exact point symmetries given by (2.12). The BGI approximate point symmetries are computed and classified following the procedure described in Section 2.2.2, and the FS approximate local symmetries are obtained following Section 2.2.3. In particular, the Fushchich–Shtelen system (2.25) for (3.7) is given by

(3.9) \begin{equation} v_{tt}-v_xv_{xx} = 0,\qquad w_{tt}+F_1(v,v_t)-v_xw_{xx}-w_xv_{xx}=0, \end{equation}

where $u(x,t)=v(x,t)+\epsilon w(x,t)$ . For the Fushchich–Shtelen system (3.9), point equivalence transformations are given by (3.8) for x, t and $F_1$ , with dependent variable transformations

(3.10) \begin{equation}v=\tilde{v}, \qquad w=\tilde{w}+\dfrac{1}{2}\,\left(1-C_1^{-2}\right)\, C_2 \,\tilde{t}^2.\end{equation}

The resulting classification is presented in Table 1. In the table, $Q_i$ denote arbitrary functions of their arguments and $a_j$ arbitrary constants. The table is organized as follows. The first column lists evolutionary forms of the six point symmetry generators $X^0_k$ (2.12), $k=1,\ldots, 6$ , of the unperturbed PDE (2.9). The second column lists the forms of the arbitrary function $F_1(u, u_t)$ for which the corresponding $X^0_k$ is stable in the BGI sense, and the corresponding BGI approximate point symmetry of the perturbed wave equation (3.7). The third column contains the same information for the FS approximate point symmetries of the perturbed equation (3.7). Trivial approximate symmetries are not listed (see Section 3.1).

In Table 1, the generators of BGI and FS symmetries were sought in the evolutionary forms, depending on at most first derivatives of the fields. As a result, approximate BGI and FS point symmetries were obtained, as well as some local first-order approximate FS symmetries. We also note that for the point symmetry defined by $\hat{X}^0_5$ to be stable in the FS sense, the arbitrary function $F_1=F_1(v, v_t)$ must satisfy three coupled nonlinear PDEs. Because of the complexity of the latter, a general form of $F_1$ could not be obtained explicitly; however, it was verified that the form of $F_1$ as listed satisfies those equations. The same takes place for the point symmetry defined by $\hat{X}^0_6$ .

Table 1 illustrates differences between BGI and FS frameworks as there are stable symmetries in one framework and unstable in the other framework. Some specific examples are listed below. The following observations can be made.

• • When $F_1(u,u_t)=u_t$ or $F_1={\rm const}$ , all point symmetries (2.12) of the unperturbed PDE (2.9) are stable as BGI and FS symmetries.

• • For all forms $F_1(u,u_t)$ , the PDE (3.7) and the FS system (3.9) are invariant under t- and x-translations. Consequently, the exact symmetries $\hat{X}^0_3$ and $\hat{X}^0_4$ reappear as BGI and FS approximate symmetries without change.

• • As an example, as seen in the first row corresponding to $\hat{X}^0_1$ , for $F_1=uu_t$ , this symmetry is stable as a point BGI approximate symmetry and as a local (but not point) first-order FS approximate symmetry.

• • It can be separately proven that for every form of $F_1$ where a symmetry $\hat{X}^0_i$ is stable as a BGI point approximate symmetry, it is also stable as a local (point or first-order) FS approximate symmetry.

• • The converse is not true; for example, the scaling symmetry $\hat{X}^0_6$ is stable as a local FS approximate symmetry for $F_1= v^{a_3}Q_6(v_t)+a_2v_t + a_1$ , but is only stable as a BGI approximate symmetry when $a_3=-1$ .

We also note that genuine BGI approximate symmetries are given by $\hat{X}_1$ , $\hat{X}_2$ , $\hat{X}_5$ and $\hat{X}_6$ ; similarly, genuine FS approximate symmetries correspond to the local generators $\hat{Z}_1$ , $\hat{Z}_2$ , $\hat{Z}_5$ and $\hat{Z}_6$ .

3.4. Instability of local symmetries of unperturbed PDEs in terms of higher-order approximate symmetries: An example

For an ordinary differential equation (ODE), all local symmetries are stable in the BGI sense: each local symmetry of a given ODE corresponds to an approximate BGI local, often higher-order, symmetry of its perturbed version [40]. For a PDE, in general, this is not the case. For the BGI framework, differential functions ${\zeta_1}_{i_1 i_2 \ldots i_p }^{(p)} $ in the determining equation (2.23) contain derivatives of u of orders higher than those in the differential function $\zeta^1$ . It follows that the left-hand side of equation (2.23) splits into a system of linear PDEs in $\zeta_1$ . On the other side, the function H may contain derivatives of u with respect to other variables different than those in the left-hand side of equation (2.23). This can lead to some constraints on the unperturbed symmetry component $\zeta_0$ ; in that case, an exact local symmetry of the unperturbed PDE (2.1) may not correspond to a local approximate BGI symmetry of the perturbed PDE (2.13). A similar argument holds for FS approximate symmetries. The main reason, as it can be seen in the example below, is the existence of multiple kinds of derivatives in PDEs and thus more restrictive conditions that arise for $\zeta_0$ when the determining equations are being split with respect to higher-order derivatives.

As an illustration, consider the PDE (3.7) with $F_1(u,u_t)=uu_t$ :

(3.11) \begin{equation} u_{tt}+\epsilon uu_t=u_xu_{xx}\end{equation}

and the related Fushchich–Shtelen system (3.9)

(3.12) \begin{equation}v_{tt}-v_xv_{xx} = 0,\qquad w_{tt}+vv_t-v_xw_{xx}-w_xv_{xx}=0.\end{equation}

From Table 1, one can see that $\hat{X}^0_2=t \,\partial/\partial u$ is unstable as a point symmetry in both BGI and FS frameworks; that is, the point symmetry $\hat{X}^0_2$ admitted by the PDE (3.11) with $\epsilon=0$ corresponds to no approximate point or first-order symmetry arising from BGI or FS approaches. First, we examine whether or not it is possible to construct a local, possibly higher-order, BGI approximate symmetry of (3.11) that would correspond to $\hat{X}^0_2$ . The generator of such a symmetry would have the form

(3.13) \begin{equation}\hat{X}_2 =C_2\hat{X}_2^0+\epsilon \hat{X}_2^1=\left( C_2 t +\epsilon \zeta_1[u]\right)\dfrac{\partial}{\partial u},\end{equation}

where $C_2={\rm const}\ne 0$ . The determining equation (2.23) for BGI local symmetries reads

(3.14) \begin{equation} \left(D_t^2\zeta_1-u_xD_x^2\zeta_1-u_{xx}D_x\zeta_1\right)\bigg|_{u_{tt}=u_xu_{xx}}={H},\end{equation}

where one readily finds

(3.15) \begin{equation}H=C_2(tu_t+u).\end{equation}

One can show by a direct computation that whatever the dependence of $\zeta_1$ on partial derivatives of u is chosen to be, higher-order derivatives of u that arise in (3.14) lead to constraints on $C_2$ that result in $C_2=0$ , which means that no nontrivial BGI point symmetry (3.13) corresponding to $\hat{X}^0_2$ exists.

Second, we seek a local, possibly higher-order, approximate FS symmetry of the PDE (3.11) corresponding to $\hat{X}^0_2$ . Such a symmetry would arise as an exact local symmetry of Fushchich–Shtelen system (3.12). The corresponding evolutionary generator (2.29) has the form

(3.16) \begin{equation} Z= \psi_1[v,w] \dfrac{\partial}{\partial v}+\psi_2[v,w] \dfrac{\partial}{\partial w}.\end{equation}

As noted in Section 2.2.3, the determining equation for the first equation of the system (3.12) is satisfied when $\psi_1=C_2t$ as in $\hat{X}^0_2$ . Now the determining equation for the second PDE of (3.12) leads to

(3.17) \begin{equation} \left(D_t^2\psi_2-v_xD_x^2\psi_2-v_{xx}D_x\psi_2\right)\bigg|_{v_{tt}=v_xv_{xx},w_{tt}=-vv_t+v_xw_{xx}+w_xv_{xx} }=C_2(tv_t+v).\end{equation}

It can be shown that for any dependence $\psi_2[v,w]$ , constraints on $C_2$ exist, leading to $C_2=0$ . Consequently, there is no higher-order FS symmetry corresponding to the unstable point symmetry $\hat{X}^0_2$ admitted by the wave equation (3.11) with $\epsilon=0$ .

3.5. A relation between BGI and FS approximate symmetries

The computational example of Section 3.3 above illustrated the fact that BGI and FS frameworks can yield rather different approximate point symmetry classifications for the same PDE with a small parameter. However, in certain situations, the two approaches can lead to related results. We now show that for a specific class of (1+1)-dimensional PDEs, a stable BGI approximate point symmetry always correspond to a stable FS approximate local symmetry.

Consider the following class of PDEs on u(x,t), written in the Kovalevskaya form with respect to an independent variable t:

(3.18) \begin{equation} \dfrac{\partial ^n u}{\partial t^n} =F_0[u],\qquad F_0[u]\equiv F_0(x,t,u,\partial u, \partial^2 u,...,\partial^k u),\end{equation}

and its perturbed version with a small parameter $\epsilon$ :

(3.19) \begin{equation} \dfrac{\partial ^n u}{\partial t^n}=F_0[u]+\epsilon F_1[u],\qquad F_1[u]\equiv F_1(x,t,u,\partial u, \partial^2 u,...,\partial^\ell u).\end{equation}

A local BGI approximate symmetry of a PDE (3.19) has the form (2.22)

(3.20) \begin{equation}\hat{X} = \hat{X}^0+\epsilon \hat{X}^1 = \left(\zeta_0[u]+\epsilon \zeta_1[u]\right)\dfrac{\partial}{\partial u}\, .\end{equation}

As per Theorem 2.1 and Remark 2.2, the O(1) term in (3.20) corresponds to a local symmetry

(3.21) \begin{equation}\hat{X}^0 = \zeta_0[u]\,\dfrac{\partial}{\partial u}\end{equation}

of the unperturbed equation (3.18).

In order to compute FS approximate symmetries of a PDE (3.19), we substitute $u(x,t)=v(x,t)+\epsilon w(x,t)+o(\epsilon)$ into (3.19) and split the orders of $\epsilon$ to get the Fushchich–Shtelen system

(3.22) \begin{align} \dfrac{\partial ^n v}{\partial t^n} & =F_0[v],\nonumber \\[4pt]\dfrac{\partial ^n w}{\partial t^n} & ={F_0}_vw+{F_0}_{v_i}w_i+{F_0}_{v_{ij}}w_{ij}+...+{F_0}_{v_{i_1i_2...i_k}}w_{i_1i_2...i_k}+F_1[v]. \end{align}

The evolutionary generator of a local FS approximate symmetry has the form (2.29). The determining equations (2.4) for exact local symmetries of (3.22) are

(3.23a) \begin{align} &\hat{Z}^{(n)}\left(\dfrac{\partial ^n v}{\partial t^n}- F_0\right)=0,\\[-10pt]\nonumber \end{align}

(3.23b) \begin{align} &\hat{Z}^{(n)}\left(\dfrac{\partial ^n w}{\partial t^n}-{F_0}_vw-{F_0}_{v_i}w_i-{F_0}_{v_{ij}}w_{ij}-...-{F_0}_{v_{i_1i_2...i_k}}w_{i_1i_2...i_k}-F_1\right)=0,\end{align}

holding on solutions of (3.22).

Theorem 3.1. If (3.20) is a BGI approximate local symmetry generator of a PDE (3.19) having the specific form

(3.24) \begin{equation}\hat{X}=\left(\zeta_0(x,t)+\epsilon \zeta_1(x,t,u,u_x,u_t)\right)\dfrac{\partial}{\partial u}\end{equation}

and additionally, $F_0[v]$ in (3.19) satisfies the following system of equations

(3.25) \begin{align} & \zeta^0{F_0}_{uu}+\zeta_i^{0^{(1)}}{F_0}_{uu_i}+\zeta_{i_1i_2}^{0^{(2)}}{F_0}_{uu_{i_1i_2}}+...+\zeta_{i_1i_2...i_k}^{0^{(k)}}{F_0}_{uu_{i_1i_2...i_k}}=0,\nonumber\\[3pt] & \zeta^0{F_0}_{uu_i}+\zeta_i^{0^{(1)}}{F_0}_{u_iu_i}+\zeta_{i_1i_2}^{0^{(2)}}{F_0}_{u_iu_{i_1i_2}}+...+\zeta_{i_1i_2...i_k}^{0^{(k)}}{F_0}_{u_iu_{i_1i_2...i_k}}=0,\nonumber\\[3pt]& \vdots \\[3pt]& \zeta^0{F_0}_{uu_{i_1i_2...i_k}}+\zeta_i^{0^{(1)}}{F_0}_{u_iu_{i_1i_2...i_k}}+\zeta_{i_1i_2}^{0^{(2)}}{F_0}_{u_{i_1i_2}u_{i_1i_2...i_k}}\nonumber\\[3pt] & \qquad \qquad +...+\zeta_{i_1i_2...i_k}^{0^{(k)}}{F_0}_{u_{i_1i_2...i_k}u_{i_1i_2...i_k}}=0,\nonumber\end{align}

then

(3.26) \begin{equation}\hat{Z}=\zeta^0(x,t)\dfrac{\partial}{\partial v}+\zeta^1(x,t,v,v_x,v_t) \dfrac{\partial}{\partial w}\end{equation}

is a FS approximate local symmetry of the perturbed PDE (3.19) to the point symmetry generator $\hat{X}^{0}=\zeta^0\,\partial/\partial v$ of the unperturbed PDE (3.18).

Proof. We need to show that under the stated conditions, the determining equations (2.20) for BGI approximate symmetries of (3.19) are equivalent to the determining equations for FS approximate symmetries. Since the first PDE of the Fushchich–Shtelen system (3.22) is the same as the unperturbed equation (3.18), the first FS determining equation (3.23a) is satisfied for any $\zeta^0$ and $\zeta^1$ as long as $\zeta^0$ is an exact point symmetry component of (3.18).

The second FS determining equation (3.23b) with $\psi_1=\zeta^0$ , $\psi_2=\zeta^1$ can be rewritten as

(3.27) \begin{equation} \left(\zeta_t^{1^{(n)}}-\zeta^1{F_0}_{v}-\zeta_i^{1^{(1)}}{F_0}_{v_i}-\zeta_{i_1i_2}^{1^{(2)}}{F_0}_{v_{i_1i_2}}-...-\zeta_{i_1i_2...i_k}^{1^{(k)}}{F_0}_{v_{i_1i_2...i_k}}\right)\bigg|_{{\partial ^n v}/{\partial t^n}=F_0}=G,\end{equation}

where

\begin{align*} G&=w\left(\zeta^0{F_0}_{vv}+\zeta_i^{0^{(1)}}{F_0}_{vv_i}+\zeta_{i_1i_2}^{0^{(2)}}{F_0}_{vv_{i_1i_2}}+...+\zeta_{i_1i_2...i_k}^{0^{(k)}}{F_0}_{vv_{i_1i_2...i_k}}\right)\nonumber\\[4pt] &\quad +w_i\left(\zeta^0{F_0}_{vv_i}+\zeta_i^{0^{(1)}}{F_0}_{v_iv_i}+\zeta_{i_1i_2}^{0^{(2)}}{F_0}_{v_iv_{i_1i_2}}+...+\zeta_{i_1i_2...i_k}^{0^{(k)}}{F_0}_{v_iv_{i_1i_2...i_k}}\right) +...\nonumber\\[4pt] &\quad +w_{i_1i_2...i_k}\left(\zeta^0{F_0}_{vv_{i_1i_2...i_k}}+\zeta_i^{0^{(1)}}{F_0}_{v_iv_{i_1i_2...i_k}}+\zeta_{i_1i_2}^{0^{(2)}}{F_0}_{v_{i_1i_2}v_{i_1i_2...i_k}}+... \right.\\[4pt]& \quad \left. +\,\zeta_{i_1i_2...i_k}^{0^{(k)}}{F_0}_{v_{i_1i_2...i_k}v_{i_1i_2...i_k}}\right) + \zeta^0{F_1}_v+\zeta_i^{0^{(1)}}{F_1}_{v_i}+...+\zeta_{i_1i_2...i_\ell}^{0^{(\ell)}}{F_1}_{v_{i_1i_2...i_\ell}}.\nonumber\end{align*}

As $\zeta^0$ and $F_0[v]$ satisfy (3.25), G reduces to

(3.28) \begin{equation} G= \zeta^0{F_1}_v+\zeta_i^{0^{(1)}}{F_1}_{v_i}+...+\zeta_{i_1i_2...i_\ell}^{0^{(\ell)}}{F_1}_{v_{i_1i_2...i_\ell}}.\end{equation}

Now, we proceed to check the determining equation (2.20) of BGI approximate symmetries for (3.19). The left-hand side of (2.20) simplifies to

\begin{equation*} \left(\zeta_t^{1^{(n)}}-\zeta^1{F_0}_{u}-\zeta_i^{1^{(1)}}{F_0}_{u_i}-\zeta_{i_1i_2}^{1^{(2)}}{F_0}_{u_{i_1i_2}}-...-\zeta_{i_1i_2...i_k}^{1^{(k)}}{F_0}_{u_{i_1i_2...i_k}}\right)\bigg|_{{\partial ^n u}/{\partial t^n}=F_0}\end{equation*}

which is equivalent to the left-hand side of (3.27). Now, the right-hand side of (2.20), the function H, is the coefficient of $\epsilon$ in

(3.29) \begin{equation} -\hat{X}^{0^{(n)}}\left(\dfrac{\partial ^n u}{\partial t^n}-F_0-\epsilon F_1\right)\bigg|_{{\partial ^n u}/{\partial t^n}=F_0+\epsilon F_1}.\end{equation}

Since $\zeta^0=\zeta^0(x,t)$ , none of the terms in (3.29) contains $\partial^n u/\partial t^n$ . Hence, the coefficient of $\epsilon$ in (3.29) is

(3.30) \begin{equation} H=\hat{X}^{0^{(n)}} F_1= \zeta^0{F_1}_u+\zeta_i^{0^{(1)}}{F_1}_{u_i}+...+\zeta_{i_1i_2...i_\ell}^{0^{(\ell)}}{F_0}_{u_{i_1i_2...i_\ell}}.\end{equation}

The latter is equivalent to G (3.28). It follows that the determining equation (3.5) of FS symmetries for the system (3.22) and the determining equation (2.20) of BGI approximate symmetries for the PDE (3.19) are equivalent. Hence, $\hat{Z}$ (3.26) is a FS approximate local symmetry of the system (3.22).

The above theorem states that when a point symmetry of an unperturbed PDE yields a BGI approximate point symmetry but not a FS approximate point symmetry of the perturbed PDE, under the conditions of the theorem, there exists a corresponding higher-order FS approximate symmetry of the perturbed PDE instead. This is illustrated by the following example.

Example 3.1 Consider again the PDE (3.11) $u_{tt}+\epsilon uu_t=u_xu_{xx}$ . Using Table 1, we observe that $X^0_1=\partial/\partial u$ is unstable as a FS approximate point symmetry but it is a stable point symmetry in the sense of BGI; the corresponding BGI generator is given by

\begin{equation*}\hat{X}_1=\left(1-\dfrac{\epsilon}{10}t(tu_t+4u)\right)\dfrac{\partial}{\partial u},\end{equation*}

with H (3.30) given by (3.15) with $C_1=1$ . Consider now the Fushchich–Shtelen system (3.12) for the PDE (3.11). Using determining equation (3.5) for exact local symmetries of (3.12), one can find that

\begin{equation*}\hat{Z}_1=\dfrac{\partial}{\partial v}-\dfrac{1}{10}t(tv_t+4v)\dfrac{\partial}{\partial w}\end{equation*}

is a higher-order FS approximate symmetry generator of the PDE (3.11).

Remark 3.1. The conditions of Theorem 3.1 are not satisfied when $\zeta^0_u\neq 0$ , $\zeta^0_{u_x}\neq 0$ or $\zeta^0_{u_t}\neq 0.$

Example 3.2. The perturbed wave equation

(3.31) \begin{equation} u_{tt}+\epsilon uu_t=e^{u}u_{xx},\end{equation}

admits an approximate point symmetry with the evolutionary form

(3.32) \begin{equation} \hat{X}=\left(2-xu_x-\epsilon t\left(tu_t+4\right)\right)\dfrac{\partial}{\partial u}\end{equation}

corresponding to the stable point symmetry $\hat{X}^0=\zeta^0\,{\partial}/{\partial u}=(2-xu_x)\,{\partial}/{\partial u}$ . Here, $\zeta^0=2-xu_x$ does not satisfy the conditions of Theorem 3.1 since it involves $u_x$ . It turns out that $\hat{X}^0$ is unstable as a FS approximate point symmetry of (3.31). Indeed, it is easy to check that

(3.33) \begin{equation}\hat{Z}=\left(2-xv_x\right)\dfrac{\partial}{\partial v}-t\left(tv_t+4\right)\dfrac{\partial}{\partial w}\end{equation}

is not a local symmetry of the Fushchich–Shtelen system of the PDE (3.31) given by

(3.34) \begin{equation} v_{tt}-v_xv_{xx} =0, \qquad w_{tt}+vv_t-e^vw_{xx}-e^vwv_{xx} =0.\end{equation}

However, a direct computation shows that a slightly modified version of (3.33) given by

\begin{equation*}\hat{Z}=\left(2-xv_x\right)\dfrac{\partial}{\partial v}-\left(t(tv_t+4)+xw_x\right)\dfrac{\partial}{\partial w}\end{equation*}

yields the corresponding first-order FS approximate symmetry.

4. Exact and approximate point symmetry classification of a one-dimensional perturbed wave model in a fiber-reinforced solid

One-dimensional nonlinear wave equations

(4.1) \begin{equation}u_{tt}=K(u_x)u_{xx}\end{equation}

on the unknown u(x,t) and various forms of $K(u_x)$ arise in multiple physical contexts, in particular, in nonlinear mechanics [32]. The point symmetry classification of the PDE family (4.1) has been performed by Oron and Rosenau [34]. If $K(u_x)=c^2={\rm const}$ , the PDE (4.1) becomes linear:

(4.2) \begin{equation}u_{tt}=c^2 u_{xx},\end{equation}

and consequently, its Lie symmetry algebra is infinite-dimensional and consists of the vector fields of the form

(4.3) \begin{equation} X^0_\infty =(\alpha_1+\alpha_2)\dfrac{\partial}{\partial t}+c(\alpha_1-\alpha_2)\dfrac{\partial}{\partial x}+(C_1u+\beta_1+\beta_2)\dfrac{\partial}{\partial u},\end{equation}

parameterized by an arbitrary constant $C_1$ and four arbitrary functions $\alpha_1(x+ct)$ , $\beta_1(x+ct)$ , $\alpha_2(x-ct)$ , and $\beta_2(x-ct)$ .

In the current section, we consider a special form of the arbitrary function $K(u_x)=c^2+\epsilon Q(u_x)$ in (4.1), which yields a PDE family

(4.4) \begin{equation}u_{tt}=(c^2+\epsilon Q(u_x))u_{xx}\end{equation}

with a small parameter $\epsilon$ . It is assumed that $Q(u_x)\ne {\rm const}$ . Such models arise, for example, in the analysis of wave propagation in fiber-reinforced elastic solids [3, 13] with small fiber strengths. The PDEs (4.4) are nonlinear perturbed versions of the linear PDE (4.2) and therefore have a reduced set of symmetries compared to that of the linear wave equation. It is of interest to follow the algorithms presented in Section 2 to compare the exact point symmetry classification of the PDE family (4.4) as it stands with approximate (BGI and FS) point symmetries of the PDEs (4.4) viewed as perturbations of the linear wave equation (4.2).

We classify exact and approximate (BGI and FS) point symmetries for (4.4). The classification is performed with respect to the forms of the arbitrary function $Q(u_x)$ , with each classification case holding for an arbitrary $\epsilon$ . In the classifications, cases are simplified using the equivalence transformations of the perturbed equation (4.4), given by

(4.5) \begin{align}t & =C_1\tilde{t}+C_2,\quad x=C_3\tilde{x}+C_4,\quad u=C_5\tilde{u}+C_6\tilde{x}+C_7\tilde{t}+C_8,\nonumber\\c^2 & = \dfrac{C_3^2}{C_1^2}\tilde{c}^2,\quad Q(u_x) = \dfrac{C_3^2}{C_1^2}\tilde{Q}(\tilde{u}_{\tilde{x}}),\end{align}

involving arbitrary constants $C_i$ , $i=1,\ldots, 8$ . It follows that by taking $C_1=1/c$ , $C_3=C_5=1$ , and other constants zero, upon dropping tildes, one obtains the PDE (4.4) with $c^2=1$ :

(4.6) \begin{equation}u_{tt}=(1+\epsilon Q(u_x))u_{xx},\end{equation}

which will be considered below.

The results below are presented modulo the equivalence transformations (4.5), usually without obvious trivial approximate symmetries (see Section 3.1); some trivial approximate symmetries will be pointed out.

4.1. Exact point symmetries of the wave equation family (4.6)

The exact symmetry generator for the PDE (4.6) has the form

(4.7) \begin{equation} Y=\xi^1(x,t,u;\,\epsilon)\dfrac{\partial}{\partial t}+\xi^2(x,t,u;\,\epsilon)\dfrac{\partial}{\partial x}+\eta(x,t,u;\,\epsilon)\dfrac{\partial}{\partial u}.\end{equation}

The following cases arise, holding for an arbitrary $\epsilon$ and non-constant Q.

1. 1. In the general case of arbitrary $Q(u_x)$ , one has the five-dimensional Lie group of point symmetries generated by

   (4.8) \begin{equation}Y_1=\dfrac{\partial}{\partial t},\quad Y_2=\dfrac{\partial}{\partial x},\quad Y_3=\dfrac{\partial}{\partial u},\quad Y_4 = t\dfrac{\partial}{\partial u},\quad Y_5=t\dfrac{\partial}{\partial t}+x\dfrac{\partial}{\partial x}+u\dfrac{\partial}{\partial u},\end{equation}
   describing respectively translations in t, x, u, the Galilei transformation in the direction of the displacement u, and a homogeneous space-time scaling.

2. 2. In the case when $Q(u_x)=u_x$ , the Lie algebra (4.8) is extended by a point symmetry generator

   (4.9) \begin{equation} Y_6=x\dfrac{\partial}{\partial u}-\epsilon\left(x\dfrac{\partial}{\partial x}+\dfrac{3t}{2}\dfrac{\partial}{\partial t}\right). \end{equation}

4.2. BGI approximate point symmetry classification of the wave equation family (4.6)

The BGI approximate point symmetry generator for the PDE (4.6) has the form

(4.10) \begin{equation}X= X^0+\epsilon X^1=X^0 + \epsilon \left(\xi_1^1(x,t,u)\dfrac{\partial}{\partial t}+\xi_1^2(x,t,u)\dfrac{\partial}{\partial x}+\eta_1(x,t,u)\dfrac{\partial}{\partial u}\right),\end{equation}

where, according to Theorem 2.1, the freedom in $X^0$ does not exceed that in $X^0_\infty$ (4.3). From the determining equations (2.19) for BGI point symmetries, the following cases arise.

1. 1. $Q(u_x)$ arbitrary: the O(1) and $O(\epsilon)$ components of the generator (4.10) are given by

   (4.11a) \begin{align} X^0=&(C_1t+C_2)\dfrac{\partial}{\partial t}+(C_1x+C_3)\dfrac{\partial}{\partial x}+(C_1u+C_4t+C_5)\dfrac{\partial}{\partial u},\\[-10pt]\nonumber \end{align}
   (4.11b) \begin{align} X^1=&(\lambda_1+\lambda_2)\dfrac{\partial}{\partial t}+ \left(\lambda_1-\lambda_2\right)\dfrac{\partial}{\partial x}+(C_6u+\lambda_3+\lambda_4)\dfrac{\partial}{\partial u}, \end{align}

   where $C_i={\rm const}$ , $\lambda_1$ and $\lambda_3$ are arbitrary functions of $x+t$ , and $\lambda_2$ and $\lambda_4$ are arbitrary functions of $x-t$ . Consequently, the nonlinear wave equation (4.6) for an arbitrary $Q(u_x)$ admits the approximate point symmetries

   (4.12) \begin{align} X_1 & =t\dfrac{\partial}{\partial t}+x\dfrac{\partial}{\partial x}+u\dfrac{\partial}{\partial u},\quad X_2=\dfrac{\partial}{\partial t},\quad X_3=\dfrac{\partial}{\partial x},\quad X_4=t\dfrac{\partial}{\partial u},\quad X_5 = \dfrac{\partial}{\partial u},\nonumber\\[5pt]X_\infty & =\epsilon\left[\left(\lambda_1+\lambda_2\right)\dfrac{\partial}{\partial t}+\left(\lambda_1-\lambda_2\right)\dfrac{\partial}{\partial x}+\left(C_6u+\lambda_3+\lambda_4\right)\dfrac{\partial}{\partial u}\right],\end{align}
   which are, respectively, the re-numbered exact point symmetries (4.8), and a trivial approximate symmetry $X_\infty$ corresponding to the infinite symmetry set (4.3) of the linear wave equation (4.2). The difference between the freedom in (4.3) and the space spanned by $X_1,\ldots,X_5$ in (4.12) corresponds to unstable point symmetries of the linear wave equation.

2. 2. $Q(u_x)=A\ln(u_x+B)+C$ , where A, B and C are arbitrary constants: here, the nonlinear wave equation (4.6) admits an additional genuine approximate symmetry given by

   (4.13) \begin{equation} X_g=(u+Bx)\dfrac{\partial}{\partial u}-\epsilon \, \dfrac{At}{2}\dfrac{\partial}{\partial t}.\end{equation}

3. 3. $Q(u_x)=u_x$ : the exact symmetry generator (4.3) of the linear wave equation (4.2) reduces to

   (4.14) \begin{equation}X^0=(C_1t+C_2)\dfrac{\partial}{\partial t}+(C_1x+C_3)\dfrac{\partial}{\partial x}+(C_1u+\beta_1+\beta_2)\dfrac{\partial}{\partial u},\end{equation}
   and the $O(\epsilon)$ approximate symmetry components have the form
   (4.15) \begin{align}\xi_1^1 & = \lambda_1(x+t)+\lambda_2(x-t)-\dfrac{1}{4}\displaystyle\int_{}^{t}\left(\beta_2'(x-t)+2z\beta_1''\left(\left(t-2z\right)+x\right)\right)dz, \nonumber\\[5pt]\xi_1^2 & =H(x,t), \qquad \eta_1 =\left( C_4+\dfrac{\beta_2'(x-t)-\beta_1'(x+t)}{4}\right)u+\lambda_3(x+t)+\lambda_4(x-t),\end{align}
   where H(x,t) is an arbitrary solution of the PDEs: $H_t={\xi^1_1}_x$ , $H_x={\xi^1_1}_t+\dfrac{1}{2}\left({\beta_1}_x+{\beta_2}_x\right)$ . In this third case, the point symmetries of the unperturbed linear wave equation (4.2) with arbitrary $\beta_1(x+t)$ and $\beta_2(x-t)$ remain stable and yield genuine approximate symmetries with $O(\epsilon)$ components given by the terms in (4.15) that contain $\beta_1$ and $\beta_2$ . In particular, a genuine approximate symmetry similar to (4.13) is admitted:
   (4.16) \begin{equation} Y_g=x\dfrac{\partial}{\partial u}-\epsilon \, \dfrac{t}{2}\dfrac{\partial}{\partial t}.\end{equation}

4.3. FS approximate point symmetry classification of the wave equation family (4.6)

For the perturbed PDE (4.6) with the solution form

(4.17) \begin{equation}u(x,t)=v(x,t)+\epsilon w(x,t)+o(\epsilon),\end{equation}

the Fuschich-Shtelen system (2.25) reads

(4.18) \begin{equation} v_{tt}=v_{xx},\qquad w_{tt} =w_{xx}+ Q(v_x)v_{xx}.\end{equation}

The basic group of point equivalence transformations of the PDE system (4.18) is given by

(4.19) \begin{align}t & =C_1\tilde{t}+C_2,\qquad\qquad\qquad\ \ x=C_3\tilde{x}+C_4,\nonumber\\v & =C_5\tilde{v}+C_6\tilde{x}+C_7\tilde{t}+C_8,\quad w=C_5\tilde{w}+C_9\tilde{x}+C_{10}\tilde{t}+C_{11},\\c^2 & = \dfrac{C_3^2}{C_1^2},\qquad\qquad\qquad\qquad\ \ Q(u_x) = \dfrac{C_3^2}{C_1^2}\,\tilde{Q}(\tilde{u}_{\tilde{x}}),\nonumber\end{align}

in terms of arbitrary constants $C_1, \ldots, C_{11}$ . We now find exact point symmetries of the system (4.18) that correspond to FS approximate point symmetries of the PDE (4.6). The infinitesimal generator of such symmetries has the form

(4.20) \begin{equation} Z= \lambda^1(x,t,v,w)\dfrac{\partial}{\partial x}+ \lambda^2(x,t,v,w)\dfrac{\partial}{\partial t}+\phi_1(x,t,v,w) \dfrac{\partial}{\partial v}+ \phi_2(x,t,v,w) \dfrac{\partial}{\partial w}. \end{equation}

The solution of the determining equations (2.4) leads to the following classification.

1. 1. $Q(v_x)$ arbitrary:

   (4.21) \begin{align} Z_1 & =\dfrac{\partial}{\partial t},\qquad Z_2=\dfrac{\partial}{\partial x},\qquad Z_3=\dfrac{\partial}{\partial v},\qquad Z_4=t\dfrac{\partial}{\partial v},\nonumber\\Z_5 & = t\dfrac{\partial}{\partial t}+x\dfrac{\partial}{\partial x}+v\dfrac{\partial}{\partial v}+w\dfrac{\partial}{\partial w},\qquad Z_6 = v\dfrac{\partial}{\partial w},\\Z_\infty & =\left(\beta_1(x+t)+\beta_2(x-t)\right)\dfrac{\partial}{\partial w}, \nonumber \end{align}
   where $\beta_1=\beta_1(x+t)$ and $\beta_2=\beta_2(x-t)$ are arbitrary functions. In this general case, no genuine FS approximate symmetries arise. Indeed, the generators $Z_1,\ldots,Z_5$ mimic the exact point symmetry generators (4.8), and $Z_6, Z_\infty$ are trivial FS symmetries. Including the above general symmetries, and modulo the equivalence transformations (4.19), the system (4.18) admits additional point symmetries in the following cases.

2. 2. $Q(v_x)=v_x^s$ , $s\neq 0$ :

   (4.22) \begin{equation} Z_7=v\dfrac{\partial}{\partial v}+(s+1)w\dfrac{\partial}{\partial w}. \end{equation}

3. 3. $Q(v_x)=e^{v_x}:$

   (4.23) \begin{equation} Z_7'=x\dfrac{\partial}{\partial v}+w\dfrac{\partial}{\partial w}. \end{equation}

The symmetries given by $Z_7$ and $Z_7'$ are genuine FS approximate point symmetries of the perturbed PDE (4.6). In particular, $Z_7$ corresponds to a scaling which is different for O(1) and $O(\epsilon)$ components of the approximate solution (4.17).

4.4. Summary

For an arbitrary Q, the perturbed one-dimensional wave equation (4.6) admits five exact symmetries given by (4.8) and it has these five symmetries and a trivial approximate symmetry as BGI approximate symmetries (4.12). For $Q=u_x$ , the equation (4.6) admits (4.8) and an additional exact symmetry given by (4.9). For BGI classification with $Q(u_x)=u_x$ , the PDE (4.6) has an infinite set of BGI approximate symmetries with approximate symmetry components given by (4.15). Note that the exact symmetry generator $Y_6$ in (4.9)

\begin{equation*}Y_6=x\dfrac{\partial}{\partial u}-\epsilon\left(\dfrac{3t}{2}\dfrac{\partial}{\partial t}+x\dfrac{\partial}{\partial x}\right)\end{equation*}

can be obtained from the BGI approximate components (4.15) by taking

\begin{equation*}\beta_1=\dfrac{x+t}{2},\quad \beta_2=\dfrac{x-t}{2},\quad \lambda_1=-\dfrac{11}{16}(x+t),\quad \lambda_2=\dfrac{11}{16}(x-t).\end{equation*}

It follows that the BGI approximate symmetry classification of the wave equation (4.6) includes the exact symmetry classification of (4.6) but corresponds to a subset of exact point symmetries (4.3) of the unperturbed (linear) wave equation (4.2). An additional case appears in the BGI approximate symmetry classification when $Q=A\ln(u_x+B)+C$ , with a corresponding additional approximate symmetry given by (4.13). This case does not arise in the FS symmetry classification.

For an arbitrary Q, the PDE (4.6) admits exact point symmetry generators (4.8) and trivial approximate symmetries given by (4.21). In comparison with the exact and BGI symmetry classifications of (4.6), two different cases appear in FS approximate symmetry classification: $Q=v_x^s$ , $s\neq 0$ and $Q=e^{v_x}$ . For $Q=v_x^s$ , the PDE (4.6) admits an additional FS approximate symmetry given by (4.22). (By contrast, in the exact and BGI symmetry classifications of (4.6), this case appears only when $s=1$ .) For $Q=e^{v_x}$ , a stable exact symmetry $x \,\partial/\partial u$ of the linear wave equation (4.2) yields a genuine FS approximate symmetry of (4.6) given by (4.23).

5. Approximate and numerical solutions modeling breaking waves in fiber-reinforced materials

The displacements in shear waves propagating in an incompressible hyperelastic material with a single family of fibers directed along the wave propagation are governed by a nonlinear one dimensional wave equation

(5.1) \begin{equation} u_{tt}=(\alpha+3\beta u_x^2)u_{xx},\quad u=u(x,t),\end{equation}

where the constants $\alpha, \beta >0$ are the material parameters [13]. In this section, we consider wave equations (4.4) with $Q(u_x)=Bu_x^s$ , $B>0$ , $s\neq 0$ , which include the model (5.1). By a re-scaling of x, t and u, these PDEs can be brought into a simpler form

(5.2) \begin{equation} u_{tt}=(1+\epsilon u_x^s)u_{xx}. \end{equation}

5.1. Approximate solutions of the PDE (5.2) in the FS framework

Here we use Fuschich-Shtelen approximate symmetries to construct an approximate solution for the PDE (5.2) in the usual FS form

(5.3) \begin{equation}u(x,t)=v(x,t)+\epsilon w(x,t)+o(\epsilon).\end{equation}

In the first-order of precision in $\epsilon$ , the equation (5.2) is equivalent to the Fuschich-Shtelen system (4.18) with $Q(v_x)=v_x^s$ :

(5.4) \begin{equation} v_{tt}=v_{xx},\qquad w_{tt} =w_{xx}+ v_x^sv_{xx},\end{equation}

which admits the symmetry generator (4.22). The corresponding characteristic equations are given by

(5.5) \begin{equation}\dfrac{dt}{0} =\dfrac{dx}{0} =\dfrac{dv}{v}=\dfrac{dw}{(s+1)w}.\end{equation}

Consequently, if v(x,t) is any solution for the first equation of the system (5.4), then the invariant solution following from the characteristic equations (5.5) is given by $w(x,t)=v^{s+1}\phi(x,t)$ . Consider traveling wave solutions of the first equation in (5.4):

(5.6) \begin{equation} v=g\left( x \pm t\right).\end{equation}

Substituting (5.6) and $w= g^{s+1}\phi$ into the second PDE of (5.4), one gets to the PDE in $\phi$

(5.7) \begin{equation} g^{s+1}\left(\phi_{tt}-\phi_{xx}\right)+2(s+1)g^sg'\left(\pm\phi_{t}-\phi_{x}\right)-(g')^{s}g'' =0.\end{equation}

When $s\ne -1$ , the PDE (5.7) has a general solution

\begin{equation*} \phi=g^{-s-1}\left[h \pm \dfrac{t(g')^{s+1}- \int_{}^{t}(g')^{s+1}\left(\pm (2r-t)+x\right)dr}{2(s+1)}\right],\end{equation*}

where $h=h(x,t)$ satisfies $h_{tt}=h_{xx}$ . Similarly, when $s=-1$ , the solution form changes to

\begin{equation*} \phi=h\pm \dfrac{1}{2}t\ln\left(|g'|\right)\,.\end{equation*}

In light of the above results, the higher-order solution part w has the form

(5.8) \begin{equation} w= \begin{cases} h\pm \dfrac{t(g')^{s+1}- \int_{}^{t}(g')^{s+1}\left(\pm (2r-t)+x\right)dr}{2(s+1)}, & s\neq -1, \\[2ex] h\pm \dfrac{t\ln\left(|g'|\right)}{2}, & s=-1. \end{cases}\end{equation}

Finally, when $s\neq -1$ , the perturbed equation (5.2) has the approximate solution (5.3) given by

(5.9a) \begin{align} u(x,t) = g(x\pm t)+\epsilon\left[ h(x,t)\pm \dfrac{ t(g')^{s+1}- \int_{}^{t}(g')^{s+1}\left(\pm (2r-t)+x\right)dr}{2(s+1)}\right]+o(\epsilon)\,.\end{align}

When $s=-1$ , the approximate solution takes the form

(5.9b) \begin{equation} u(x,t)=g(x\pm t)+\epsilon \left(h(x,t)\pm \dfrac{t\ln\left(|g'\left( x\pm t\right)\right|)}{2}\right)+o(\epsilon).\end{equation}

Example 5.1. As a specific example that will be used below, we consider the PDE (5.1) describing shear waves in a fiber-reinforced solid, re-scaled to the form (5.2) with $s=2$ :

(5.10) \begin{equation} u_{tt} = \left(1 + \epsilon \left(u_x\right)^2\right) u_{xx},\end{equation}

and assume in (5.1) that $\beta/\alpha\sim\epsilon\ll 1$ , which corresponds to weak fiber effects. We also choose $v(x,t)=\exp({-(x-t)^2})$ . Then, the solution (5.9a) of the PDE (5.2) with $h=0$ reduces to

(5.11) \begin{align} u(x,t)&=e^{-(x-t)^2}+\dfrac{\epsilon}{6} \biggl[8t\left( x-t\right) ^{3}e^{-3\left( x-t\right) ^{2}}+\dfrac{1}{9}\bigg(\left( 12tx-6t^{2}-6x^{2}-2\right) e^{-3\left( x-t\right) ^{2}}\nonumber\\[5pt] &+\left( 12tx+6t^{2}+6x^{2}+2\right) e^{-3\left( x+t\right)^{2}}\bigg)\biggl] + o(\epsilon). \end{align}

Note that for any fixed t, the approximate solution (5.11) approaches zero as $x \rightarrow \infty$ . Also, for any $x \in ({-}\infty,\infty)$ , the solution (5.11) is bounded as $t \rightarrow \infty$ . The solution (5.11) is not a purely right-traveling wave solution but describes an evolving wave form (see Section 5.2 below). In particular, the PDE (5.10) is known to have breaking wave-type solutions [13].

5.2. Numerical simulations of (5.10) and breaking waves

We now compute numerical solutions of the wave equation (5.10) in Example 5.1 in order to model its breaking wave behavior (see Section 5.3 below for details) and provide a reference for comparison with the approximate solutions developed in Section 5.1. Gaussian-type initial conditions corresponding to a right-traveling wave and periodic boundary conditions

(5.12) \begin{equation} u(x,0)=e^{-x^2},\quad u_t(x,0)=2xe^{-x^2},\quad u({-}L,t)=u(L,t),\quad L>0,\quad t \in [0,T],\end{equation}

are posed in the space-time domain $x\in[{-}L,L]$ , $t\geq 0$ , and the equation (5.10) is solved using an explicit finite-difference cross-stencil scheme with constant spatial and temporal steps $\tilde{h}$ , $\tilde{\tau}$ . Following [12], we use a conservative finite-difference scheme developed for the PDE (5.10):

(5.13) \begin{align} U_{t\check{t}} & - U_{x\bar{x}} - \epsilon\,\dfrac{U_x^3 - U_{\bar{x}}^3}{\tilde{h}} =0, \quad\tilde{h} = \frac{2L}{M},\quad \tilde{\tau} = \frac{T}{N}, \nonumber\\[3pt]x_m & = -L + m\tilde{h},\qquad\qquad\qquad\! m = 0,\ldots,M, \\[3pt]t_n & = 0 + n\tilde{\tau}, \qquad\qquad\quad\qquad\! n = 0,\ldots,N,\nonumber\end{align}

with $U=U^n_m$ approximating the value of u(x, t) at the mesh node $(x_m, t_n)$ . Here, $U_{t\check{t}}$ and $U_{x\bar{x}}$ represent the second-order central differences, $U_x$ the first-order forward difference and $U_{\bar{x}}$ the first-order backward difference:

(5.14) \begin{align} U_{t\check{t}} & = \dfrac{U^{n+1}_m - 2U^n_m + U^{n-1}_m}{\tilde{\tau}^2}, U_{x\bar{x}} = \dfrac{U^n_{m+1} - 2U^n_m + U^n_{m-1}}{\tilde{h}^2},\nonumber\\U_x & = \dfrac{U^n_{m+1} - U^n_m}{\tilde{h}}, U_{\bar{x}} = \dfrac{U^n_m - U^n_{m-1}}{\tilde{h}}.\end{align}

The numerical solutions provide a good agreement with approximate solutions (5.11), for a broad range of values of $\epsilon$ , from the initial dimensionless time to the time close to the wave breaking. The breaking time increases approximately as $\epsilon^{-1}$ as $\epsilon$ decreases (see Section 5.3 below). Here we present sample computation results for a relatively large value of the small parameter, $\epsilon=0.5$ . The computation is performed from $t=0$ to $t=4$ close to the wave breaking time. A comparison of the numerical solution and the approximate solution (5.11) of PDE (5.10) with initial and boundary conditions (5.12) at several time snapshots is presented in Figure 1(a). The relative difference at the time step $t_n$ between the approximate and numerical solutions is calculated using 2-norms according to the formula

(5.15) \begin{equation} E_n = E(t_n) = \frac{||u_{approx} - u_{num}||_2}{||u_{approx}||_2}\end{equation}

and is shown in Figure 1(b).

Figure 1. (a) Numerical and approximate profiles of u according to the PDE (5.10) ( $\epsilon = 0.5$ ) with initial conditions (5.12) for $L=10$ , $h=0.01$ , $\tau =h/8$ , and $t = 0, 2, 4, 6, 8$ . (b) Relative difference (5.15) between numerical and approximate solutions.

5.3. Estimates of wave breaking times using approximate and numerical solutions

As discussed in [13], the nonlinearity $(1 + \epsilon \left(u_x\right)^2)$ leads to greater characteristic speed values at points where $|u_x|$ is larger. This can lead to the intersection of characteristic curves, which corresponds to the wave breaking and the formation of a shock. This behavior can be studied using the method of characteristics. While for linear hyperbolic PDEs, such as the constant-coefficient wave equation $u_{tt} = c^2u_{xx}$ in the simplest case, characteristic curves can be found in terms of explicit formulas such as $x=x_0\pm ct$ and lead to explicit exact solutions, the situation is significantly more complex for nonlinear hyperbolic PDEs. Using the method described in [25, 49], one can show that the equation (5.10) can be reduced to the first-order characteristic form

(5.16) \begin{equation} u_t = \pm \dfrac{1}{2\sqrt{\epsilon}}\left(\sqrt{\epsilon}\, u_x\sqrt{1+\epsilon\left(u_x\right)^2} + \ln{\left(\sqrt{\epsilon}\, u_x + \sqrt{1+\epsilon\left(u_x\right)^2}\right)}\right) \end{equation}

on the characteristic curves

(5.17) \begin{equation} \dfrac{dx}{dt} = \pm \sqrt{1+\epsilon\left(u_x\right)^2}\,. \end{equation}

In physical terms, the part of the wave with a time derivative given by (5.16) moves at a finite velocity given by (5.17). The integration of (5.17) yields a constant of integration $x_0$ that corresponds to the point on the characteristic curve where $t = t_0$ is some initial time. Thus, $x = x(x_0, t)$ , from which $u_t(x(x_0, t), t)$ in (5.16) is a function of $x_0$ and t. Similarly, $u_x(x(x_0, t), t)$ in (5.17) is a function of $x_0$ and t.

The formation of a shock where the solution becomes multi-valued takes place when characteristic curves intersect. Without explicit knowledge of $u_x(x, t)$ , no explicit solution $x(x_0,t)$ of (5.17) is available. To estimate the breaking time $T_{b}$ , we use time-progressing linear approximations to characteristic curves, in conjunction with the finite-difference numerical solution of the PDE (5.10) described in Section 5.2. At each time layer $t_n$ in (5.13), linearized characteristics are launched forward in time from each grid point $(x_m, t_n)$ . The smallest time of the intersection of such characteristics estimates the wave breaking time.

For example, when the numerical computation has reached the time layer $t=t_n$ , the linearized characteristics are launched from each spatial grid point $x_m$ , $m = 1,\ldots,M-1$ . In particular, the right-traveling characteristics are approximated by the lines

(5.18) \begin{equation} x = x_m + t\sqrt{1+\epsilon\left(u_x(x_m, t_n)\right)^2}, \end{equation}

where $u_x=U_x$ is the first-order forward finite difference (5.14). To approximate the breaking time $t = \tau$ numerically, we solve (5.18) for the time $\tau$ when two different characteristics intersect. Given two starting points, $x_{m_1}$ and $x_{m_2}$ , we get the system

(5.19) \begin{equation} x = x_{m_1} + \tau\sqrt{1+\epsilon\left(u_x(x_{m_1}, t_n)\right)^2},\quad x = x_{m_2} + \tau\sqrt{1+\epsilon\left(u_x(x_{m_2}, t_n)\right)^2}.\end{equation}

Solving for $\tau$ yields

(5.20) \begin{equation} \tau = \dfrac{x_{m_1} - x_{m_2}}{m_2-m_1}, \end{equation}

where

(5.21) \begin{equation}m_i = \sqrt{1+\epsilon\left(u_x(x_{m_i}, t_n)\right)^2},\qquad i=1, 2\end{equation}

are the slopes of the characteristic lines. We choose $x_{m_1}$ and $x_{m_2}$ to be adjacent grid points, $x_{m_1}+\tilde{h} = x_{m_2}$ . The numerator of (5.20) is constant, so the approximate breaking time corresponds to the largest denominator of (5.20). We determine $x_{m_1}$ corresponding to the largest difference between the slopes $m_1$ and $m_2$ , then solve for $\tau$ .

The meaning of $\tau$ is thus the estimated time from $t_n$ to the wave breaking time $T_{b}$ ; one consequently has an estimate

(5.22) \begin{equation}T_{b}\sim t_n + \tau.\end{equation}

As the wave evolves, the slopes of the linearized characteristics will change and therefore so will $\tau$ . To account for this, at each time step, we repeat the calculation for $\tau$ . We use the first-order forward finite-difference approximation to compute $u_x$ in (5.21) at each time step. Figure 2(a) shows a plot of the value of the time to the break $\tau$ versus the time at which it was calculated, for several values of $\epsilon$ .

Figure 2. (a) Estimates (solid color) and numerical (black dashed) values for the time to break $\tau$ (5.20) for $\epsilon=0.1$ , $0.2$ , … , 1 (right to left). (b) Numerical wave profiles of $u_{xx}$ ( $\epsilon = 0.5$ ) for $t = 0, 1, 2$ .

Alternatively, one can numerically estimate the wave breaking time by defining it as the time when $\text{min}(u_{xx}) \leq \delta$ for some negative number $\delta$ . Choosing for example $\delta = -30$ , and using the second-order central difference approximation to the derivative for $u_{xx}$ in the numerical solution, we calculate the numerical break time for each $\epsilon$ . Using Richardson extrapolation of the break time values found with spatial step sizes $\tilde{h} = 0.01$ and $\tilde{h} = 0.005$ and temporal step size $\tilde{\tau} = 0.00125$ , we found the break times for each $\epsilon$ in the limit as $\tilde{h} \rightarrow 0$ as shown in Figure 2(a). This approach uses the fact that as the wave approaches the breaking time, the second spatial derivative $u_{xx}$ tends to negative infinity (see Figure 2(b)).

We determined the actual, non-linearized characteristic curves by numerically integrating (5.17). The curves are shown in Figure 3.

Figure 3. (a) Characteristic curves found by numerical integration of (5.17) with $\epsilon = 0.5$ . The blue horizontal line is the breaking time as determined by Richardson extrapolation. (b) The same plot as (a) with fewer characteristic lines shown. The thick red characteristic lines correspond to the earliest intersection.

It is interesting that one can also approximately determine the wave breaking time from the approximate solution (5.11) by finding the time when the second spatial derivative $u_{xx}$ of the approximate solution develops an additional root, as shown in Figure 4. This corresponds to an additional inflection point in the wave itself, which can be observed in Figure 1(a) closer to the wave breaking.

The numerically determined break time values ( $\tau_{num}$ ) and the approximate-determined break time values ( $\tau_{approx}$ ) are given in Table 2. Both sets of data are also plotted in Figure 5. The results suggest that the wave breaking time behaves like $\tau \sim \epsilon^{-1}$ .

Table 2. Numerical and approximate wave breaking time estimates for the PDE (5.10) vs. the small parameter values ( $\epsilon$ )

Figure 4. Wave profiles of the approximate solution $u_{xx}$ ( $\epsilon = 0.5$ ) for $t = 0, 1, 3, 3.5$ . Note the development of extra roots as time increases.

Figure 5. Numerical and approximate-derived break times $\tau$ as a function of $\epsilon$ .

6. Discussion

In this paper, we considered exact and approximate symmetries of partial differential equations with small parameter, in the Baikov–Gazizov–Ibragimov (BGI) [4–6] and the Fushchich–Shtelen (FS) [21] frameworks. The goals of this paper were to study relations, properties and applications of approximate symmetries of perturbed PDEs in comparison with exact symmetries of the unperturbed and perturbed PDEs and to use the approximate symmetries to construct approximate solutions of nonlinear PDEs.

For ordinary differential equations with a small parameter, all BGI approximate symmetries are stable: each symmetry of an unperturbed equation reappears as a point or higher-order approximate BGI symmetry [40]. For PDEs, the situation is different: we found that the point symmetry of the unperturbed PDE, in general, does not correspond to a point or a higher-order approximate symmetry of the perturbed model. The main reason of this symmetry instability is that the determining equations (2.23) for BGI local symmetries of the perturbed PDE (2.13), whatever the dependence of the approximate symmetry components $\zeta^1$ (2.22), always contains derivatives of u higher than those in $\zeta^1$ . This leads to a split system of PDEs in $\zeta^1$ with some restrictions on the unperturbed symmetry components $\zeta^0$ . A similar argument holds in FS framework. As an illustration, in Section 3.4, we showed that there was no higher-order (BGI and FS) approximate symmetry for the perturbed PDE (3.11) corresponding to the point symmetry $t \partial/\partial u$ of the unperturbed wave equation (2.9).

A detailed stability classification of point symmetries of the nonlinear wave equation (2.9) in terms of approximate BGI and FS frameworks was performed in Section 3.3. The results illustrate that the two approximate symmetry are significantly different. Yet as shown in Section 3.5, there exists a connection between the BGI and FS approximate symmetries; in particular, each stable BGI point symmetry of the form $(\zeta^0(x,t)+\epsilon \zeta^1(x,t,u,u_x,u_t)) \partial/\partial u$ yields a higher-order approximate FS symmetry in the form $\zeta^0(x,t) \partial/\partial v+ \zeta^1(x,t,v,v_x,v_t)\partial/\partial w$ (Theorem 3.1).

In the FS framework, we found that there exist FS approximate point symmetries of perturbed PDEs that do not correspond to the stable point symmetries of the unperturbed models and also cannot appear in the BGI framework (Remark 2.3). In such FS symmetries, the O(1) and $O(\epsilon)$ solution and symmetry components ‘mix’ together; in particular, the O(1) symmetry component depends on $O(\epsilon)$ solution component. An example is given by (2.30).

The classifications of exact and approximate (BGI and FS) symmetries of a perturbed PDE not only illustrates differences between these approximate symmetry approaches but also leads to different types of approximate symmetries that can be used to construct approximate solutions for the given PDE. In Section 4, exact and approximate point symmetries for the perturbed one-dimensional wave equation (4.6) arising in the context of nonlinear waves in elastic materials were classified. It was shown that the exact symmetry classification of the perturbed PDE (4.6) is included in the BGI approximate symmetry classification of (4.6). Moreover, the classifications of BGI and Fushchich–Shtelen approximate symmetries of (4.6) yielded different classification cases in terms of the arbitrary constitutive function. While the BGI classification included a logarithmic case that did not appear in the FS classification, two new cases appeared in the FS classification: $Q(v_x)=v_x^s$ and $Q(v_x)=e^{v_x}$ that did not arise in the BGI framework except for $s=1$ . All these cases yielded genuine BGI and FS approximate point symmetries for the nonlinear wave equation (4.6).

In Section 5, for the nonlinear wave PDE family (5.2) with a small power nonlinearity, an FS approximate symmetry generator was used to derive explicit families of approximately invariant solutions (5.9). For the case of quadratic nonlinearity (5.10), these approximate solutions were compared with finite-difference numerical solutions based on the method of characteristics. The latter indicate the appearance of breaking waves from smooth initial conditions in finite time. Though the approximate solutions remain smooth, their behavior can be used to predict the breaking times; these breaking times showed a good agreement with those estimated from the numerical computations. Both approaches indicate that the breaking time is inversely proportional to the magnitude of the small parameter that controls the nonlinearity.

The current work employed GeM module [15–17] for Maple symbolic manipulation software to generate and simplify determining equations for approximate symmetry and equivalence transformation components, the determining equations which were subsequently solved with standard Maple routines.

The approximate symmetry area still contains many open questions, including the general question of whether or not there exists some ‘unified’ approach that would include, and perhaps extend, both BGI and FS approximate symmetry frameworks. Another related open question is that of systematic computation of approximate conservation laws and approximately conserved quantities of dynamic ODEs and PDEs, and relations between such conservation laws and approximate symmetries. An obvious area of application of approximate conservation laws is the development and validation of specialized numerical methods for nonlinear DE models with small parameters.

A further extension of interest is the computation of approximate symmetries of DEs that include higher-order perturbation terms:

(6.1) \begin{equation} F[u;\epsilon]=F_0[u]+\epsilon F_1[u]+\ldots +\epsilon^N F_N[u]=o(\epsilon^N),\qquad N>1\end{equation}

(cf. (2.13)). Examples of such situation arise, for instance, in weakly nonlinear shallow water-type equations with very weak dispersion. For higher-order perturbations (6.1), it is natural to generalize the BGI symmetry generator (2.22)

(6.2) \begin{equation}\hat{X} = \hat{X}^0+\epsilon \hat{X}^1+\ldots+\epsilon^N \hat{X}^N,\end{equation}

and the FS solution representation (6.3) to

(6.3) \begin{equation} u(x) = v_1(x)+\epsilon v_2(x)+\ldots+\epsilon^N v_{N+1}(x) + o(\epsilon).\end{equation}

It is of interest to work out details and examples of approximate symmetry computations for such higher-order perturbations. An even further extension would include the consideration of singularly perturbed DEs where setting the small parameter $\epsilon=0$ changes the order of the differential equation or a PDE class it belongs to.

Another open area in the domain of approximate symmetries involves DE systems with more than one independent small parameters. Such models arise in multiple contexts. For instance, the dimensionless Serre–Su–Gardner–Green–Naghdi equations [39]

\begin{align*} {u}_{t}+ \varepsilon {u}{u}_{x}+h_{x} & = \dfrac{\delta^2}{3h}\left((h)^3\left({u}_{xt}+\varepsilon u{u}_{xx}-\varepsilon (u_{x})^2\right)\right)_{x},\\h_{t}+(h{u})_{x} & =0\end{align*}

involve two independent parameters $\varepsilon$ and $\delta$ ; a dimensionless viscoelastic wave model [13]

\begin{align*} u_{tt} &=\left (\alpha+3\beta u_x^2\right)u_{xx}+\eta u_x \left[\left(8u_x^2+2\right)u_{xx}u_{tx}+\left(2u_x^2+1\right)u_xu_{txx}\right] \nonumber\\&+ \zeta u_x^3 \left[\left(24u_x^2+4\right)u_{xx}u_{tx}+\left(4u_x^2+1\right)u_xu_{txx}\right]\end{align*}

includes four real parameters, three of which may be small: $\beta, \eta, \zeta \ll \alpha$ .