The class of second order quasilinear equations: models, solutions and background of classification
O. Makarenko A. Popov S. Skurativskyi




The class of second order quasilinear  equations: models, solutions and background of classification  


 O. Makarenko^†,[e-mail:
<makalex51@gmail.com>], A. Popov^,[e-mail:
<popovanton567@gmail.com>], S. Skurativskyi^,[e-mail: <
skurserg@gmail.com>] 


^†Institute for Applied System Analysis NTU ”KPI”

Politekhnichna st., 14, 14B,  Kyiv, Ukraine, 03056

^Institute of Physics and Technology NTU ”KPI”

Prosp. Peremohy, 37, Kyiv, Ukraine, 03056

^Subbotin institute of geophysics, Nat. Acad. of Sci. of Ukraine

     Bohdan Khmelnytskyi str. 63-G, Kyiv, Ukraine





 Abstract.The paper is concerned with the unsteady solutions to the model of mutually penetrating continua and quasilinear hyperbolic modification of the Burgers equation (QHMB). The studies were focused on the peculiar solutions of models in question. 
On the base of  these models and their solutions, the ideas of second order quasilinear models classification were developed.  



 Keyword: hyperbolic equations, attractors, multivaluedness
 







§ INTRODUCTION


Nonlinear models for phenomena and systems are the cornerstones in modern physics. The examples of these models are well known, namely the Burgers, KdV, Liouville, nonlinear Shredinger  equations and etc. The derivation of analytical solutions for these equations is a challenge for scientists. Therefore, the numerical treatments of such models   are developed intensively. 

Among the models mentioned above it is worth  accentuating the second order in time nonlinear hyperbolic differential equations <cit.> having a broad applications recently. Note that the partial solutions for these equations had been derived <cit.>.

Due to the significance of considered equation and their solutions, this paper deals with the classification of these equations, definition of general expression of quasilinear models. Models' solutions obtained via the numerical modelling are analyzed in detail. We also discuss the possible ways of investigations, in particular,  combining the concept of dynamical systems (attractors) and artificial intelligence methods (neural networks). The aspects related to the multivaluedness of solutions including symmetries are considered as well. 
     

	
	

§ WAVE REGIMES IN MEDIA WITH OSCILLATING INCLUSIONS

	
To begin with, let us note that the class of second order quasilinear models is not empty  and covers  many different  models originated from the physics and biology. In particular, consider the model of mutually penetrating continua which uses for the description of physical processes in complex media <cit.>. It turned out that this model possesses the specific wave solutions. Consider these solutions in more detail. 

The model we are going to deal with has the following form  

    ρ∂ ^2 u/∂ t^2  = ∂σ/∂ x - m ρ∂ ^2 w
    /∂ t^2 ,     ∂ ^2 w /∂ t^2
     + Φ( w - u) = 0,

where ρ is  medium's density, u and w are the displacements of carrying medium and oscillator from the rest state,  mρ is the density of oscillating continuum.  We also use  the cubic constitutive equation for the carrying medium
 σ  = e_1 u_x  + e_3 u_x^3, where e_1, e_3 are the  elastic moduli, and relation for applied force Φ(x)=ω ^2 x+δ x^3, where ω  denotes the natural frequency of oscillator, whereas  the parameter δ appears  due to accounting for the cubic term in the expansion of restoring force in a power series. 


The  traveling wave solutions of model (<ref>) have the following form

    u=U(s),    w=W(s),    s=x-Dt,

 where the parameter D
 is a  constant velocity of the wave front. 

 Inserting (<ref>) into  model (<ref>), it is easy to see that  the   functions U and W
 satisfy the dynamical system 

    D^2 U' = ρ^-1σ( U') - mD^2 W^',
       W” + Ω ^2 ( W - U) +δ D^-2(W-U)^3= 0,

where Ω  = ω D^-1.


This system can be written in the form 

    W' =α _1 R + α _3 R^3,     
      U' =R,      (α_1+3α_3 R^2)R'+Ω^2(W-U)+δ D^-2(W-U)^3=0,

 where α _1  =
