On nonlinear population waves
Applied Mathematics and Computation 215 (2009) 2950–2964
Contents lists available at ScienceDirect
Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc
On nonlinear population waves Nikolay K. Vitanov a,*, Ivan P. Jordanov a, Zlatinka I. Dimitrova b a b
Institute of Mechanics, Bulgarian Academy of Sciences, Akad. G. Bonchev Str., Bl. 4, 1113 Sofia, Bulgaria Institute of Solid State Physics, Bulgarian Academy of Sciences, Blvd. Tzarigradsko Chaussee 72, 1784 Sofia, Bulgaria
a r t i c l e Keywords: Nonlinear dynamics Population dynamics Migration Solitary waves Kinks Asymptotic stability
i n f o
a b s t r a c t We discuss a model system of partial differential equations for description of the spatiotemporal dynamics of interacting populations. We are interested in the waves caused by migration of the populations. We assume that the migration is a diffusion process influenced by the changing values of the birth rates and coefficients of interaction among the populations. For the particular case of one population and one spatial dimension the general model is reduced to analytically tractable PDE with polynomial nonlinearity up to 4th order. We investigate this particular case and obtain two kinds of solutions: (i) approximate solution for small value of the ratio between the coefficient of diffusion and the wave velocity and (ii) exact solutions which describe nonlinear kink and solitary waves. In an appropriate phase space the kinks correspond to a connection between two states represented by a saddle point and a stable node. Finally we derive conditions for the asymptotic stability of the obtained solutions. Ó 2009 Elsevier Inc. All rights reserved.
1. Introduction The dynamics of many systems in physics, chemistry, biology, fluid dynamics, economics or society is nonlinear [1–9]. If the nonlinearity compensates for the dispersion nonlinear traveling waves can propagate through the corresponding system [10–14]. A classical example for nonlinear dynamics is the dynamics of interacting populations [15,16]. One important aspect of this dynamics is the movement of humans or animals and this movement can lead to large solitary migration waves [17,18]. The solitary waves are very popular and frequently obtained in various areas of natural sciences [19–23]. In the last years many authors have obtained exact solutions which describe solitary waves, solitons, foldons, etc. in various equations and system of equations. Just several examples are the Korteweg–de Vries (KdV) system [24–26], Ito equations [27], Boussinesq equation [28], Korteweg–de Vries equation [29–33], Kuramoto–Sivashinsky equation and its generalizations [34–38], sine-Gordon and sinh-Gordon equations [39–44], Kadomtsev–Petviashvili equation [45,46], Nizhnik–Novikov–Veselov system [47,48], Broer–Kaup equation [49–51], nonlinear Schrödinger equation [52], Maxwell–Bloch system [53] or Camassa– Holm equation [54]. In this paper we shall concentrate our attention on the analytical investigation of one-dimensional migration waves caused by difference in the spatial density of the members of one population. The general model for n interacting populations on two-dimensional surface is described in Section 2. In Section 3 we discuss the case of one population and one spatial dimension and obtain exact particular solutions of the studied model equation and its particular cases. Section 4 is devoted to approximate solution of the model equation for the case of large values of the ratio between the wave velocity and diffusion coefficient. In this section we derive conditions for asymptotic stability of the obtained solutions. Several concluding remarks are summarized in Section 5. * Corresponding author. E-mail addresses: [email protected], [email protected] (N.K. Vitanov). 0096-3003/$ - see front matter Ó 2009 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2009.09.041
2951
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
2. Formulation of the problem We consider a two-dimensional plane populated by members of n populations. Let us select a small region Ds ¼ DxDy in this plane. In this region we observe DN i members from the ith population ði ¼ 1; 2; . . . ; nÞ. Then if Ds is sufficiently small we can define the density of the ith population
qi ð~r; tÞ ¼
DN i : Ds
ð1Þ
We assume that the members of the ith population can move through the borders of the small area under consideration. r; tÞ of the current of this movement. The net number of the members of the ith population, Then we can define the density~ ji ð~ ~ Þdl. Within the observed area the number of the members of the ith ~ is ð~ crossing a small line dl with normal vector m ji m population can change by births, deaths or interactions with other populations (one example is the interaction among predators and preys). These processes will be summarized mathematically by the function R below. The total change of the number of members from the ith population in some large area is a consequence of migration through the boundaries of the area and the birth, death and interaction processes within the area. As a mathematical expression this reads
@ @t
Z S
dsqi ð~ r; tÞ ¼
Z
~þ dl~ji m
Z
l
dsRi :
ð2Þ
S
The first term in (2) describes the net rate of increase of the density of the ith population, the second term describes the net rate of immigration into the area and the third term describes the net rate of increase exclusive of immigration. The differential form of Eq. (2) is
@ qi þ r ~ji ¼ Ri : @t
ð3Þ
Now we have to specify ~ ji and Ri . Ri can be taken in the general form analogous to the form proposed by Dimitrova and Vitanov [55–59]
" Ri ¼ r 0i qi 1
n X
#
a0ij qi þ Ai ;
ð4Þ
j¼1
where
Ai ¼
n X
r ij qj
j¼1
n X
a0ij ðaijk þ ril Þqi qj
j;l¼1
n X
a0ij rik aijl qj qk ql
ð5Þ
j;k:l¼1
and n is the number of the interacting populations. Eq. (4) is obtained as follows: (1.) Similar to the generalized Lotka–Volterra equations we assume that the changes of the density because of reproduction and interaction between the populations are
" Ri ¼ r i qi 1
n X
#
aij qj ;
j¼1
where ri is the growth ratio of the ith population and aij is the interaction coefficient measuring to what extent the growth of the ith population is influenced by the jth population. (2.) Growth ratios and interaction coefficients depend on the density of members of the populations
" ri ¼
r 0i
1þ
n X
# rik qk ;
k¼1
" 0 ij
aij ¼ a 1 þ
n X
#
aijl ql :
l¼1
We shall assume that ~ j has a form similar to the general form of the linear multicomponent diffusion [60]
~ji ¼
n X
Dik rqk ;
ð6Þ
k¼1
where Dik ¼ Dik ðqi ; qj ;~ r; tÞ. We note that in the population dynamics it is not necessary that Dik ¼ Dki . Below we shall treat Dik as constants. After the substitution of (4) and (6) in (3) we obtain the following model equations
" # n n n n X X X @ qi X 0 0 0 0 Dij Dqj ¼ r i qi 1 ðaij r ij Þqj aij ðaijk þ rik Þqj qk aij rik aijl qj qk ql ; @t j¼1 j¼1 j;k¼1 j;k;l¼1
D¼
@2 @2 þ 2: 2 @x @y
ð7Þ
The fixed points of the system (7) correspond to stationary densities q01 ; q02 ; . . . ; q0n . The stationary densities qi0 are qi0 ¼ 0 as well as the solutions of the system of equations
2952
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
1
n X
a0ij rij qj 0
j¼1
n X
a0ij ðaijk þ rik Þqj 0 qk0
j;k¼1
n X
a0ij rik aijl qj 0 qk0 ql0 ¼ 0:
ð8Þ
j;k;l¼1
Let us consider a deviation Q i of the density qi of the ith population from the corresponding stationary density q0i as follows:
qi ðx; y; tÞ ¼ Q i ðx; y; tÞ þ q0i ; i ¼ 1; 2; . . . ; n:
ð9Þ
The substitution of (9) in the right-hand side of (7) leads us to
" n n n X X @ qi X 0 Dij Dqj r i ðQ i þ q0i Þ a0ij rij Q j þ a0ij ðaijk þ rik ÞðQ j Q k þ Q j q0k þ Q k q0j Þ @t j¼1 j¼1 j;k¼1 þ
n X
#
a0ij rik aijl ðQ j Q k Q l þ Q j Q k q0l þ Q k Q l q0j þ Q l Q j q0k þ Q j q0k q0l þ Q k q0l q0j þ Q l q0j q0k Þ :
ð10Þ
j;k;l¼1
For the case of one population and one spatial dimension (7) reduces to the following model equation:
@ q1 DDq1 ¼ r0 q1 1 ða0 r 11 Þq1 a0 ða111 þ r 11 Þq21 a0 r 11 a111 q31 : @t
ð11Þ
The fixed points of (11) are
q10 ¼ 0; q10
1 ¼ ; r 11
q10
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 1 a0 a0 þ 4a0 a111 ¼ : 2 a0 a111
ð12Þ
We are interesting in the behavior of perturbations (linear and nonlinear) around the basic state q0 , i.e.
q1 ðx; tÞ ¼ q0 þ Q ðx; tÞ:
ð13Þ
The substitution of (13) in (11) leads to the following nonlinear equation for Q
@Q @2Q D 2 ¼ EQ 4 þ FQ 3 þ GQ 2 þ HQ; @t @x
ð14Þ
where
E ¼ r 0 a0 r 11 a111 ; 0
F ¼ r 0 a0 ða111 þ r11 þ 4q0 r 11 a111 Þ;
0
G ¼ r ½r11 a 3q0 a0 ða111 þ r11 þ 2r 11 a111 q0 Þ; H ¼ r 0 ½1 q0 ð3a0 r 11 q0 þ 4a0 r 11 a111 q20 2r11 þ 2a0 þ 3a0 a111 q0 Þ: Below we shall treat (14) as general form of a model equation which for the particular case when E; F; G; H satisfy the above relationships describes population density waves. In the following section we shall obtain and discuss: (i) traveling wave solutions of (14) and (ii) traveling wave solutions of important equations connected to (14). 3. Traveling wave solutions of (14) and its particular cases 3.1. Case E – 0; F – 0; G – 0; H – 0 We consider Eq. (14) and search for solutions of the kind of traveling waves: Q ðx; tÞ ¼ Q ðnÞ ¼ Q ðx v tÞ where velocity of the wave. We have to solve the equation 2
D
d Q dn2
þv
dQ þ EQ 4 þ FQ 3 þ GQ 2 þ HQ ¼ 0: dn
v is the ð15Þ
Following [61] we assume that Q ðnÞ has the form
QðnÞ ¼
n X i¼0
ai /i ;
vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi uX u r d/ cj /j ; ¼ t dn j¼0
ð16Þ
where ¼ 1, and aj and cj are parameters that we shall determine below. (For more results and methods for obtaining exact solutions of nonlinear PDE see for an example [62–65].) We note that the methodology based on (16) is closely connected to the methodology of Kudryashov [65,66] which will be much used below. Substituting (16) in (15) we obtain an equation that contains powers of /. We have to balance the highest power arising from the second derivative in (15) with the highest power arising in the term containing Q 4 in the same equation. The resulting balance equation is
r ¼ 3n þ 2;
n ¼ 2; 3; . . .
ð17Þ
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
2953
If n ¼ 2 then r ¼ 8. Assuming that c0 ¼ c1 ¼ c3 ¼ c4 ¼ c6 ¼ c7 ¼ 0, c2 ¼ p2 ; c5 ¼ 2pq; c8 ¼ q2 we search for a solution of the form
Q ðnÞ ¼ a0 þ a1 / þ a2 /2 ;
d/ ¼ ðp/ þ q/4 Þ: dn
ð18Þ
Substituting this in (15) we obtain the following system of nine equations for the parameters a0 ; a1 ; a2 ; p; q; D; E; F; G; H
a2 10Dq2 þ Ea32 ¼ 0; a1 Dq2 þ Ea32 ¼ 0; a22 Fa2 þ 4Ea2 a0 þ 6Ea21 ¼ 0; a2 2v q þ 14Dpq þ 3Fa1 a2 þ 12Ea1 a2 a0 þ 4Ea31 ¼ 0; 5Da1 pq þ Ga22 þ v a1 q þ 3Fa22 a0 þ 3Fa21 a2 þ 6Ea20 a22 þ 12Ea2 a0 a21 þ Ea41 ¼ 0;
v;
ð19Þ
a1 ð2Ga2 þ 6Fa0 a2 þ Fa21 þ 12Ea20 a2 þ 4Ea0 a21 Þ ¼ 0; 3Fa20 a2 þ 3Fa0 a21 þ 2Ga0 a2 þ Ga21 þ 4Ea30 a2 þ 6Ea20 a21 þ Ha2 þ 4Da2 p2 þ 2v a2 p ¼ 0; a1 ð3Fa20 þ H þ 2Ga0 þ 4Ea30 þ Dp2 þ v pÞ ¼ 0; Ga0 þ H þ Fa20 þ Ea30 ¼ 0: Eq. (19) implies that a1 ¼ 0. Let in addition
¼ 1. The solution of the system is
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 10Dq2 E2
F a0 ¼ ; a2 ¼ 4Epffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 7 F 10DF v¼ ; 80 E
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 F 10DF ; p¼ ; E 80 DE 2 3 3F 3F G¼ ; H¼ 8E 64E2
ð20Þ
and q; F; E; D are free parameters. The expression for the solitary wave depends on the solution of the differential equation in (18) and it is given by
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 3 p exp½3pðn þ cÞ Q ðnÞ ¼ a0 þ a2 1 q exp½3pðn þ cÞ
ð21Þ
for the case ½/3 =ðp þ q/3 Þ > 0 and
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 3 p exp½3pðn þ cÞ Q ðnÞ ¼ a0 þ a2 1 þ q exp½3pðn þ cÞ
ð22Þ
for the case ½/3 =ðp þ q/3 Þ < 0. Above c is a constant of integration. Fig. 1 shows a kink of kind (21) as well as its derivative which is a well formed solitary wave. For the discussed case of modeling of the movement of population members and for q0 ¼ 0 we obtain
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 10Dq2 a2 ¼ ; r 0 a0 r 11 a111 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 r 11 þ a111 p¼ 10Dr0 a0 ðr11 þ a111 Þ; 80 Dr 11 a111 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 7 a111 þ r 11 v¼ 10Dr 0 a0 ða111 þ r 11 Þ; 80 a111 r 11 r11 þ a111 a0 ¼ ; 4r 11 a111
where the parameters r 11 ; a0 and a111 have to satisfy the requirements
64r 11 a11 ¼ 3a0 ðr 11 þ a11 Þ3 ;
8ðr 11 a0 Þr 11 a111 ¼ 3a0 ðr11 þ a111 Þ2 :
We have two requirements for the parameters a0 ; r0 ; a111 ; r11 . Thus we can choose a0 and r 0 as free parameters. Then a111 and r 11 are determined by the above two relationships. Thus for prescribed values of the growth rate and the interaction coefficient the kink can exist only for selected values of the adaptation parameters a111 and r 11 .
2954
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
(a)
1 0.8 0.6
Q 0.4 0.2 0 -0.2 -40
-20
0 ξ
20
40
(b)
0.2 0.1
Q*
0 -0.1 -0.2 -40
-20
0
20
40
ξ Fig. 1. (a) An example for solitary wave Q ðnÞ of kind (21). The wave parameters are as follows: a0 ¼ 1; a2 ¼ 2:519; p ¼ 0:316; q ¼ 1:265; 0:221; G ¼ 0:06; H ¼ 0:03, q0 ¼ 0. (b) Q ¼ dQ =dn, where Q ðnÞ is this one from (a).
v¼
3.2. Case E ¼ 0; F – 0; G – 0; H – 0 Let us discuss the case when the number of the members of the population does not influence the interaction between the members. In this case a111 ¼ 0 and from (14) we obtain the equation
@Q @2Q D 2 ¼ FQ 3 þ GQ 2 þ HQ: @t @x
ð23Þ
This equation is closely connected to the Burgers–Huxley equation which has many applications as well as interesting exact solutions (see for an example [67–71]). In our short analysis of some solutions of (23) we shall follow Kudryashov [66]. He proposed to look for exact solutions of nonlinear differential equations in the form
QðnÞ ¼
n X k¼0
aðkÞ/ðnÞk þ
l n X 1 d/ bl ; /ðnÞ dn l¼1
ð24Þ
where /ðnÞ is known solution of a nonlinear ordinary differential equation in the polynomial form and n can be determined by the first step of the Painleve test [72]. Below we shall consider the cases n ¼ 2 and n ¼ 3. Let first n ¼ 2. We search for solution of in the form (16) where ¼ 1 and aj and cj are parameters that we shall determine below. The balance of the powers leads to the relationship r ¼ 2n þ 2; n ¼ 2; 3; . . . If n ¼ 2 then r ¼ 6. Assuming that c0 ¼ c1 ¼ c3 ¼ c5 ¼ 0; c2 ¼ p2 ; c4 ¼ 2pq; c6 ¼ q2 we shall search for a solution in the form
QðnÞ ¼ a0 þ a1 / þ a2 /2 ;
d/ ¼ ðp/ þ q/3 Þ: dn
ð25Þ
In the method of Kudyashov (25) corresponds to bl ¼ 0; l ¼ 1; 2. Let ¼ 1. The solution of the corresponding system of nonlinear algebraic equations leads to a1 ¼ 0 and to the following relationships for the other wave parameters (below l ¼ 1)
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 G 2FD þ v F þ 6DpF pffiffiffiffiffiffiffiffiffiffiffiffiffi ; 3 F 2FD 1 ðG2 D þ v 2 F þ 6D2 p2 FÞ ; H¼ 3 FD pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 3=2 2 2 2FD 1 2GF Dv 2 2FDv F þ G p¼ : p ffiffiffiffiffiffiffiffiffiffiffiffi ffi 12 F 2 D GD þ v 2FD a2 ¼ 2l
pffiffiffiffiffiffiffiffiffiffiffiffiffi q 2FD ; F
2955
a0 ¼
ð26Þ
The solitary wave solution depends on the solution of the differential equation in (25). It is
Q ðnÞ ¼ a0 þ a2
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi p exp½2pðx þ cÞ 1 q exp½2pðx þ cÞ
ð27Þ
for the case ½/2 =ðp þ q/2 Þ > 0 and
Q ðnÞ ¼ a0 þ a2
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi p exp½2pðx þ cÞ 1 þ q exp½2pðx þ cÞ
ð28Þ
for the case ½/2 =ðp þ q/2 Þ < 0. An example for a kink of kind (27) is shown in Fig. 2. Let us now apply the method of Kudryashov and search for a solution of (23) of the form (24) with bl – 0; l ¼ 1; 2 and with /ðnÞ which is a solution of the Ricatti equation
d/ ¼ /2 þ 2a/ þ b; dn
ð29Þ
where a and b are parameters. After some algebraic calculations we obtain the solution
pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffipffiffiffiffiffiffiffiffiffiffi 2FD 1 d/ 1 7v 2FD 8GD 6v 2 2 FD Q ðnÞ ¼ /ðnÞ þ ; þ F 24 F /ðnÞ dn FD
ð30Þ
where /ðnÞ is a solution of the Ricatti equation with
pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 2GD þ 6a0 FD þ v 2FD pffiffiffiffiffiffiffiffiffiffiffiffiffi a¼ ; 6 D 2FD
pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 7v 2FD 8GD 6v ; a0 ¼ 24 FD
ð31Þ
1 2G2 D F v 2 þ 6FHD b¼ 48 FD2 The form of the solution of the Riccati equation is
/ðnÞ ¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ða2 þ bÞ tanh ða2 þ bÞn þ b þ a
þ
cþ
2 4aþ
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi expð4anÞ cosh ða2 þ bÞn þ b h i h i pffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 ffiffiffiffiffiffiffiffi 1 ffiffiffiffiffiffiffiffi p exp 4a þ a2 þ b n þ bÞ þ p exp 4a a2 þ b n þ bÞ 2 2 2 4a
a þb
ð32Þ
a þb
(a) 150
Q
100
50 -40
-20
0
20
ξ Fig. 2. (a) An example for solitary wave Q ðnÞ of kind (27). The wave parameters are as follows: a111 ¼ 0; D ¼ 0:1; G ¼ 0:1; r 0 ¼ 100; a0 ¼ 0:02; q0 ¼ 50; q ¼ 1000.
2956
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
and b and c are constants of integration. Let now n ¼ 3 in (24) and /ðnÞ is solution of the equation
d/ ¼ q/3 þ s/2 þ p/ þ q; dn
ð33Þ
where q; s; p; u are parameters. In this case the solution of (23) is
QðnÞ ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi 2q 2FD 2 2FDð2s2 D 18pqD 3qv Þ 4s 2FD þ /ðnÞ þ / ðnÞ: 18qFD 3F F
ð34Þ
In order to obtain solution of (33) we set /ðnÞ ¼ uðnÞ þ k where here k is a parameter. It is determined by the equation qk3 þ sk2 þ pk þ u ¼ 0 and we have to solve the following equation for uðnÞ
du ¼ qu3 þ ð3qk þ sÞu2 þ ð3qk2 þ 2sk þ pÞu ¼ a u3 þ b u2 þ c u: dn
ð35Þ
The implicit solution of this equation is
1
1=ð2c Þ 2
u C B
a u2 þb uþc C B C; n þ c ¼ ln B b
C B p ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi p ffiffiffiffiffiffiffiffiffiffiffiffiffi ffi
@ 2a uþb b 2 4a c 2c b 2 4a c A
p ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
2a uþb þ b 2 4a c 0
ð36Þ
where c is a constant of integration. The explicit solution can be obtained for selected values of the parameters. Let for an example b ¼ 0. Then for uðnÞ we have the solution
uðnÞ ¼
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi c expð2c ðn þ cÞÞ ; 1 a expð2c ðn þ cÞÞ
ðu2 =ðc þ a u2 ÞÞ > 0
ð37Þ
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi c expð2c ðn þ cÞÞ ; 1 þ a expð2c ðn þ cÞÞ
ðu2 =ðc þ a u2 ÞÞ < 0:
ð38Þ
and
uðnÞ ¼
Let us go back to the model Eq. (14) and let us set E ¼ 0, D ¼ 1; r0 ¼ a; 1 ¼ aa0 r11 as well as 1 þ a ¼ aðr11 a0 Þ. Eq. (14) is reduced to the equation
@Q @ 2 Q ¼ 2 Q 3 þ ð1 þ aÞQ 2 aQ; @t @x
0 0 and
Q ðnÞ ¼ a0 þ a2
2 p exp½pðx þ cÞ 1 þ q exp½pðx þ cÞ
ð48Þ
for the case ½/=ðp þ q/Þ < 0. Fig. 3a shows a kink of kind (44) and Fig. 3b shows a kink of kind (47) for the same values of the parameters in Fig. 3a. The model equation has other interesting solutions. One example for such a solution is presented in Fig. 3c. We call this solution KSW-wave (kink–solitary wave) because it has properties of a kink solution as well as the shape of a solitary wave. Now let us search for solution of Kudryashov kind (24) with b1 – 0 and b2 – 0 i.e. we search for solution of the kind
Q ðnÞ ¼ a0 þ a1 /ðnÞ þ a2 /2 ðnÞ þ
2 b1 d/ b2 d/ ; þ 2 /ðnÞ dn / ðnÞ dn
ð49Þ
where /ðnÞ is a solution of the Ricatti equation (29). After some algebraic manipulations we obtain for the coefficients of the solution (49)
a0 ¼
2 1 G2 b1 24
þ 51 v 2 12 HD 35 v b1 G þ 252 25
1 12v 5b1 G ; 5 G D a2 ¼ 6 ; G pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 6 4v þ 10v 2 25DH b1 ¼ ; 5 G D b2 ¼ 6 : G
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 25G2 b1 v 2 þ 540v 4 240v 3 b1 G þ 625H2 D2 DG
;
a1 ¼
The solution of the Riccati equation is again (32) but the coefficients a and b are
1 ða1 b1 GÞ ; 24 D 2 2 1 25G b1 þ 120v b1 G 216v 2 þ 900HD þ 1800a0 DG b¼ : 14400 D2 a¼
Finally let us note two equations connected to Eq. (14) with E ¼ F ¼ 0. First of all let a0 ¼ r 0 ¼ 1 and D ¼ 1. We obtain the equation
2958
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
70
(a)
60 50 40
Q
30 20 10 0 -10 -60
-40
-20
0
20
40
60
ξ 1
(b) 0.8 0.6
Q 0.4 0.2
0 -60
-40
-20
0
20
40
60
ξ 10
(c) 5
Q
0
-5
-10 -5
0
5
ξ Fig. 3. Several waves for the particular case of model Eq. (41) with quadratic nonlinearity. (a) Solution of kind (44) with parameters a0 ¼ 0:0, q ¼ 0:01; G ¼ 0:001; H ¼ 0:0003; rho0 ¼ 0; r 0 ¼ H; D ¼ 0:001. (b) Solution of kind (47) with parameters which are the same as these for (a) with exception of a which is 0. (c) Another solitary wave solution of kind (44). D ¼ 0:1; q ¼ 0:3; G ¼ 0:1; H ¼ 0:1.
@2Q @n2
þv
@Q þ Q ð1 Q Þ ¼ 0; @n
ð50Þ
i.e. the Fisher equation [76] much discussed since the pioneering work of Kolmogorov and coauthors [77]. Ablowitz and Zeppetela [78] first obtained exact analytical solution of (50) as follows:
1 QðnÞ ¼ h pffiffiffi i2 1 r expðn= 6Þ
ð51Þ
pffiffiffi with r < 0; n ¼ x v t; v ¼ 5= 6. From more recent research we mention the work of Malfliet [22] who has used the tanhmethod to solve this equation and obtained two solutions containing tanh functions. Other solutions of the Fisher equation
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
2959
are obtained by Wazwaz [73] who also used the tanh-method to obtain solutions which contain coth-function. Finally let D ¼ 1 and H ¼ 1. For the cases G ¼ 1; E ¼ F ¼ 0; F ¼ 1; E ¼ G ¼ 0; and E ¼ 1; F ¼ G ¼ 0 we obtain the generalized Fisher equation of the kind
@Q @Q ¼ 2 ¼ Q ð1 Q a Þ @t @x
ð52Þ
with a ¼ 1; a ¼ 2 or a ¼ 3, respectively. This particular case has been discussed by Wang [79] who obtained the solution
Q ðx; tÞ ¼
2=a
1 a aþ4 b 1 ; tanh pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi x pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi t þ þ 2 2 2 2a þ 4 2 2a þ 4
ð53Þ
where b is a constant. 4. Approximate solution of (15) and stability of the solutions of (14) Above we have obtained exact but particular solution of (15). Now we shall obtain approximate solution of (15) which will be valid when v =D 1, i.e. when the value of the wave velocity is much larger than the value of the diffusion coefficient. In order to understand the approximate solution we need a phase plane analysis of (15). Eq. (15) can be reduced to a system of two ODEs by introducing the variable P ¼ Q 0 .
P ¼ Q 0; P0 ¼
v D
P
ð54Þ
E 4 F 3 G 2 H Q Q Q Q: D D D D
We can investigate the dynamics in the phase plane ðQ ; Q 0 Þ. In the general case the fixed points are Q 1 ¼ 0 as well as the solutions of the equation EQ 3 þ FQ 2 þ GQ þ H ¼ 0. For this general case and for particular case of the population dynamics application of the model and when q0 – 0 we obtain large relationships for coordinates of the fixed points as well as for the eigenvalues connected to the linear stability of the fixed points. In order to illustrate the method of analysis we shall discuss the case q0 ¼ 0. The fixed points of (54) for this case are
Q 1 ¼ 0;
1 Q2 ¼ ; r 11
Q 3;4
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 1 a0 a0 þ 4a0 a111 ¼ : 2 a0 a111
The eigenvalues determining the linear stability properties of the fixed points are as follows:
1 h v pffiffiffiffiffiffii X1 ; 2 D h 1 v pffiffiffiffiffiffii ¼ X2 ; 2 D h 1 v pffiffiffiffiffiffii ¼ X3 ; 2 D
Q 1 : k1;2 ¼ Q 2 : k1;2 Q 3 : k1;2
v 2
r0 ; D D v 2 4r 0 4r 0 a ðr a Þ 0 11 111 X2 ¼ þ ; þ D D Dr211 X1 ¼
4
ð55Þ ð56Þ ð57Þ
where
X3 ¼
v 2
r 0 r 11 a0 r0 r11 r 0 r 11 2 8 þ2 D Da111 Da2111 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi r 0 a0 ða0 þ 4a111 Þ r0 2 þ8 ; Da111 D
Q 4 : k1;2 ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi a0 ða0 þ 4a111 Þ r 0 r 11 a0 ða0 þ 4a111 Þ r0 a0 þ 4 þ2 Da111 Da0 a111 Da2111
3=2
1 h v pffiffiffiffiffiffii X4 ; 2 D
ð58Þ
where
v 2
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
r 0 r 11 3=2 a0 ða0 þ 4a111 Þ r 0 r 11 a0 ða0 þ 4a111 Þ a0 r0 r11 r 0 r 11 r0 a0 X4 ¼ 2 8 2 4 þ2 2 2 D Da111 Da111 Da0 a111 Da111 Da111 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi r 0 a0 ða0 þ 4a111 Þ r0 þ2 þ8 : D Da111 The diffusion coefficient D must be positive and the wave velocity v must be different from 0. Let first X i ¼ 0 where i can be 1, 2, 3 or 4. Then the correspondent fixed point Q i is a degenerate node. If v > 0 the degenerate node is stable and if v < 0 the degenerate node is unstable. If X i is positive the eigenvalues k1;2 for the fixed point Q i are real. If in addition X i > ðv =DÞ2 Q i is a
2960
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
saddle point. If X i < ðv =DÞ2 Q i is a stable node. Finally if X i is negative Q i is a spiral point. For v > 0 the spiral point is stable and for v < 0 the spiral point is unstable. We are interested in finding solutions of the model equations which tend to reach constant values at n ! 1. In the phase plane ðQ ; Q 0 Þ this corresponds to a path connecting a stable and unstable fixed points. One possibility is a connection of saddle point and stable node. We shall investigate this possibility below. In order to obtain approximate solution of the model equation we shall use expansion in term of power series in appropriate small parameter. We have to solve (15) in presence of the boundary conditions Q ð1Þ ¼ l ¼ const; ðl – 0Þ and Q ð1Þ ¼ 0. Let us rewrite Eq. (15) as follows. We assume Q ¼ Q =l; v ¼ v =D; E ¼ El3 =D; F ¼ F l2 =D; G ¼ Gl=D; H ¼ H=D. (15) becomes
Q 00 þ v Q 0 þ E Q 4 þ F Q 3 þ G Q 2 þ H Q ¼ 0;
ð59Þ
where denotes d=dn. Let us introduce the parameter ¼ 1=v . When is small the dominant member in (59) is the term containing the first derivative of Q . Then for n ! 1 Q tends to constant values. Let x ¼ n1=2 and gðxÞ ¼ Q ðx=1=2 Þ. We expand gðxÞ in power series of 0
2
gðxÞ ¼ g 0 ðxÞ þ g 1 ðxÞ þ 2 g 2 ðxÞ þ . . .
ð60Þ
The substitution of (60) in (59) leads to following equations for
g 00 g 01
þ þ
E g 40 þ F g 30 þ G g 20 þ H g 0 ¼ 0; 4E g 30 g 1 þ 3F g 20 g 1 þ 2G g 0 g 1 þ
0
and
1
:
ð61Þ H g 1 þ g 000 ¼ 0:
ð62Þ
It is not possible to obtain general analytical solution of the above system of Eqs. (61) and (62). In order to illustrate the wave profile we shall obtain analytical solution for a particular case. Let q0 ¼ 0 and a111 ¼ r11 ¼ a0 . Then F ¼ G ¼ 0. E ¼ r 0 a30 and H ¼ r0 . The solution of (61) and (62) is
g 0 ðxÞ ¼
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi r0 3
; 3 expð3xr0 Þ r 0 a0 Z
Z
Z dxpðxÞ ; g 1 ðxÞ ¼ exp dxpðxÞ c þ dxqðxÞ exp
ð63Þ ð64Þ
where 2
pðxÞ ¼ 4
r 0 a0
3
expð3xr0 Þ r 0 a0
3
þ r0 ;
qðxÞ ¼ g 000 ðxÞ
and c is a constant of integration. Let c ¼ 0; r0 ¼ a0 ¼ 1. The approximate solution for g is
h i expð3xÞþ1 1 4 expð3xÞ ln expð3xÞ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi gðxÞ ¼ g 0 þ g 1 þ Oð Þ ¼ p þ Oð2 Þ 3 expð3xÞ þ 1 3 ½expð3xÞ þ 14=3 2
ð65Þ
is small when v =D 1 i.e. the value of the wave velocity must be much larger than the value of the diffusion coefficient. The graphics for g 0 ; g 1 and g ¼ g 0 þ g 1 are shown in Fig. 4. We mention that in the phase plane ðQ ; Q 0 Þ the solution (65) can correspond to connection between the stable node Q ¼ 0 and the saddle point Q ¼ 1=a0 . The conditions for such a connection are as follows. Q 1 is a stable node when r 0 < 0 and Q is a saddle point when r 211 þ a0 ðr 11 a111 Þ > 0. The above results give us another point of view on the solution (21). Let us set the constant of integration c ¼ 0 in (21). The solution connects the constant states Q ¼ 0 at n ! 1 and Q ¼ ðF=ð4EÞÞððq2 =ð1 q2 ÞÞ1=3 1Þ. The last state however is a fixed point only for selected values of the free parameter q. Let us illustrate this for the particular case of the population dynamics model and for q0 ¼ 0. In this case if we have to connect the fixed point Q ¼ 0 to the fixed point Q ¼ 1=r11 the value of q has to be sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1 4a111 =ðr11 þ a111 ÞÞ3 q¼ : 1 þ ð1 4a111 =ðr11 þ a111 ÞÞ3 Q 1 is a stable node when r0 < 0 and Q is a saddle point when r211 þ a0 ðr11 a111 Þ > 0. We note that the conditions for the connection stable node – saddle point for the exact solution are the same as the conditions for the connection stable node – saddle point for the approximate solution above. Finally let us investigate the asymptotic stability of the solutions of (14). Let us denote these solutions as UðnÞ. UðnÞ is asymptotically stable if a small perturbation imposed on the solution decays and the system returns to its original state described by the solution UðnÞ. We shall work in new coordinates: t and n. In these coordinates Eq. (14) becomes
@Q @Q @2Q v D 2 ¼ EQ 4 þ FQ 3 þ GQ 2 þ HQ : @t @n @n
ð66Þ
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
Fig. 4. Approximate solution of the model equation. (a) g 0 . (b) g 1 . (c) g ¼ g 0 þ g 1 ;
2961
¼ 0:01. r 0 ¼ a0 ¼ 1.
Now we search for perturbed solution in the form
Q ðn; tÞ ¼ UðnÞ þ Vðn; tÞ;
ð67Þ
where Vðn; tÞ is a small perturbation of the solution UðnÞ. In addition we assume that Vðn; tÞ ¼ 0 for jnj P L for some value of L. This means that we have a perturbation of the solution in a bounded interval and the perturbation vanishes out of this bounded interval. We assume that the perturbation V is small and because of this we can omit the nonlinear terms in the perturbation equation for V. Thus we obtain the linearized perturbation equation
@V @V @2V V D 2 ¼ Vð4EU 3 þ 3FU 2 þ 2GU þ HÞ: @t @n @n
ð68Þ
We search for solution of (68) in the form
Vðn; tÞ ¼ wðnÞ expðktÞ:
ð69Þ
2962
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
In addition we assume v ¼ v =D; k ¼ k=D; X ðUÞ ¼ XðUÞ=D ¼ ð4EU 3 þ 3FU 2 þ 2GU þ HÞ=D. The substitution of all this in (68) leads to the equation 2
d w dn2
þ v
dw þ wðk þ X ðUÞÞ ¼ 0: dn
ð70Þ
We eliminate the first derivative in (70) by means of the Liouville–Green transformation
vn wðnÞ ¼ hðnÞ exp 2
ð71Þ
and by means of the substitution X ðUÞ ¼ 1 XðUÞ Eq. (70) is reduced to the standard form 2
d h ðnÞh ¼ 0; þ ½k q dn
ð72Þ
ðnÞ ¼ XðUðnÞÞ þ v 2 =4 1. Now if q ðnÞ > 0 and q is continuous the eigenvalues k of the boundary wave problem (72) where q and the boundary conditions wðLÞ ¼ wðLÞ ¼ 0 on the interval L; n < L are all positive and the corresponding state is asymptotically stable. ðnÞ > 0. As we consider a problem from the Next we obtain the relationships which are consequence of the condition q ðnÞ > 0 population dynamics we expect that UðnÞ P 0 for all n. Let U max be the maximum value of U. Then the condition q means that there exists a critical velocity v c and waves with velocity smaller than this velocity are unstable. The stability condition for Eq. (14) is
v 2 > v 2c ¼ 4D
4EU 3max þ 3FU 2max þ 2GU max þ H :
ð73Þ
Let us write the stability conditions for the case when the model equation describes population density waves and q0 ¼ 0. For the case of solution (65) and for the parameters of the kink from Fig. 1
v 2 > v 2c ¼ 4D
1 4U 3max :
For the case of solution (21)
v 2 > v 2c ¼ 4Dr0
h i 1 4a0 r 11 a111 U 3max 3a0 ðr 11 þ a111 ÞU 2max þ 2ðr11 a0 ÞU max :
Finally we note that we can obtain approximate solutions and can investigate the stability of the solutions for all particular cases discussed in Section 3. The corresponding calculations are straightforward and analogous to the calculations from this section. Because of this we shall leave these calculations to the interested reader. 5. Concluding remarks In this paper we discuss a nonlinear model of the spatio-temporal dynamics of interacting populations. We consider in more details the particular case of 1D migration of 1 population. For this case the model system is reduced to a single (1 + 1)-dimensional nonlinear PDE which can be treated analytically. This equation describes the evolution of the density deviation from the states of stationary spatial density of the population members. If the amplitude of such a deviation becomes large solitary waves can travel through the system. We investigate such waves corresponding to connection of a saddle point and stable node in an appropriate phase plane. By means of power series expansion we obtain approximate solution for the case of large velocity/diffusion ratio. In addition by means of appropriate ansatz we obtain exact particular analytical solutions of the model equation and several of its particular cases. These solutions describe nonlinear kink and solitary waves. Finally we investigate the stability of the solutions and obtain conditions for asymptotic stability of the obtained kink wave solution of the model equation. The result is that there exists critical wave velocity and for velocities smaller than the critical velocity the population waves become unstable. The obtained solutions are examples of organized motion in a nonlinear population system. The motion of the population members is modeled by a diffusion process and one would not expect arising of collective motion. But because of the nonlinearity a self-organization can occur and then large number of population members begin to move synchronously. The result of this collective motion can be a stable nonlinear wave. Acknowledgment We would like thank to the ESF Action COST MP0801 ”Physics of Competition and Conflict” for the support of our research. References [1] W.A. Brock, D.A. Hsieh, B. LeBaron, Nonlinear Dynamics, Chaos and Instability, The MIT Press, Cambridge, MA, 1981. [2] J. Epstein, Nonlinear Dynamics Mathematical Biology and Social Science, Addison-Wesley, Reading, MA, 1997.
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964 [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52] [53] [54] [55] [56] [57] [58] [59] [60]
2963
L. Lam (Ed.), Introduction to Nonlinear Physics, Springer, New York, 1997. A. Oron, S.H. Davis, S.G. Bankoff, Long-scale evolution of thin liquid films, Rev. Mod. Phys. 69 (1997) 931–980. N.K. Vitanov, F.H. Busse, Bounds on the heat transport in a horizontal fluid layer with stress-free boundaries, ZAMP 48 (1997) 310–324. N.K. Vitanov, Upper bound on the heat transport in a layer of fluid of infinite Prandtl number, rigid lower boundary, and stress-free upper boundary, Phys. Rev. E 61 (2000) 956–959. H. Kantz, D. Holstein, M. Razwitz, N.K. Vitanov, Markov chain model for turbulent wind speed data, Physica A 342 (2004) 315–321. P. Ch. Ivanov, L.A.N. Amaral, A.L. Goldberger, S. Havlin, M.G. Rosenblum, S. Struzik, H.E. Stahley, Multifractality in human heartbeat dynamics, Nature 399 (1999) 461–465. T. Nagatani, The physics of traffic jams, Rep. Prog. Phys. 65 (2002) 1331–1386. A.C. Newell, Solitons in Mathematics and Physics, SIAM, Philadelphia, PA, 1985. R.K. Bullough, P.J. Caudrey (Eds.), Solitons, Springer, Berlin, 1980. M.J. Ablowitz, P.A. Clarkson, Solitons, Nonlinear Evolution Equations and Inverse Scattering, Cambridge University Press, Cambridge, 1991. P.G. Drazin, R.S. Johnson, Solitons: An Introduction, Cambridge University Press, Cambridge, 2002. T. Dauxois, M. Peyrard, Physics of Solitons, Cambridge University Press, Cambridge, 2006. D. Tilman, The ecological consequences of changes in biodiversity: a search for general principles, Ecology 80 (1999) 1455–1474. J. Henrich, R. McElveath, The evolution of cultural evolution, Evol. Anthropol. 12 (2003) 123–135. N.C. Stenseth, A. Mysterud, G. Ottersen, J.W. Hurrell, K.S. Chan, W. Lima, Ecological effects of climate fluctuations, Science 297 (2002) 1292–1296. P.P. Majumder, Ethnic populations of India as seen from evolutionary perspective, J. Biosci. 26 (2001) 533–545. M.L. Wang, Solitary wave solutions for variant Boussinesq equations, Phys. Lett. A 199 (1995) 169–172. Y.S. Kivsar, B. Luther-Davies, Dark optical solitons: physics and applications, Phys. Rep. 298 (1998) 81–197. K.S. Liu, Z.T. Fu, S.D. Liu, Q. Zhan, Jacobi elliptic function expansion method and periodic wave solutions of nonlinear wave equations, Phys. Lett. A 289 (2001) 69–74. W. Malfliet, Solitary wave solutions of nonlinear wave equations, Am. J. Phys. 60 (1992) 650–654. N.K. Vitanov, N.K. Martinov, On the solitary waves in the sine-Gordon model of the two-dimensional Josephson junction, Z. Phys. B 100 (1996) 129– 135. C. Dai, G. Zhou, J. Zhang, Exotic localized structures based on variable separation solution of (2 + 1)-dimensional KdV equation via the extended tanhmethod, Chaos Solitons Fract. 33 (2007) 1458–1467. C.-L. Zhang, G.-P. Cai, J.-Y. Qiang, Chaos solitons and fractals in (2 + 1)-dimensional KdV system derived from a periodic wave solution, Chaos Solitons Fract. 34 (2007) 1575–1583. C.-L. Zheng, H.-P. Zhu, Localized excitations with periodic and chaotic behaviors in (1 + 1)-dimensional Korteweg–de Vries type system, Chaos Solitons Fract. 34 (2007) 496–502. A.-M. Wazwaz, Multiple-soliton solutions for the generalized (1 + 1)-dimensional and the generalized (2 + 1)-dimensional Ito equations, Appl. Math. Comput. 202 (2008) 840–849. A.-M. Wazwaz, Multiple soliton solutions for the Boussinesq equation, Appl. Math. Comput. 192 (2007) 479–486. D. Korteweg, G. de Vries, On the change of the form of long waves advancing in a rectangular canal and on a new type of long stationary waves, Philos. Mag. 36 (5) (1895) 422–443. N.J. Zabusky, M.D. Kruskal, Interaction of ‘‘solitons” in a collisionless plasma and the recurrence of initial states, Phys. Rev. Lett. 15 (1965) 240–243. C.S. Gardner, J.M. Greene, M.D. Kruskal, R.M. Miura, Method for solving Korteweg–de Vries equation, Phys. Rev. Lett. 19 (1967) 1095–1097. T.R. Marchant, Asymptotic solutions for a third-order Korteveg-de-Vries equation, Chaos Solitons Fract. 22 (2004) 261–270. T.R. Marchant, Asymptotic solutions on a non-zero mean level, Chaos Solitons Fract. 32 (2007) 1328–1336. Y. Kuramoto, T. Tsuzuki, Persistent propagation of concentration waves in dissipative media far from thermal equilibrium, Prog. Theor. Phys. 55 (1976) 356–369. G.I. Sivashinsky, Large cells in nonlinear Marangoni convection, Physica D 4 (1982) 227–235. N.A. Kudryashov, Exact soliton solutions of the generalized evolution equation of wave dynamics, PMM USSR 52 (1988) 361–365. N.A. Kudryashov, Exact solutions of the generalized Kuramoto–Sivashinsky equation, Phys. Lett. A 147 (1990) 287–291. A.-M. Wazwaz, New solitary wave solutions to the Kuramoto–Sivashinsky and the Kawahara equations, Appl. Math. Comput. 182 (2006) 1642–1650. A.-M. Wazwaz, The tanh method: exact solution of the sine-Gordon and sinh-Gordon equations, Appl. Math. Comput. 167 (2005) 1196–1210. N. Martinov, N. Vitanov, On the correspondence between the self-consistent 2D Poisson–Boltzmann structures and the sine-Gordon waves, J. Phys. A: Math. Gen. 25 (1992) L51–L56. N.K. Martinov, N.K. Vitanov, On the self-consistent thermal equilibrium structures in two-dimensional negative-temperature systems, Canad. J. Phys. 72 (1994) 618–624. N. Martinov, N. Vitanov, On some solutions of the two-dimensional sine-Gordon equation, J. Phys. A: Math. Gen. 25 (1992) L419–L426. N.K. Martinov, N.K. Vitanov, New class of running-wave solutions of the 2 + 1-dimensional sine-Gordon equation, J. Phys. A: Math. Gen. 27 (1994) 4611–4618. N.K. Vitanov, On traveling waves and double-periodic structures in two-dimensional sine-Gordon systems, J. Phys. A: Math. Gen. 29 (1996) 5195– 5207. V. Druma, On an analytical solution of the two-dimensional Korteweg–de Vries equation, JETP Lett. 19 (12) (1974) 219–225. Z.-Y. Ma, X.-F. Wu, J.-M. Zhu, Multisoliton excitations for the Kadomtsev–Petviashvili equation and the coupled Burgers equation, Chaos Solitons Fract. 31 (2007) 648–657. C. Dai, F. Liu, J. Zhang, Novel types of interactions between solitons in the (2 + 1)-dimensional Nizhnik–Novikov–Veselov system, Chaos Solitons Fract. 36 (2008) 437–445. L.B. Cheng, New soliton structures with non-propagatory behavior in three-dimensional system, Chaos Solitons Fract. 36 (2008) 253–262. J. Zhu, Z. Ma, New exact solutions to the (2 + 1)-dimensional Broer–Kaup equation, Chaos Solitons Fract. 34 (2007) 476–481. W. Huang, Y. Liu, Z. Lu, Doubly periodic wave and folded solitary wave solutions for (2 + 1)-dimensional and higher order Broer–Kaup equation, Chaos Solitons Fract. 31 (2007) 54–63. C.-L. Zheng, H.-P. Zhu, L.-Q. Chen, Exact solution and semifolded structures of generalized Broer–Kaup system in (2 + 1)-dimensions, Chaos Solitons Fract. 26 (2005) 187–194. H.-Y. Ruan, Y.-H. Chen, Abundant soliton structures of (2 + 1)-dimensional NLS equation, Chaos Solitons Fract. 15 (2003) 947–955. S.B. Leble, N.V. Ustinov, On soliton and periodic solutions of Maxwell–Bloch system for two-level medium with degeneracy, Chaos Solitons Fract. 26 (2005) 187–194. A. Parker, Wave dynamics for peaked solitons of the Camassa–Holm equation, Chaos Solitons Fract. 35 (2008) 220–237. Z.I. Dimitrova, N.K. Vitanov, Influence of adaptation on the nonlinear dynamics of a system of competing populations, Phys. Lett. A 272 (2000) 368– 380. Z.I. Dimitrova, N.K. Vitanov, Dynamical consequences of adaptation of growth rates in a system of three competing populations, J. Phys. A: Math. Gen. 34 (2001) 7459–7473. Z.I. Dimitrova, N.K. Vitanov, Adaptation and its impact on the dynamics of a system of three competing populations, Physica A 300 (2001) 91–115. Z.I. Dimitrova, N.K. Vitanov, Chaotic pairwise competition, Theor. Popul. Biol. 66 (2004) 1–12. N.K. Vitanov, Z.I. Dimitrova, H. Kantz, On the trap of extinction and its elimination, Phys. Lett. A 349 (2006) 350–355. E.H. Kerner, Further considerations on the statistical mechanics of biological associations, Bull. Math. Biophys. 21 (1959) 217–253.
2964
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
[61] E. Fan, Multiple traveling wave solutions of nonlinear evolution equations using a unified algebraic method, J. Phys. A: Math. Gen. 35 (2002) 6853– 6872. [62] M.L. Wang, Solitary wave solutions for variant Boussinesq equations, Phys. Lett. A 199 (1999) 169–172. [63] L. Yang, J. Liu, K. Yang, Exact solutions of nonlinear PDE, nonlinear transformations and reduction of nonlinear PDE to a quadrature, Phys. Lett. A 278 (2001) 267–270. [64] E.J. Parker, B.R. Duffy, Travelling solitary wave solutions to a compound KdV-Burgers equation, Phys. Lett. A 229 (1997) 200–217. [65] N.A. Kudryashov, Simplest equation method to look for exact solutions of nonlinear differential equations, Chaos Solitons Fract. 24 (2005) 1217–1231. [66] N.A. Kudyashov, Exact solitary waves of the Fisher equation, Phys. Lett. A 342 (2005) 99–106. [67] F. Carrielo, M. Tabor, Painleve expansions for nonintegrable evolution-equations, Physica D 39 (1989) 77–94. [68] N.A. Kudryashov, Exact solutions of a family of Fisher equations, Theor. Math. Phys. 94 (1993) 211–218. [69] P. Broadbridge, B.H. Bradshaw, G.R. Fulford, G.K. Aldis, Huxley and Fisher equations for gene propagation, Anziam J. 44 (2002) 11–20. [70] O.Yu. Efimova, N.A. Kudyashov, Exact solutions of the Burgers–Huxley equation, J. Appl. Math. Mech. 68 (2004) 413–420. [71] N.A. Kudryashov, N.B. Loguinova, Extended simplest equation method for nonlinear differential equations, Appl. Math. Comput. 205 (2008) 396–402. [72] J. Weiss, M. Tabor, G. Carnevalle, The Painleve property for partial differential equations, J. Math. Phys. 24 (1983) 522–526. [73] A.-M. Wazwaz, The tanh method for traveling wave solutions of nonlinear equations, Appl. Math. Comput. 154 (2004) 713–723. [74] W.X. Ma, B. Fuchssteiner, Explicit and exact solutions to a Kolmogorov–Petrovskii–Piskunov equation, Int. J. Non-linear Dyn. 31 (1996) 329–338. [75] T. Kawahara, M. Tanaka, Interactions of travelling fronts: an exact solution of a nonlinear diffusion equation, Phys. Lett. A 97 (1983) 311–314. [76] R.A. Fisher, The wave of advance of advantageous gene, Ann. Eugenics 7 (1937) 355–369. [77] A. Kolmogorov, I. Petrovskii, N. Piskunov, A study of the diffusion equation with increasing in the amount of substance, and its application to a biological problem, Bull. Moscow Univ. Math. Mech. 1 (1937) 1–26. [78] M.J. Ablowitz, A. Zeppetela, Explicit solutions of Fisher equation for a special wave speed, Bull. Math. Biol. 41 (1979) 835–840. [79] X.Y. Wang, Exact and explicit wave solutions for the generalized Fisher equation, Phys. Lett. A 131 (1988) 277–279.
Contents lists available at ScienceDirect
Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc
On nonlinear population waves Nikolay K. Vitanov a,*, Ivan P. Jordanov a, Zlatinka I. Dimitrova b a b
Institute of Mechanics, Bulgarian Academy of Sciences, Akad. G. Bonchev Str., Bl. 4, 1113 Sofia, Bulgaria Institute of Solid State Physics, Bulgarian Academy of Sciences, Blvd. Tzarigradsko Chaussee 72, 1784 Sofia, Bulgaria
a r t i c l e Keywords: Nonlinear dynamics Population dynamics Migration Solitary waves Kinks Asymptotic stability
i n f o
a b s t r a c t We discuss a model system of partial differential equations for description of the spatiotemporal dynamics of interacting populations. We are interested in the waves caused by migration of the populations. We assume that the migration is a diffusion process influenced by the changing values of the birth rates and coefficients of interaction among the populations. For the particular case of one population and one spatial dimension the general model is reduced to analytically tractable PDE with polynomial nonlinearity up to 4th order. We investigate this particular case and obtain two kinds of solutions: (i) approximate solution for small value of the ratio between the coefficient of diffusion and the wave velocity and (ii) exact solutions which describe nonlinear kink and solitary waves. In an appropriate phase space the kinks correspond to a connection between two states represented by a saddle point and a stable node. Finally we derive conditions for the asymptotic stability of the obtained solutions. Ó 2009 Elsevier Inc. All rights reserved.
1. Introduction The dynamics of many systems in physics, chemistry, biology, fluid dynamics, economics or society is nonlinear [1–9]. If the nonlinearity compensates for the dispersion nonlinear traveling waves can propagate through the corresponding system [10–14]. A classical example for nonlinear dynamics is the dynamics of interacting populations [15,16]. One important aspect of this dynamics is the movement of humans or animals and this movement can lead to large solitary migration waves [17,18]. The solitary waves are very popular and frequently obtained in various areas of natural sciences [19–23]. In the last years many authors have obtained exact solutions which describe solitary waves, solitons, foldons, etc. in various equations and system of equations. Just several examples are the Korteweg–de Vries (KdV) system [24–26], Ito equations [27], Boussinesq equation [28], Korteweg–de Vries equation [29–33], Kuramoto–Sivashinsky equation and its generalizations [34–38], sine-Gordon and sinh-Gordon equations [39–44], Kadomtsev–Petviashvili equation [45,46], Nizhnik–Novikov–Veselov system [47,48], Broer–Kaup equation [49–51], nonlinear Schrödinger equation [52], Maxwell–Bloch system [53] or Camassa– Holm equation [54]. In this paper we shall concentrate our attention on the analytical investigation of one-dimensional migration waves caused by difference in the spatial density of the members of one population. The general model for n interacting populations on two-dimensional surface is described in Section 2. In Section 3 we discuss the case of one population and one spatial dimension and obtain exact particular solutions of the studied model equation and its particular cases. Section 4 is devoted to approximate solution of the model equation for the case of large values of the ratio between the wave velocity and diffusion coefficient. In this section we derive conditions for asymptotic stability of the obtained solutions. Several concluding remarks are summarized in Section 5. * Corresponding author. E-mail addresses: [email protected], [email protected] (N.K. Vitanov). 0096-3003/$ - see front matter Ó 2009 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2009.09.041
2951
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
2. Formulation of the problem We consider a two-dimensional plane populated by members of n populations. Let us select a small region Ds ¼ DxDy in this plane. In this region we observe DN i members from the ith population ði ¼ 1; 2; . . . ; nÞ. Then if Ds is sufficiently small we can define the density of the ith population
qi ð~r; tÞ ¼
DN i : Ds
ð1Þ
We assume that the members of the ith population can move through the borders of the small area under consideration. r; tÞ of the current of this movement. The net number of the members of the ith population, Then we can define the density~ ji ð~ ~ Þdl. Within the observed area the number of the members of the ith ~ is ð~ crossing a small line dl with normal vector m ji m population can change by births, deaths or interactions with other populations (one example is the interaction among predators and preys). These processes will be summarized mathematically by the function R below. The total change of the number of members from the ith population in some large area is a consequence of migration through the boundaries of the area and the birth, death and interaction processes within the area. As a mathematical expression this reads
@ @t
Z S
dsqi ð~ r; tÞ ¼
Z
~þ dl~ji m
Z
l
dsRi :
ð2Þ
S
The first term in (2) describes the net rate of increase of the density of the ith population, the second term describes the net rate of immigration into the area and the third term describes the net rate of increase exclusive of immigration. The differential form of Eq. (2) is
@ qi þ r ~ji ¼ Ri : @t
ð3Þ
Now we have to specify ~ ji and Ri . Ri can be taken in the general form analogous to the form proposed by Dimitrova and Vitanov [55–59]
" Ri ¼ r 0i qi 1
n X
#
a0ij qi þ Ai ;
ð4Þ
j¼1
where
Ai ¼
n X
r ij qj
j¼1
n X
a0ij ðaijk þ ril Þqi qj
j;l¼1
n X
a0ij rik aijl qj qk ql
ð5Þ
j;k:l¼1
and n is the number of the interacting populations. Eq. (4) is obtained as follows: (1.) Similar to the generalized Lotka–Volterra equations we assume that the changes of the density because of reproduction and interaction between the populations are
" Ri ¼ r i qi 1
n X
#
aij qj ;
j¼1
where ri is the growth ratio of the ith population and aij is the interaction coefficient measuring to what extent the growth of the ith population is influenced by the jth population. (2.) Growth ratios and interaction coefficients depend on the density of members of the populations
" ri ¼
r 0i
1þ
n X
# rik qk ;
k¼1
" 0 ij
aij ¼ a 1 þ
n X
#
aijl ql :
l¼1
We shall assume that ~ j has a form similar to the general form of the linear multicomponent diffusion [60]
~ji ¼
n X
Dik rqk ;
ð6Þ
k¼1
where Dik ¼ Dik ðqi ; qj ;~ r; tÞ. We note that in the population dynamics it is not necessary that Dik ¼ Dki . Below we shall treat Dik as constants. After the substitution of (4) and (6) in (3) we obtain the following model equations
" # n n n n X X X @ qi X 0 0 0 0 Dij Dqj ¼ r i qi 1 ðaij r ij Þqj aij ðaijk þ rik Þqj qk aij rik aijl qj qk ql ; @t j¼1 j¼1 j;k¼1 j;k;l¼1
D¼
@2 @2 þ 2: 2 @x @y
ð7Þ
The fixed points of the system (7) correspond to stationary densities q01 ; q02 ; . . . ; q0n . The stationary densities qi0 are qi0 ¼ 0 as well as the solutions of the system of equations
2952
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
1
n X
a0ij rij qj 0
j¼1
n X
a0ij ðaijk þ rik Þqj 0 qk0
j;k¼1
n X
a0ij rik aijl qj 0 qk0 ql0 ¼ 0:
ð8Þ
j;k;l¼1
Let us consider a deviation Q i of the density qi of the ith population from the corresponding stationary density q0i as follows:
qi ðx; y; tÞ ¼ Q i ðx; y; tÞ þ q0i ; i ¼ 1; 2; . . . ; n:
ð9Þ
The substitution of (9) in the right-hand side of (7) leads us to
" n n n X X @ qi X 0 Dij Dqj r i ðQ i þ q0i Þ a0ij rij Q j þ a0ij ðaijk þ rik ÞðQ j Q k þ Q j q0k þ Q k q0j Þ @t j¼1 j¼1 j;k¼1 þ
n X
#
a0ij rik aijl ðQ j Q k Q l þ Q j Q k q0l þ Q k Q l q0j þ Q l Q j q0k þ Q j q0k q0l þ Q k q0l q0j þ Q l q0j q0k Þ :
ð10Þ
j;k;l¼1
For the case of one population and one spatial dimension (7) reduces to the following model equation:
@ q1 DDq1 ¼ r0 q1 1 ða0 r 11 Þq1 a0 ða111 þ r 11 Þq21 a0 r 11 a111 q31 : @t
ð11Þ
The fixed points of (11) are
q10 ¼ 0; q10
1 ¼ ; r 11
q10
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 1 a0 a0 þ 4a0 a111 ¼ : 2 a0 a111
ð12Þ
We are interesting in the behavior of perturbations (linear and nonlinear) around the basic state q0 , i.e.
q1 ðx; tÞ ¼ q0 þ Q ðx; tÞ:
ð13Þ
The substitution of (13) in (11) leads to the following nonlinear equation for Q
@Q @2Q D 2 ¼ EQ 4 þ FQ 3 þ GQ 2 þ HQ; @t @x
ð14Þ
where
E ¼ r 0 a0 r 11 a111 ; 0
F ¼ r 0 a0 ða111 þ r11 þ 4q0 r 11 a111 Þ;
0
G ¼ r ½r11 a 3q0 a0 ða111 þ r11 þ 2r 11 a111 q0 Þ; H ¼ r 0 ½1 q0 ð3a0 r 11 q0 þ 4a0 r 11 a111 q20 2r11 þ 2a0 þ 3a0 a111 q0 Þ: Below we shall treat (14) as general form of a model equation which for the particular case when E; F; G; H satisfy the above relationships describes population density waves. In the following section we shall obtain and discuss: (i) traveling wave solutions of (14) and (ii) traveling wave solutions of important equations connected to (14). 3. Traveling wave solutions of (14) and its particular cases 3.1. Case E – 0; F – 0; G – 0; H – 0 We consider Eq. (14) and search for solutions of the kind of traveling waves: Q ðx; tÞ ¼ Q ðnÞ ¼ Q ðx v tÞ where velocity of the wave. We have to solve the equation 2
D
d Q dn2
þv
dQ þ EQ 4 þ FQ 3 þ GQ 2 þ HQ ¼ 0: dn
v is the ð15Þ
Following [61] we assume that Q ðnÞ has the form
QðnÞ ¼
n X i¼0
ai /i ;
vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi uX u r d/ cj /j ; ¼ t dn j¼0
ð16Þ
where ¼ 1, and aj and cj are parameters that we shall determine below. (For more results and methods for obtaining exact solutions of nonlinear PDE see for an example [62–65].) We note that the methodology based on (16) is closely connected to the methodology of Kudryashov [65,66] which will be much used below. Substituting (16) in (15) we obtain an equation that contains powers of /. We have to balance the highest power arising from the second derivative in (15) with the highest power arising in the term containing Q 4 in the same equation. The resulting balance equation is
r ¼ 3n þ 2;
n ¼ 2; 3; . . .
ð17Þ
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
2953
If n ¼ 2 then r ¼ 8. Assuming that c0 ¼ c1 ¼ c3 ¼ c4 ¼ c6 ¼ c7 ¼ 0, c2 ¼ p2 ; c5 ¼ 2pq; c8 ¼ q2 we search for a solution of the form
Q ðnÞ ¼ a0 þ a1 / þ a2 /2 ;
d/ ¼ ðp/ þ q/4 Þ: dn
ð18Þ
Substituting this in (15) we obtain the following system of nine equations for the parameters a0 ; a1 ; a2 ; p; q; D; E; F; G; H
a2 10Dq2 þ Ea32 ¼ 0; a1 Dq2 þ Ea32 ¼ 0; a22 Fa2 þ 4Ea2 a0 þ 6Ea21 ¼ 0; a2 2v q þ 14Dpq þ 3Fa1 a2 þ 12Ea1 a2 a0 þ 4Ea31 ¼ 0; 5Da1 pq þ Ga22 þ v a1 q þ 3Fa22 a0 þ 3Fa21 a2 þ 6Ea20 a22 þ 12Ea2 a0 a21 þ Ea41 ¼ 0;
v;
ð19Þ
a1 ð2Ga2 þ 6Fa0 a2 þ Fa21 þ 12Ea20 a2 þ 4Ea0 a21 Þ ¼ 0; 3Fa20 a2 þ 3Fa0 a21 þ 2Ga0 a2 þ Ga21 þ 4Ea30 a2 þ 6Ea20 a21 þ Ha2 þ 4Da2 p2 þ 2v a2 p ¼ 0; a1 ð3Fa20 þ H þ 2Ga0 þ 4Ea30 þ Dp2 þ v pÞ ¼ 0; Ga0 þ H þ Fa20 þ Ea30 ¼ 0: Eq. (19) implies that a1 ¼ 0. Let in addition
¼ 1. The solution of the system is
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 10Dq2 E2
F a0 ¼ ; a2 ¼ 4Epffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 7 F 10DF v¼ ; 80 E
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 F 10DF ; p¼ ; E 80 DE 2 3 3F 3F G¼ ; H¼ 8E 64E2
ð20Þ
and q; F; E; D are free parameters. The expression for the solitary wave depends on the solution of the differential equation in (18) and it is given by
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 3 p exp½3pðn þ cÞ Q ðnÞ ¼ a0 þ a2 1 q exp½3pðn þ cÞ
ð21Þ
for the case ½/3 =ðp þ q/3 Þ > 0 and
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 3 p exp½3pðn þ cÞ Q ðnÞ ¼ a0 þ a2 1 þ q exp½3pðn þ cÞ
ð22Þ
for the case ½/3 =ðp þ q/3 Þ < 0. Above c is a constant of integration. Fig. 1 shows a kink of kind (21) as well as its derivative which is a well formed solitary wave. For the discussed case of modeling of the movement of population members and for q0 ¼ 0 we obtain
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 10Dq2 a2 ¼ ; r 0 a0 r 11 a111 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 r 11 þ a111 p¼ 10Dr0 a0 ðr11 þ a111 Þ; 80 Dr 11 a111 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 7 a111 þ r 11 v¼ 10Dr 0 a0 ða111 þ r 11 Þ; 80 a111 r 11 r11 þ a111 a0 ¼ ; 4r 11 a111
where the parameters r 11 ; a0 and a111 have to satisfy the requirements
64r 11 a11 ¼ 3a0 ðr 11 þ a11 Þ3 ;
8ðr 11 a0 Þr 11 a111 ¼ 3a0 ðr11 þ a111 Þ2 :
We have two requirements for the parameters a0 ; r0 ; a111 ; r11 . Thus we can choose a0 and r 0 as free parameters. Then a111 and r 11 are determined by the above two relationships. Thus for prescribed values of the growth rate and the interaction coefficient the kink can exist only for selected values of the adaptation parameters a111 and r 11 .
2954
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
(a)
1 0.8 0.6
Q 0.4 0.2 0 -0.2 -40
-20
0 ξ
20
40
(b)
0.2 0.1
Q*
0 -0.1 -0.2 -40
-20
0
20
40
ξ Fig. 1. (a) An example for solitary wave Q ðnÞ of kind (21). The wave parameters are as follows: a0 ¼ 1; a2 ¼ 2:519; p ¼ 0:316; q ¼ 1:265; 0:221; G ¼ 0:06; H ¼ 0:03, q0 ¼ 0. (b) Q ¼ dQ =dn, where Q ðnÞ is this one from (a).
v¼
3.2. Case E ¼ 0; F – 0; G – 0; H – 0 Let us discuss the case when the number of the members of the population does not influence the interaction between the members. In this case a111 ¼ 0 and from (14) we obtain the equation
@Q @2Q D 2 ¼ FQ 3 þ GQ 2 þ HQ: @t @x
ð23Þ
This equation is closely connected to the Burgers–Huxley equation which has many applications as well as interesting exact solutions (see for an example [67–71]). In our short analysis of some solutions of (23) we shall follow Kudryashov [66]. He proposed to look for exact solutions of nonlinear differential equations in the form
QðnÞ ¼
n X k¼0
aðkÞ/ðnÞk þ
l n X 1 d/ bl ; /ðnÞ dn l¼1
ð24Þ
where /ðnÞ is known solution of a nonlinear ordinary differential equation in the polynomial form and n can be determined by the first step of the Painleve test [72]. Below we shall consider the cases n ¼ 2 and n ¼ 3. Let first n ¼ 2. We search for solution of in the form (16) where ¼ 1 and aj and cj are parameters that we shall determine below. The balance of the powers leads to the relationship r ¼ 2n þ 2; n ¼ 2; 3; . . . If n ¼ 2 then r ¼ 6. Assuming that c0 ¼ c1 ¼ c3 ¼ c5 ¼ 0; c2 ¼ p2 ; c4 ¼ 2pq; c6 ¼ q2 we shall search for a solution in the form
QðnÞ ¼ a0 þ a1 / þ a2 /2 ;
d/ ¼ ðp/ þ q/3 Þ: dn
ð25Þ
In the method of Kudyashov (25) corresponds to bl ¼ 0; l ¼ 1; 2. Let ¼ 1. The solution of the corresponding system of nonlinear algebraic equations leads to a1 ¼ 0 and to the following relationships for the other wave parameters (below l ¼ 1)
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 G 2FD þ v F þ 6DpF pffiffiffiffiffiffiffiffiffiffiffiffiffi ; 3 F 2FD 1 ðG2 D þ v 2 F þ 6D2 p2 FÞ ; H¼ 3 FD pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 3=2 2 2 2FD 1 2GF Dv 2 2FDv F þ G p¼ : p ffiffiffiffiffiffiffiffiffiffiffiffi ffi 12 F 2 D GD þ v 2FD a2 ¼ 2l
pffiffiffiffiffiffiffiffiffiffiffiffiffi q 2FD ; F
2955
a0 ¼
ð26Þ
The solitary wave solution depends on the solution of the differential equation in (25). It is
Q ðnÞ ¼ a0 þ a2
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi p exp½2pðx þ cÞ 1 q exp½2pðx þ cÞ
ð27Þ
for the case ½/2 =ðp þ q/2 Þ > 0 and
Q ðnÞ ¼ a0 þ a2
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi p exp½2pðx þ cÞ 1 þ q exp½2pðx þ cÞ
ð28Þ
for the case ½/2 =ðp þ q/2 Þ < 0. An example for a kink of kind (27) is shown in Fig. 2. Let us now apply the method of Kudryashov and search for a solution of (23) of the form (24) with bl – 0; l ¼ 1; 2 and with /ðnÞ which is a solution of the Ricatti equation
d/ ¼ /2 þ 2a/ þ b; dn
ð29Þ
where a and b are parameters. After some algebraic calculations we obtain the solution
pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffipffiffiffiffiffiffiffiffiffiffi 2FD 1 d/ 1 7v 2FD 8GD 6v 2 2 FD Q ðnÞ ¼ /ðnÞ þ ; þ F 24 F /ðnÞ dn FD
ð30Þ
where /ðnÞ is a solution of the Ricatti equation with
pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 2GD þ 6a0 FD þ v 2FD pffiffiffiffiffiffiffiffiffiffiffiffiffi a¼ ; 6 D 2FD
pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 7v 2FD 8GD 6v ; a0 ¼ 24 FD
ð31Þ
1 2G2 D F v 2 þ 6FHD b¼ 48 FD2 The form of the solution of the Riccati equation is
/ðnÞ ¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ða2 þ bÞ tanh ða2 þ bÞn þ b þ a
þ
cþ
2 4aþ
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi expð4anÞ cosh ða2 þ bÞn þ b h i h i pffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 ffiffiffiffiffiffiffiffi 1 ffiffiffiffiffiffiffiffi p exp 4a þ a2 þ b n þ bÞ þ p exp 4a a2 þ b n þ bÞ 2 2 2 4a
a þb
ð32Þ
a þb
(a) 150
Q
100
50 -40
-20
0
20
ξ Fig. 2. (a) An example for solitary wave Q ðnÞ of kind (27). The wave parameters are as follows: a111 ¼ 0; D ¼ 0:1; G ¼ 0:1; r 0 ¼ 100; a0 ¼ 0:02; q0 ¼ 50; q ¼ 1000.
2956
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
and b and c are constants of integration. Let now n ¼ 3 in (24) and /ðnÞ is solution of the equation
d/ ¼ q/3 þ s/2 þ p/ þ q; dn
ð33Þ
where q; s; p; u are parameters. In this case the solution of (23) is
QðnÞ ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi 2q 2FD 2 2FDð2s2 D 18pqD 3qv Þ 4s 2FD þ /ðnÞ þ / ðnÞ: 18qFD 3F F
ð34Þ
In order to obtain solution of (33) we set /ðnÞ ¼ uðnÞ þ k where here k is a parameter. It is determined by the equation qk3 þ sk2 þ pk þ u ¼ 0 and we have to solve the following equation for uðnÞ
du ¼ qu3 þ ð3qk þ sÞu2 þ ð3qk2 þ 2sk þ pÞu ¼ a u3 þ b u2 þ c u: dn
ð35Þ
The implicit solution of this equation is
1
1=ð2c Þ 2
u C B
a u2 þb uþc C B C; n þ c ¼ ln B b
C B p ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi p ffiffiffiffiffiffiffiffiffiffiffiffiffi ffi
@ 2a uþb b 2 4a c 2c b 2 4a c A
p ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
2a uþb þ b 2 4a c 0
ð36Þ
where c is a constant of integration. The explicit solution can be obtained for selected values of the parameters. Let for an example b ¼ 0. Then for uðnÞ we have the solution
uðnÞ ¼
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi c expð2c ðn þ cÞÞ ; 1 a expð2c ðn þ cÞÞ
ðu2 =ðc þ a u2 ÞÞ > 0
ð37Þ
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi c expð2c ðn þ cÞÞ ; 1 þ a expð2c ðn þ cÞÞ
ðu2 =ðc þ a u2 ÞÞ < 0:
ð38Þ
and
uðnÞ ¼
Let us go back to the model Eq. (14) and let us set E ¼ 0, D ¼ 1; r0 ¼ a; 1 ¼ aa0 r11 as well as 1 þ a ¼ aðr11 a0 Þ. Eq. (14) is reduced to the equation
@Q @ 2 Q ¼ 2 Q 3 þ ð1 þ aÞQ 2 aQ; @t @x
0 0 and
Q ðnÞ ¼ a0 þ a2
2 p exp½pðx þ cÞ 1 þ q exp½pðx þ cÞ
ð48Þ
for the case ½/=ðp þ q/Þ < 0. Fig. 3a shows a kink of kind (44) and Fig. 3b shows a kink of kind (47) for the same values of the parameters in Fig. 3a. The model equation has other interesting solutions. One example for such a solution is presented in Fig. 3c. We call this solution KSW-wave (kink–solitary wave) because it has properties of a kink solution as well as the shape of a solitary wave. Now let us search for solution of Kudryashov kind (24) with b1 – 0 and b2 – 0 i.e. we search for solution of the kind
Q ðnÞ ¼ a0 þ a1 /ðnÞ þ a2 /2 ðnÞ þ
2 b1 d/ b2 d/ ; þ 2 /ðnÞ dn / ðnÞ dn
ð49Þ
where /ðnÞ is a solution of the Ricatti equation (29). After some algebraic manipulations we obtain for the coefficients of the solution (49)
a0 ¼
2 1 G2 b1 24
þ 51 v 2 12 HD 35 v b1 G þ 252 25
1 12v 5b1 G ; 5 G D a2 ¼ 6 ; G pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 6 4v þ 10v 2 25DH b1 ¼ ; 5 G D b2 ¼ 6 : G
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 25G2 b1 v 2 þ 540v 4 240v 3 b1 G þ 625H2 D2 DG
;
a1 ¼
The solution of the Riccati equation is again (32) but the coefficients a and b are
1 ða1 b1 GÞ ; 24 D 2 2 1 25G b1 þ 120v b1 G 216v 2 þ 900HD þ 1800a0 DG b¼ : 14400 D2 a¼
Finally let us note two equations connected to Eq. (14) with E ¼ F ¼ 0. First of all let a0 ¼ r 0 ¼ 1 and D ¼ 1. We obtain the equation
2958
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
70
(a)
60 50 40
Q
30 20 10 0 -10 -60
-40
-20
0
20
40
60
ξ 1
(b) 0.8 0.6
Q 0.4 0.2
0 -60
-40
-20
0
20
40
60
ξ 10
(c) 5
Q
0
-5
-10 -5
0
5
ξ Fig. 3. Several waves for the particular case of model Eq. (41) with quadratic nonlinearity. (a) Solution of kind (44) with parameters a0 ¼ 0:0, q ¼ 0:01; G ¼ 0:001; H ¼ 0:0003; rho0 ¼ 0; r 0 ¼ H; D ¼ 0:001. (b) Solution of kind (47) with parameters which are the same as these for (a) with exception of a which is 0. (c) Another solitary wave solution of kind (44). D ¼ 0:1; q ¼ 0:3; G ¼ 0:1; H ¼ 0:1.
@2Q @n2
þv
@Q þ Q ð1 Q Þ ¼ 0; @n
ð50Þ
i.e. the Fisher equation [76] much discussed since the pioneering work of Kolmogorov and coauthors [77]. Ablowitz and Zeppetela [78] first obtained exact analytical solution of (50) as follows:
1 QðnÞ ¼ h pffiffiffi i2 1 r expðn= 6Þ
ð51Þ
pffiffiffi with r < 0; n ¼ x v t; v ¼ 5= 6. From more recent research we mention the work of Malfliet [22] who has used the tanhmethod to solve this equation and obtained two solutions containing tanh functions. Other solutions of the Fisher equation
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
2959
are obtained by Wazwaz [73] who also used the tanh-method to obtain solutions which contain coth-function. Finally let D ¼ 1 and H ¼ 1. For the cases G ¼ 1; E ¼ F ¼ 0; F ¼ 1; E ¼ G ¼ 0; and E ¼ 1; F ¼ G ¼ 0 we obtain the generalized Fisher equation of the kind
@Q @Q ¼ 2 ¼ Q ð1 Q a Þ @t @x
ð52Þ
with a ¼ 1; a ¼ 2 or a ¼ 3, respectively. This particular case has been discussed by Wang [79] who obtained the solution
Q ðx; tÞ ¼
2=a
1 a aþ4 b 1 ; tanh pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi x pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi t þ þ 2 2 2 2a þ 4 2 2a þ 4
ð53Þ
where b is a constant. 4. Approximate solution of (15) and stability of the solutions of (14) Above we have obtained exact but particular solution of (15). Now we shall obtain approximate solution of (15) which will be valid when v =D 1, i.e. when the value of the wave velocity is much larger than the value of the diffusion coefficient. In order to understand the approximate solution we need a phase plane analysis of (15). Eq. (15) can be reduced to a system of two ODEs by introducing the variable P ¼ Q 0 .
P ¼ Q 0; P0 ¼
v D
P
ð54Þ
E 4 F 3 G 2 H Q Q Q Q: D D D D
We can investigate the dynamics in the phase plane ðQ ; Q 0 Þ. In the general case the fixed points are Q 1 ¼ 0 as well as the solutions of the equation EQ 3 þ FQ 2 þ GQ þ H ¼ 0. For this general case and for particular case of the population dynamics application of the model and when q0 – 0 we obtain large relationships for coordinates of the fixed points as well as for the eigenvalues connected to the linear stability of the fixed points. In order to illustrate the method of analysis we shall discuss the case q0 ¼ 0. The fixed points of (54) for this case are
Q 1 ¼ 0;
1 Q2 ¼ ; r 11
Q 3;4
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 1 a0 a0 þ 4a0 a111 ¼ : 2 a0 a111
The eigenvalues determining the linear stability properties of the fixed points are as follows:
1 h v pffiffiffiffiffiffii X1 ; 2 D h 1 v pffiffiffiffiffiffii ¼ X2 ; 2 D h 1 v pffiffiffiffiffiffii ¼ X3 ; 2 D
Q 1 : k1;2 ¼ Q 2 : k1;2 Q 3 : k1;2
v 2
r0 ; D D v 2 4r 0 4r 0 a ðr a Þ 0 11 111 X2 ¼ þ ; þ D D Dr211 X1 ¼
4
ð55Þ ð56Þ ð57Þ
where
X3 ¼
v 2
r 0 r 11 a0 r0 r11 r 0 r 11 2 8 þ2 D Da111 Da2111 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi r 0 a0 ða0 þ 4a111 Þ r0 2 þ8 ; Da111 D
Q 4 : k1;2 ¼
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi a0 ða0 þ 4a111 Þ r 0 r 11 a0 ða0 þ 4a111 Þ r0 a0 þ 4 þ2 Da111 Da0 a111 Da2111
3=2
1 h v pffiffiffiffiffiffii X4 ; 2 D
ð58Þ
where
v 2
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
r 0 r 11 3=2 a0 ða0 þ 4a111 Þ r 0 r 11 a0 ða0 þ 4a111 Þ a0 r0 r11 r 0 r 11 r0 a0 X4 ¼ 2 8 2 4 þ2 2 2 D Da111 Da111 Da0 a111 Da111 Da111 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi r 0 a0 ða0 þ 4a111 Þ r0 þ2 þ8 : D Da111 The diffusion coefficient D must be positive and the wave velocity v must be different from 0. Let first X i ¼ 0 where i can be 1, 2, 3 or 4. Then the correspondent fixed point Q i is a degenerate node. If v > 0 the degenerate node is stable and if v < 0 the degenerate node is unstable. If X i is positive the eigenvalues k1;2 for the fixed point Q i are real. If in addition X i > ðv =DÞ2 Q i is a
2960
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
saddle point. If X i < ðv =DÞ2 Q i is a stable node. Finally if X i is negative Q i is a spiral point. For v > 0 the spiral point is stable and for v < 0 the spiral point is unstable. We are interested in finding solutions of the model equations which tend to reach constant values at n ! 1. In the phase plane ðQ ; Q 0 Þ this corresponds to a path connecting a stable and unstable fixed points. One possibility is a connection of saddle point and stable node. We shall investigate this possibility below. In order to obtain approximate solution of the model equation we shall use expansion in term of power series in appropriate small parameter. We have to solve (15) in presence of the boundary conditions Q ð1Þ ¼ l ¼ const; ðl – 0Þ and Q ð1Þ ¼ 0. Let us rewrite Eq. (15) as follows. We assume Q ¼ Q =l; v ¼ v =D; E ¼ El3 =D; F ¼ F l2 =D; G ¼ Gl=D; H ¼ H=D. (15) becomes
Q 00 þ v Q 0 þ E Q 4 þ F Q 3 þ G Q 2 þ H Q ¼ 0;
ð59Þ
where denotes d=dn. Let us introduce the parameter ¼ 1=v . When is small the dominant member in (59) is the term containing the first derivative of Q . Then for n ! 1 Q tends to constant values. Let x ¼ n1=2 and gðxÞ ¼ Q ðx=1=2 Þ. We expand gðxÞ in power series of 0
2
gðxÞ ¼ g 0 ðxÞ þ g 1 ðxÞ þ 2 g 2 ðxÞ þ . . .
ð60Þ
The substitution of (60) in (59) leads to following equations for
g 00 g 01
þ þ
E g 40 þ F g 30 þ G g 20 þ H g 0 ¼ 0; 4E g 30 g 1 þ 3F g 20 g 1 þ 2G g 0 g 1 þ
0
and
1
:
ð61Þ H g 1 þ g 000 ¼ 0:
ð62Þ
It is not possible to obtain general analytical solution of the above system of Eqs. (61) and (62). In order to illustrate the wave profile we shall obtain analytical solution for a particular case. Let q0 ¼ 0 and a111 ¼ r11 ¼ a0 . Then F ¼ G ¼ 0. E ¼ r 0 a30 and H ¼ r0 . The solution of (61) and (62) is
g 0 ðxÞ ¼
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi r0 3
; 3 expð3xr0 Þ r 0 a0 Z
Z
Z dxpðxÞ ; g 1 ðxÞ ¼ exp dxpðxÞ c þ dxqðxÞ exp
ð63Þ ð64Þ
where 2
pðxÞ ¼ 4
r 0 a0
3
expð3xr0 Þ r 0 a0
3
þ r0 ;
qðxÞ ¼ g 000 ðxÞ
and c is a constant of integration. Let c ¼ 0; r0 ¼ a0 ¼ 1. The approximate solution for g is
h i expð3xÞþ1 1 4 expð3xÞ ln expð3xÞ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi gðxÞ ¼ g 0 þ g 1 þ Oð Þ ¼ p þ Oð2 Þ 3 expð3xÞ þ 1 3 ½expð3xÞ þ 14=3 2
ð65Þ
is small when v =D 1 i.e. the value of the wave velocity must be much larger than the value of the diffusion coefficient. The graphics for g 0 ; g 1 and g ¼ g 0 þ g 1 are shown in Fig. 4. We mention that in the phase plane ðQ ; Q 0 Þ the solution (65) can correspond to connection between the stable node Q ¼ 0 and the saddle point Q ¼ 1=a0 . The conditions for such a connection are as follows. Q 1 is a stable node when r 0 < 0 and Q is a saddle point when r 211 þ a0 ðr 11 a111 Þ > 0. The above results give us another point of view on the solution (21). Let us set the constant of integration c ¼ 0 in (21). The solution connects the constant states Q ¼ 0 at n ! 1 and Q ¼ ðF=ð4EÞÞððq2 =ð1 q2 ÞÞ1=3 1Þ. The last state however is a fixed point only for selected values of the free parameter q. Let us illustrate this for the particular case of the population dynamics model and for q0 ¼ 0. In this case if we have to connect the fixed point Q ¼ 0 to the fixed point Q ¼ 1=r11 the value of q has to be sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1 4a111 =ðr11 þ a111 ÞÞ3 q¼ : 1 þ ð1 4a111 =ðr11 þ a111 ÞÞ3 Q 1 is a stable node when r0 < 0 and Q is a saddle point when r211 þ a0 ðr11 a111 Þ > 0. We note that the conditions for the connection stable node – saddle point for the exact solution are the same as the conditions for the connection stable node – saddle point for the approximate solution above. Finally let us investigate the asymptotic stability of the solutions of (14). Let us denote these solutions as UðnÞ. UðnÞ is asymptotically stable if a small perturbation imposed on the solution decays and the system returns to its original state described by the solution UðnÞ. We shall work in new coordinates: t and n. In these coordinates Eq. (14) becomes
@Q @Q @2Q v D 2 ¼ EQ 4 þ FQ 3 þ GQ 2 þ HQ : @t @n @n
ð66Þ
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
Fig. 4. Approximate solution of the model equation. (a) g 0 . (b) g 1 . (c) g ¼ g 0 þ g 1 ;
2961
¼ 0:01. r 0 ¼ a0 ¼ 1.
Now we search for perturbed solution in the form
Q ðn; tÞ ¼ UðnÞ þ Vðn; tÞ;
ð67Þ
where Vðn; tÞ is a small perturbation of the solution UðnÞ. In addition we assume that Vðn; tÞ ¼ 0 for jnj P L for some value of L. This means that we have a perturbation of the solution in a bounded interval and the perturbation vanishes out of this bounded interval. We assume that the perturbation V is small and because of this we can omit the nonlinear terms in the perturbation equation for V. Thus we obtain the linearized perturbation equation
@V @V @2V V D 2 ¼ Vð4EU 3 þ 3FU 2 þ 2GU þ HÞ: @t @n @n
ð68Þ
We search for solution of (68) in the form
Vðn; tÞ ¼ wðnÞ expðktÞ:
ð69Þ
2962
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
In addition we assume v ¼ v =D; k ¼ k=D; X ðUÞ ¼ XðUÞ=D ¼ ð4EU 3 þ 3FU 2 þ 2GU þ HÞ=D. The substitution of all this in (68) leads to the equation 2
d w dn2
þ v
dw þ wðk þ X ðUÞÞ ¼ 0: dn
ð70Þ
We eliminate the first derivative in (70) by means of the Liouville–Green transformation
vn wðnÞ ¼ hðnÞ exp 2
ð71Þ
and by means of the substitution X ðUÞ ¼ 1 XðUÞ Eq. (70) is reduced to the standard form 2
d h ðnÞh ¼ 0; þ ½k q dn
ð72Þ
ðnÞ ¼ XðUðnÞÞ þ v 2 =4 1. Now if q ðnÞ > 0 and q is continuous the eigenvalues k of the boundary wave problem (72) where q and the boundary conditions wðLÞ ¼ wðLÞ ¼ 0 on the interval L; n < L are all positive and the corresponding state is asymptotically stable. ðnÞ > 0. As we consider a problem from the Next we obtain the relationships which are consequence of the condition q ðnÞ > 0 population dynamics we expect that UðnÞ P 0 for all n. Let U max be the maximum value of U. Then the condition q means that there exists a critical velocity v c and waves with velocity smaller than this velocity are unstable. The stability condition for Eq. (14) is
v 2 > v 2c ¼ 4D
4EU 3max þ 3FU 2max þ 2GU max þ H :
ð73Þ
Let us write the stability conditions for the case when the model equation describes population density waves and q0 ¼ 0. For the case of solution (65) and for the parameters of the kink from Fig. 1
v 2 > v 2c ¼ 4D
1 4U 3max :
For the case of solution (21)
v 2 > v 2c ¼ 4Dr0
h i 1 4a0 r 11 a111 U 3max 3a0 ðr 11 þ a111 ÞU 2max þ 2ðr11 a0 ÞU max :
Finally we note that we can obtain approximate solutions and can investigate the stability of the solutions for all particular cases discussed in Section 3. The corresponding calculations are straightforward and analogous to the calculations from this section. Because of this we shall leave these calculations to the interested reader. 5. Concluding remarks In this paper we discuss a nonlinear model of the spatio-temporal dynamics of interacting populations. We consider in more details the particular case of 1D migration of 1 population. For this case the model system is reduced to a single (1 + 1)-dimensional nonlinear PDE which can be treated analytically. This equation describes the evolution of the density deviation from the states of stationary spatial density of the population members. If the amplitude of such a deviation becomes large solitary waves can travel through the system. We investigate such waves corresponding to connection of a saddle point and stable node in an appropriate phase plane. By means of power series expansion we obtain approximate solution for the case of large velocity/diffusion ratio. In addition by means of appropriate ansatz we obtain exact particular analytical solutions of the model equation and several of its particular cases. These solutions describe nonlinear kink and solitary waves. Finally we investigate the stability of the solutions and obtain conditions for asymptotic stability of the obtained kink wave solution of the model equation. The result is that there exists critical wave velocity and for velocities smaller than the critical velocity the population waves become unstable. The obtained solutions are examples of organized motion in a nonlinear population system. The motion of the population members is modeled by a diffusion process and one would not expect arising of collective motion. But because of the nonlinearity a self-organization can occur and then large number of population members begin to move synchronously. The result of this collective motion can be a stable nonlinear wave. Acknowledgment We would like thank to the ESF Action COST MP0801 ”Physics of Competition and Conflict” for the support of our research. References [1] W.A. Brock, D.A. Hsieh, B. LeBaron, Nonlinear Dynamics, Chaos and Instability, The MIT Press, Cambridge, MA, 1981. [2] J. Epstein, Nonlinear Dynamics Mathematical Biology and Social Science, Addison-Wesley, Reading, MA, 1997.
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964 [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52] [53] [54] [55] [56] [57] [58] [59] [60]
2963
L. Lam (Ed.), Introduction to Nonlinear Physics, Springer, New York, 1997. A. Oron, S.H. Davis, S.G. Bankoff, Long-scale evolution of thin liquid films, Rev. Mod. Phys. 69 (1997) 931–980. N.K. Vitanov, F.H. Busse, Bounds on the heat transport in a horizontal fluid layer with stress-free boundaries, ZAMP 48 (1997) 310–324. N.K. Vitanov, Upper bound on the heat transport in a layer of fluid of infinite Prandtl number, rigid lower boundary, and stress-free upper boundary, Phys. Rev. E 61 (2000) 956–959. H. Kantz, D. Holstein, M. Razwitz, N.K. Vitanov, Markov chain model for turbulent wind speed data, Physica A 342 (2004) 315–321. P. Ch. Ivanov, L.A.N. Amaral, A.L. Goldberger, S. Havlin, M.G. Rosenblum, S. Struzik, H.E. Stahley, Multifractality in human heartbeat dynamics, Nature 399 (1999) 461–465. T. Nagatani, The physics of traffic jams, Rep. Prog. Phys. 65 (2002) 1331–1386. A.C. Newell, Solitons in Mathematics and Physics, SIAM, Philadelphia, PA, 1985. R.K. Bullough, P.J. Caudrey (Eds.), Solitons, Springer, Berlin, 1980. M.J. Ablowitz, P.A. Clarkson, Solitons, Nonlinear Evolution Equations and Inverse Scattering, Cambridge University Press, Cambridge, 1991. P.G. Drazin, R.S. Johnson, Solitons: An Introduction, Cambridge University Press, Cambridge, 2002. T. Dauxois, M. Peyrard, Physics of Solitons, Cambridge University Press, Cambridge, 2006. D. Tilman, The ecological consequences of changes in biodiversity: a search for general principles, Ecology 80 (1999) 1455–1474. J. Henrich, R. McElveath, The evolution of cultural evolution, Evol. Anthropol. 12 (2003) 123–135. N.C. Stenseth, A. Mysterud, G. Ottersen, J.W. Hurrell, K.S. Chan, W. Lima, Ecological effects of climate fluctuations, Science 297 (2002) 1292–1296. P.P. Majumder, Ethnic populations of India as seen from evolutionary perspective, J. Biosci. 26 (2001) 533–545. M.L. Wang, Solitary wave solutions for variant Boussinesq equations, Phys. Lett. A 199 (1995) 169–172. Y.S. Kivsar, B. Luther-Davies, Dark optical solitons: physics and applications, Phys. Rep. 298 (1998) 81–197. K.S. Liu, Z.T. Fu, S.D. Liu, Q. Zhan, Jacobi elliptic function expansion method and periodic wave solutions of nonlinear wave equations, Phys. Lett. A 289 (2001) 69–74. W. Malfliet, Solitary wave solutions of nonlinear wave equations, Am. J. Phys. 60 (1992) 650–654. N.K. Vitanov, N.K. Martinov, On the solitary waves in the sine-Gordon model of the two-dimensional Josephson junction, Z. Phys. B 100 (1996) 129– 135. C. Dai, G. Zhou, J. Zhang, Exotic localized structures based on variable separation solution of (2 + 1)-dimensional KdV equation via the extended tanhmethod, Chaos Solitons Fract. 33 (2007) 1458–1467. C.-L. Zhang, G.-P. Cai, J.-Y. Qiang, Chaos solitons and fractals in (2 + 1)-dimensional KdV system derived from a periodic wave solution, Chaos Solitons Fract. 34 (2007) 1575–1583. C.-L. Zheng, H.-P. Zhu, Localized excitations with periodic and chaotic behaviors in (1 + 1)-dimensional Korteweg–de Vries type system, Chaos Solitons Fract. 34 (2007) 496–502. A.-M. Wazwaz, Multiple-soliton solutions for the generalized (1 + 1)-dimensional and the generalized (2 + 1)-dimensional Ito equations, Appl. Math. Comput. 202 (2008) 840–849. A.-M. Wazwaz, Multiple soliton solutions for the Boussinesq equation, Appl. Math. Comput. 192 (2007) 479–486. D. Korteweg, G. de Vries, On the change of the form of long waves advancing in a rectangular canal and on a new type of long stationary waves, Philos. Mag. 36 (5) (1895) 422–443. N.J. Zabusky, M.D. Kruskal, Interaction of ‘‘solitons” in a collisionless plasma and the recurrence of initial states, Phys. Rev. Lett. 15 (1965) 240–243. C.S. Gardner, J.M. Greene, M.D. Kruskal, R.M. Miura, Method for solving Korteweg–de Vries equation, Phys. Rev. Lett. 19 (1967) 1095–1097. T.R. Marchant, Asymptotic solutions for a third-order Korteveg-de-Vries equation, Chaos Solitons Fract. 22 (2004) 261–270. T.R. Marchant, Asymptotic solutions on a non-zero mean level, Chaos Solitons Fract. 32 (2007) 1328–1336. Y. Kuramoto, T. Tsuzuki, Persistent propagation of concentration waves in dissipative media far from thermal equilibrium, Prog. Theor. Phys. 55 (1976) 356–369. G.I. Sivashinsky, Large cells in nonlinear Marangoni convection, Physica D 4 (1982) 227–235. N.A. Kudryashov, Exact soliton solutions of the generalized evolution equation of wave dynamics, PMM USSR 52 (1988) 361–365. N.A. Kudryashov, Exact solutions of the generalized Kuramoto–Sivashinsky equation, Phys. Lett. A 147 (1990) 287–291. A.-M. Wazwaz, New solitary wave solutions to the Kuramoto–Sivashinsky and the Kawahara equations, Appl. Math. Comput. 182 (2006) 1642–1650. A.-M. Wazwaz, The tanh method: exact solution of the sine-Gordon and sinh-Gordon equations, Appl. Math. Comput. 167 (2005) 1196–1210. N. Martinov, N. Vitanov, On the correspondence between the self-consistent 2D Poisson–Boltzmann structures and the sine-Gordon waves, J. Phys. A: Math. Gen. 25 (1992) L51–L56. N.K. Martinov, N.K. Vitanov, On the self-consistent thermal equilibrium structures in two-dimensional negative-temperature systems, Canad. J. Phys. 72 (1994) 618–624. N. Martinov, N. Vitanov, On some solutions of the two-dimensional sine-Gordon equation, J. Phys. A: Math. Gen. 25 (1992) L419–L426. N.K. Martinov, N.K. Vitanov, New class of running-wave solutions of the 2 + 1-dimensional sine-Gordon equation, J. Phys. A: Math. Gen. 27 (1994) 4611–4618. N.K. Vitanov, On traveling waves and double-periodic structures in two-dimensional sine-Gordon systems, J. Phys. A: Math. Gen. 29 (1996) 5195– 5207. V. Druma, On an analytical solution of the two-dimensional Korteweg–de Vries equation, JETP Lett. 19 (12) (1974) 219–225. Z.-Y. Ma, X.-F. Wu, J.-M. Zhu, Multisoliton excitations for the Kadomtsev–Petviashvili equation and the coupled Burgers equation, Chaos Solitons Fract. 31 (2007) 648–657. C. Dai, F. Liu, J. Zhang, Novel types of interactions between solitons in the (2 + 1)-dimensional Nizhnik–Novikov–Veselov system, Chaos Solitons Fract. 36 (2008) 437–445. L.B. Cheng, New soliton structures with non-propagatory behavior in three-dimensional system, Chaos Solitons Fract. 36 (2008) 253–262. J. Zhu, Z. Ma, New exact solutions to the (2 + 1)-dimensional Broer–Kaup equation, Chaos Solitons Fract. 34 (2007) 476–481. W. Huang, Y. Liu, Z. Lu, Doubly periodic wave and folded solitary wave solutions for (2 + 1)-dimensional and higher order Broer–Kaup equation, Chaos Solitons Fract. 31 (2007) 54–63. C.-L. Zheng, H.-P. Zhu, L.-Q. Chen, Exact solution and semifolded structures of generalized Broer–Kaup system in (2 + 1)-dimensions, Chaos Solitons Fract. 26 (2005) 187–194. H.-Y. Ruan, Y.-H. Chen, Abundant soliton structures of (2 + 1)-dimensional NLS equation, Chaos Solitons Fract. 15 (2003) 947–955. S.B. Leble, N.V. Ustinov, On soliton and periodic solutions of Maxwell–Bloch system for two-level medium with degeneracy, Chaos Solitons Fract. 26 (2005) 187–194. A. Parker, Wave dynamics for peaked solitons of the Camassa–Holm equation, Chaos Solitons Fract. 35 (2008) 220–237. Z.I. Dimitrova, N.K. Vitanov, Influence of adaptation on the nonlinear dynamics of a system of competing populations, Phys. Lett. A 272 (2000) 368– 380. Z.I. Dimitrova, N.K. Vitanov, Dynamical consequences of adaptation of growth rates in a system of three competing populations, J. Phys. A: Math. Gen. 34 (2001) 7459–7473. Z.I. Dimitrova, N.K. Vitanov, Adaptation and its impact on the dynamics of a system of three competing populations, Physica A 300 (2001) 91–115. Z.I. Dimitrova, N.K. Vitanov, Chaotic pairwise competition, Theor. Popul. Biol. 66 (2004) 1–12. N.K. Vitanov, Z.I. Dimitrova, H. Kantz, On the trap of extinction and its elimination, Phys. Lett. A 349 (2006) 350–355. E.H. Kerner, Further considerations on the statistical mechanics of biological associations, Bull. Math. Biophys. 21 (1959) 217–253.
2964
N.K. Vitanov et al. / Applied Mathematics and Computation 215 (2009) 2950–2964
[61] E. Fan, Multiple traveling wave solutions of nonlinear evolution equations using a unified algebraic method, J. Phys. A: Math. Gen. 35 (2002) 6853– 6872. [62] M.L. Wang, Solitary wave solutions for variant Boussinesq equations, Phys. Lett. A 199 (1999) 169–172. [63] L. Yang, J. Liu, K. Yang, Exact solutions of nonlinear PDE, nonlinear transformations and reduction of nonlinear PDE to a quadrature, Phys. Lett. A 278 (2001) 267–270. [64] E.J. Parker, B.R. Duffy, Travelling solitary wave solutions to a compound KdV-Burgers equation, Phys. Lett. A 229 (1997) 200–217. [65] N.A. Kudryashov, Simplest equation method to look for exact solutions of nonlinear differential equations, Chaos Solitons Fract. 24 (2005) 1217–1231. [66] N.A. Kudyashov, Exact solitary waves of the Fisher equation, Phys. Lett. A 342 (2005) 99–106. [67] F. Carrielo, M. Tabor, Painleve expansions for nonintegrable evolution-equations, Physica D 39 (1989) 77–94. [68] N.A. Kudryashov, Exact solutions of a family of Fisher equations, Theor. Math. Phys. 94 (1993) 211–218. [69] P. Broadbridge, B.H. Bradshaw, G.R. Fulford, G.K. Aldis, Huxley and Fisher equations for gene propagation, Anziam J. 44 (2002) 11–20. [70] O.Yu. Efimova, N.A. Kudyashov, Exact solutions of the Burgers–Huxley equation, J. Appl. Math. Mech. 68 (2004) 413–420. [71] N.A. Kudryashov, N.B. Loguinova, Extended simplest equation method for nonlinear differential equations, Appl. Math. Comput. 205 (2008) 396–402. [72] J. Weiss, M. Tabor, G. Carnevalle, The Painleve property for partial differential equations, J. Math. Phys. 24 (1983) 522–526. [73] A.-M. Wazwaz, The tanh method for traveling wave solutions of nonlinear equations, Appl. Math. Comput. 154 (2004) 713–723. [74] W.X. Ma, B. Fuchssteiner, Explicit and exact solutions to a Kolmogorov–Petrovskii–Piskunov equation, Int. J. Non-linear Dyn. 31 (1996) 329–338. [75] T. Kawahara, M. Tanaka, Interactions of travelling fronts: an exact solution of a nonlinear diffusion equation, Phys. Lett. A 97 (1983) 311–314. [76] R.A. Fisher, The wave of advance of advantageous gene, Ann. Eugenics 7 (1937) 355–369. [77] A. Kolmogorov, I. Petrovskii, N. Piskunov, A study of the diffusion equation with increasing in the amount of substance, and its application to a biological problem, Bull. Moscow Univ. Math. Mech. 1 (1937) 1–26. [78] M.J. Ablowitz, A. Zeppetela, Explicit solutions of Fisher equation for a special wave speed, Bull. Math. Biol. 41 (1979) 835–840. [79] X.Y. Wang, Exact and explicit wave solutions for the generalized Fisher equation, Phys. Lett. A 131 (1988) 277–279.