e_1 - D^2 ρ/m ρ D^2, α _3
= e_3 /m ρ D^2 >0. Through the report we fix 
    e_1=ρ=1,  e_3=0.5,  m=0.6,  ω=0.9
 in numerical treatments.


At first, consider system (<ref>) at δ=0. 
The detail description  of  phase plane of dynamical system had been done in the paper <cit.>, we thus summarize 
the main results only.

At α_1<1 the phase plane contains three fixed points, whereas at α>1  the only one fixed point (center) remains. 
For  α_1<0, when  D=1.2 is fixed for definiteness,  a typical phase portrait is depicted in the Fig. <ref>a. In this case, all fixed points are centers surrounded
by periodic orbits. There are  separatrices separated the regions with periodic and unbounded trajectories and two lines corresponding to discontinuity of system. 



When 0<α<1, at D=0.9 for instance, in the phase portrait (Fig. <ref>b) one can distinguish   the homoclinic trajectories that go through 
the origin. The homoclinic loop can be written in the explicit form 

    [                                                                                         s - s_0  = 3 /2Ωarcsin(
    4α _3 R^2 - 3 + 4α _1 /√(9 - 8α _1 )) - 1/2Ω√(α _1 /1 - α _1 )×;                                                                                                                                                                                ; ln(
     1/R^2 + 3α _3  - 4α _1 α _3 /4( α _1  - α _1^2 ) +
     √({1/R^2 + 3α _3  - 4α _1 α _3 /4( α _1  - α _1^2 )}^2  - α _3^2 ( 9 - 8α _1 )/16( α _1  - α _1^2 )^2 ))  |_R_0^2 ^R^2 ]


 
 Solution (<ref>) corresponds to the solitary wave solution
with infinite support.

When α_1 tends to zero, the angles between separatrices of saddle point O are growing. As a result,   at  α _1 = 0  
we observe  the transformation of  solitary wave into the compacton, i.e.,  solutions with finite support <cit.>. 
These orbits are described by the following expressions

    U^' _s=√(3/2α _3 )sinΩ s/3,      Ω s/3∈[ 0;π]  
    
        0,      Ω s/3∉[ 0;π]  
    
        
    U=
    0,      Ω s/3∈[ -∞; 0  ],
    3/Ω√(3/2α_3)(1-cosΩ s/3),      Ω s/3∈(0; π],
    6/Ω√(3/2α_3),      Ω s/3∈(π; ∞).




As above, there are  the periodic orbits enclosed in the homoclinic loops  and  periodic trajectories l lying  beyond the homoclinic contour.



 §.§ Phase diagrams in the model with cubic nonlinearity in the equation of motion for oscillating inclusions


If δ≠ 0, then system (<ref>) does not reduce to the dynamical system in the plane (R;R'). But the first integral  for  (<ref>)  can still be derived in the form

    I=μ_1/2(W-U)^4 + μ_2(W-U)^2+ f(R),

where μ_1=δ D^-2, μ_2=ω^2 D^-2, f(R) = α _3^2 R^6  + 4α
_1 - 3/2α _3 R^4  + ( α _1^2  - α _1 )R^2. Since dI/ds=0 on the trajectories of system (<ref>), then I ≡.

It is easy to see that expression (<ref>) can be used for splitting  system  (<ref>). Indeed, solving (<ref>)
with respect W-U we obtain


    W-U=±√(-μ_2 ±√(μ_2^2-2μ_1 f(R)+2μ_1 I)/μ_1).


This allows us to separate the last equation of (<ref>) from other ones. Unfortunately, the resulting equation cannot be 
integrated for general case, therefore,  let us consider its phase plane structure, which is equivalent to the structure of  level curves for function (<ref>).





Now consider the position and the form of homoclinic trajectories when the parameter δ is varied. 
Starting from the loop 
at δ=0 which coincides with the orbits of Fig. <ref>b, we see that increasing δ causes the attenuation
 of  loop's size along vertical axis.  Actually, the level curve consists of the closed curve (homoclinic loops) and unbounded trajectories. 
If δ decreases,  loop's size grows, but at δ_0 the additional heterocycle connecting four new saddle points appears.
 The bifurcational value δ_0  can be derived via  analyzing the function (<ref>). Namely, δ_0 corresponds to the moment when different branches of  (<ref>)
are tangent. This happens when μ_2^2-2μ_1 f(R)=0. Thus, the condition of contact for two branches leads us to a cubic equation with respect to R^2  with zero discriminant. 
Then δ_0=-α_3ω^4 /D^2(α_1-1)^2 or δ_0=27 α_3ω^4 /D^2α_1^2 (9-8α_1). The last value of δ_0
is not interesting because four branches of (<ref>) degenerate into two ones forming homoclinic loops.

Considering the first value δ_0, we put D=0.9 and derive  δ_0=-2.24. For δ>δ_0 we have homoclinic loops placed along horizontal axis accompanied by  
appearing  the unbounded curves  in the upper and lower parts of diagram. When δ<δ_0, the homoclinic loops are placed in the vertical quarters of the phase plane, whereas the unbounded orbits appear 
at the left and at the right sides of diagram (Fig. <ref>a).
Note that the profiles of the resulting solitary waves are different (Fig. <ref>b), namely, at δ<δ_0 the U profile looks like a bell-shape curve, 
but at δ>δ_0 it is a kink-like regime.  

The homoclinic orbits divide the phase plane into parts filled by closed curves corresponding to the periodic regimes.  If we choose  δ=-1.5>δ_0, we get  the typical phase 
portrait of system (<ref>) plotted in 
the Fig. <ref>a. In the portrait  two pairs of nontrivial fixed points  A_± (± Q;0) and B_± (0; ±ω/√(-δ)) can be distinguished. Inserting the coordinates of these points into 
relation (<ref>),  we obtain the values of I_1=-(α_1-1)^2/2α_3 and I_2=ω^4/2D^2|δ| which allows us to state the conditions of periodic regimes existence. 
For fixed δ=-1.5, I_1=-0.18 and I_2=0.27. We thus get that periodic regimes exist if I_1<I<I_2 only. 

Now let us choose  δ=-2.7<δ_0. In this case the homoclinic loop going through the origin is placed along vertical axis and  the phase portrait looks like  Fig. <ref>b. As above, we have I_1=-0.18, but I_2=0.15.
 





 §.§ Wave dynamics of   model (<ref>)


To model the wave dynamics, we used  the three level
finite-difference numerical scheme for model (<ref>). 



Solitary waves. To consider the evolution of solitary waves, let us construct the initial data v_i, q_i, G_i, F_i for numerical simulation on the base of homoclinic contour. 
To do this, we integrate dynamical system (<ref>) with initial data R(0)=10^-8, Z(0)=0, s ∈ [0;L]  and choose the right  homoclinic loop in the phase portrait (Fig. <ref>b). 
Then the profiles of W(s), U(s), and R(s) can be derived. Joining the proper arrays, we can build the profile in the form of  arch:

    v=U(ih)∪ U(L-ih),  q=u(x+τ D)=[U(ih)+τ D R(ih)]∪[U(L-ih)+τ D R(L-ih)].



The arrays G and F are formed in similar manner. 
Combining two arches and continuing the steady solutions 
at the ends of graph, we get more complicated profile.
 We apply the fixed boundary conditions, i.e. u(x=0,t)=v_1,   u(x=Kh,t)=v_K, where K is the length of an array.




Starting from the two-arch initial data, we see (Fig. <ref>) that solitary waves move to each other, vanish during approaching, and appear with negative
 amplitude and shift of phases. After collision in the zones between  waves some ripples are revealed.  Secondary collisions of waves are watched also. Note that the
 simulation of compacton solutions displays similar properties. 






Propagation of solitary waves at δ≠0 depends on the sign of δ. Numerical simulations  show that the 
collision of waves at  δ>0 is similar to the collision at  δ=0. Behavior of waves after interaction does not change essentially  when  δ<0 and close to zero (Fig. <ref>a). 
But if  δ is not small,  after collision the amplitude of solution is increasing in the place of soliton's intersection  and  after a while   the solution  is  destroyed <cit.>. 

 This suggests that we encounter  the unstable  interaction of solitary waves or the numerical scheme we used possesses spurious solutions. But if we take 
half spatial step and increase the scheme  parameter r up to 0.8, the scenario of solitary waves collision is not changed qualitatively. Therefore, the assertion on the unstable nature of collision is more preferable.

Periodic waves. To simulate the evolution of periodic waves, we take the initial profile corresponding to a periodic orbit (red curve in Fig. <ref>b)
surrounding the homoclinic loop at the phase portrait of Fig. <ref>a   at δ=-0.2 and D=0.9.  Due to the periodicity of the problem, 
the boundary conditions for numerical scheme  should be modified <cit.>.
The resulting profile derived after time interval 1200τ  has passed is depicted in Fig. <ref>b with dark green curve. It is obvious that this wave is shifted  
a distance 1200τ D to the right in a self-similar manner.



§ UNSTEADY SOLUTIONS TO QHMB


	Let us consider the evolution of localized perturbation within the framework of 	QHMB:
	
    τ∂ ^2 u/∂ t^2+α∂ u/∂ t+βφ(u)∂ u/∂ x=μ k(u)∂^2 u/∂ x^2+νψ (u) (∂ u/∂ x)^2+θ f(u),

	where τ, α, β, μ, ν, θ are constants.
	
	As an initial data for equation (<ref>)б we chose the profile 
	u(x,0)=exp(x+a/b)+exp(x-a/b) (fig. <ref>a), where a=5, b=4.5. The evolution of this two-hump perturbation is shown in figs. <ref>b,c,d.
	
	
	
From the analysis of these figures it follows that  the number of oscillations increases during evolution. The oscillations have non-regular character but remain localized in spatial domain. Thus, one can call them the “pre-turbulent oscillations” or, due to the localization, the “turbulons”. Let us also note that the appearance of this type solutions was discovered in other numerical simulations of 2D and 3D problems <cit.>.
	


§ THE CLASS OF QUASILINEAR EQUATIONS AND THEIR FORMAL DESCRIPTION



The presented examples of modelling the carrying processes manifest the variety of solutions of models considered. Analyzing these models in details, one  can be convinced that they belong to the second order quasilinear equations, i.e., the relation is a linear with respect to higher derivatives whereas it is nonlinear with respect to lower derivatives, unknown functions,  and independent variables.    It is obvious that this class of equations is very broad and needs to systematize.

Let us start from the defining the general expression of quasilinear second order models, namely 

    τ r(u) ∂ ^2 u/∂ t^2   +α s(u) ∂  u/∂ t+βφ (u) ∂ u/∂ x = μ k(u)  ∂ ^2 u/∂ x^2 +νψ (u) ( ∂  u/∂ x)^2 + γ h(u)  ∂ ^2 u/∂ x∂ t +
       +ξ b(u) ( ∂  u/∂ t)^2+θ f(u)+ χ I(u,x,t,±Δ t, ±Δ x;...),

where I is the source (integral form addmited), 
τ, α, β, μ, ν, θ are constants, r, s, k, ψ, h, b, f are the specified functions.

Note that, depending on the coefficients and functions, class (<ref>) covers the classical linear equations (elliptic, parabolic and hyperbolic) and well-known nonlinear models (sin-Gordon, Burgers, Liouville,  Hopf, hyperbolic modification for the Burgers equation, equations with blow-up solutions).

For future handling the variety of models and their solutions, we need a convenient method for their identification. We propose the following descriptor for such objects:
 
    EQ(   τ,α,β,μ,ν,γ,ξ,θ,χ; 
    {   NONLINEAR_FUNCTION: r,s,φ,k,ψ,h,b,f,I};
    {   TITLE_OF_EQUATIONS};{TYPES_OF_SOLUTIONS})


For instance, the descriptor for QHMB can be chosen in the  form
 
    EQ(   1, 1, 1, 1, 1, 0, 0, 1, 0; 
    {   NONLINEAR_FUNCTION: φ, k,ψ, f };
    {   QHMB};{compactons;blow-up;oscillations}).

When the system of equations is considered, the coefficients and functions in (<ref>) should be assumed as  matrices. In particular, rewrite system (<ref>) as follows

    ∂ ^2 u/∂ t^2    =ρ^-1(e_1+3 e_3(∂ u/∂ x)^2) ∂
    ^2 u /∂ x^2 + m (ω^2(w-u)+δ(w-u)^3), 
    ∂ ^2 w /∂ t^2
       =-ω^2(w-u) - δ(w-u)^3.

Then  its descriptor   is
 
    EQ(   ([ 1 0; 0 1 ]),
    0̂,0̂, ([ ρ^-1    0;    0    1 ]),0̂ ,0̂,0̂,([  m  0;  0 -1 ]) ,0̂; 
    {   NONLINEAR_FUNCTION: 
       k=([ e_1+3 e_3(∂ u/∂ x)^2;                    0 ]), f=([ ω^2(w-u)+δ(w-u)^3; ω^2(w-u)+δ(w-u)^3 ])};
    {   EXAMPLE_ONE};{solitary  waves; compactones}).


Remark that the construction of descriptor helps us not only to classify the equations but at the statement of  problems. In particular, the Hopf equation can be considered as a limit case of hyperbolic Burgers equation.  

Let us also outline the types of solutions of such models. For generality, it is worth to mention that  certain equations from class (<ref>) possess interesting solutions, namely, autowave solutions to the Burgers equation,  singular solutions to the  Liouville equation, blow-up solutions to parabolic and hyperbolic equations with nonlinear sources, solitons for the sin-Gordon equation. Among the solutions to the  hyperbolic modification of the Burgers equation we encounter the packets of oscillations, compactons, autowaves.     

Another type of solutions is related to the multivaluedness. This means that the solutions have  several values in   given spatial point at fixed moment of time <cit.>. In particular, the Hopf equation  admits the solutions when their profiles  “overturn"  and becomes multivalued. Similar behavior of wave solutions is observed in the Vakhnenko equation which admits the loop solutions. Studies of such solutions cause the introduction of new class of solutions named the “foldons”. Foldons are the multivalued   autowaves. It is important that both the Hopf equation and Vakhnenko equation belong to class (<ref>) and their  multivalued solutions have analytical expressions. 

Up today we have hardly any methods to treat multivalued models. In this case the group analysis methods seem to be useful. To apply group methods, the multivalued differential equations should be considered as multivalued surfaces (geometric objects). The promising approach of multivalued model studies is concerned with the asymptotic transition from the singlevalued to multivalued models.  We encounter this  when transform the hyperbolic Burgers equation to the Hopf equation. Another way to  treat the multivaluedness is the expanding the inverse scattering transform.

It is worth to mention the development of novel approach dealing with the specific type of  inverse problems. The examples considered in the paper are concerned with analyzing the solutions  and their dependence on the parameters when the equations subjected to the initial and boundary data are available. Traditionally, such problems have been solved via the trial and error methods to find the self-similar solutions. But the inverse statement of the problem is possible. Indeed, we can take the functions ũ (x,t) in such a way that they satisfy equation (<ref>), when the components in  (<ref>)  are chosen correctly. 
In general case we would like to have the method when the set of functions ũ_1(x,t), ũ_2(x,t), ... , ũ_n(x,t) satisfies  equation  (<ref>) subjected to proper set of initial conditions. This problem is similar to approaches in the theory of neural networks. The initial and boundary conditions are regarded as inputs of neural networks, whereas the functions ũ_1(x,t), ũ_2(x,t), ... , ũ_n(x,t) considered as potential solutions are identified as outputs of neural networks. It is obvious that for arbitrary equation with fixed initial conditions the solutions    u_1(x,t), u_2(x,t), ... , u_n(x,t) differ from the  “desired” set of function ũ_1(x,t), ũ_2(x,t), ... , ũ_n(x,t). Then one can  find the values of equation's coefficients providing the small deviation Δ=ũ_i(x,t)-u_i (x,t). To realize this procedure, the iteration process of approaching to the equation's coefficients can be applied. This process is an analog of studying process in neural networks. We can use multivalued neural networks in these problems as well.
  



§ CONCLUSION


Summarizing, we presented the specific solutions occurred in quasilinear models for the carrying processes in nature. The results obtained were encouraged us to develop the systematic approach to the studies of nonlinear dynamical models within the framework of quasilinear second order equations and their solutions. We outlined  new ways to analyze quasi-linear models including the multivalued cases. 

00


monog
Danylenko V.A., Danevych T.B.,  Makarenko O.S.,  Skurativskyi S.I.,  Vladimirov V.A., Self-organization in nonlocal
non-equilibrium media,  Kyiv, Subbotin in-t of geophysics NAS of Ukraine,
2011.

Makarenko_2012
 Makarenko A. S., Model equations and formation of structures in media with memory. Ukr. J. Phys. 57(4) (2012) 408–421.

Popov_Mak_2015
Popov A.S.,  Makarenko O.S., To the numerical solution of quasilinear  modification for the Burgers equation. Analysis, modelling, management.  Vol.2: Kyiv, NTTU”KPI”, 2015.  – P. 197-203.

Mak_2015_Austria
Makarenko O. S., Toward multivaluedness aspects in self-organization, complexity and computations investigations. Forth Int. Workshop on Nonlinear Dynamics and Sinchronization INDS’15, July 31, Klagenfurt, Austria, Alpen-Adria University, 2015. –  Pp. 84-93





nasu
  Danylenko V.A.,  Skurativskyi S.I, Resonance regimes of the spreading of nonlinear wave fields in media with oscillating inclusions,
  Reports of NAS of Ukraine,
    11  (2008) 108–112.

nonlin_dym
 Danylenko V.A., Skurativskyi  S.I.,  Travelling wave solutions of  nonlocal models  for media with oscillating inclusions,
  Nonlinear Dynamics and Systems Theory,
4(12) (2012) 365–374.



dan_lodz_16
Danylenko V.A., Skurativskyi S.I.,  On the dynamics of solitary wave solutions supported by
the model of mutually penetrating continua,   Dynamical systems. Mechatronics and life sciences, 2 (2015) 453–460. (arXiv:1512.05226v1 [nlin.PS] 15 Dec 2015)


dan_skur_ufj
Danylenko V.A.,  Skurativskyi S.I., Skurativska  I.A.,  Asymptotic 	wave	solutions	for	the	model of a medium 	with	Van	Der	Pol	oscillators, Ukrainian	Journal	of	Physics,	
59(9) (2014) 932–938.

skur_akust
 Skuratovskii S.I., Skuratovskaya I.A., Localized autowave solutions of the nonlinear model of complex medium,
  Electronic Journal ”Technical Acoustics”,
  6 (2010)  http: ejta.org.




Hyman-Rosenau
Hyman J. and  Rosenau Ph., Compactons: solitons with finite wavelength, Phys. Rev. Lett.,
   70  (1993) 564.


vsan_lodz2016
Vladimirov V. A. and Skurativskyi S. I., Solitary  waves in  one-dimensional pre-stressed lattice and its 
continual analog, Dynamical systems. Mechatronics and life sciences, 2  (2015) 531–542. (arXiv:1512.06125v1 [nlin.PS] 18 Dec 2015)


Dan_Skur_2016
 Danylenko V.A. and  Skurativskyi S.I., Peculiarities of wave dynamics in media with oscillating inclusions, 
International Journal of Non-Linear Mechanics, 84 (2016)  31–38.




