New classes of stable exact solutions for a nonlinear rotational DNA model

Vol. 41 (2010) ACTA PHYSICA POLONICA B No 9 NEW CLASSES OF STABLE EXACT SOLUTIONS FOR A NONLINEAR ROTATIONAL DNA MODEL∗ J. Ramírez, J.L. Romero Departamento de Matemáticas, Universidad de Cádiz P.O. Box 40, 11510 Puerto Real, Cádiz, Spain F.R. Romero†, A. Álvarez, J.F.R. Archilla Grupo de Física No Lineal, Facultad de Física, Universidad de Sevilla Avda. Reina Mercedes, s/n, 41012 Sevilla, Spain (Received April 16, 2010) We consider a system of two coupled nonlinear partial differential equations for describing the rotational motions of bases in both polynucleotide chains of the DNA molecule. The model was proposed by L.V. Yakushevich and it is well known that the model supports, for some operating regimes, traveling wave solutions as kink–(antikink) soliton solutions. We have tried to make some progress by performing an analysis of the classical symmetries of this model. Our study shows that the model does not have enough symmetries as to reduce the equations to ordinary differential equations. Nevertheless, the known symmetries have been useful for finding several classes of exact solutions, by imposing adequate Ansätze. Some of them are kink–(antikink) like solutions, but other ones are not traveling wave solutions. For some of the new solutions, we have carried out a qualitative study and analyzed some stability properties. We think that they could be significant for the description of the DNA molecule as well as for some other applications. PACS numbers: 02.30.Jr, 02.20.–a, 87.14.gk, 05.45.Yv 1. Introduction Nonlinear physics of DNA is a part of nonlinear science which has attained an enormous development since the beginning of the last eighties. It is generally accepted that it started in 1980 with the article by Englander ∗We acknowledge financial support from the DGYCYT project MTM2006-05031, and Junta de Andalucía FQM 201. P06-FQM-01448. †[email protected] (2001) 2002 J. Ramírez et al. et al. [1] entitled Nature of the Open State in Long Polynucleotide Double Helices: Possibility of Soliton Excitations, where they proposed the first nonlinear Hamiltonian model of DNA. The DNA molecule is a discrete system consisting of many atoms having a quasi-one-dimensional structure. There exist different kinds of interactions between the main atomic groups. One of them is the stacking interaction between neighboring bases along the DNA axis, these are short-range forces which stabilize the DNA structure and hold one base over the next one forming a stack of bases. The DNA can be considered as a complex dynamical system because DNA interacts with its environment in the cell, as many kind of proteins and thermal fluctuations can impose different functioning regimes. A first approximation to a theoretical description of DNA begins by ignoring all the external influences and assuming that the amplitudes of the internal motions in DNA are small. This gives rise to the linear theory of DNA which is based on linear mathematical models. But without this last assumption, necessarily, the mathematical models are nonlinear, the mathematical difficulties increase and much more interesting phenomena appear. There are many experimental data and theoretical results published about the nonlinear properties of DNA (for a review see, e.g., [2]). Many authors have contributed to the nonlinear theory of DNA by proposing improved Hamiltonian models and also working in the study of possible solutions to the dynamical equations. Among them, it is worth remarking the Peyrard– Bishop model [3] introduced for the study of DNA thermal denaturation. This model describes the motions of the displacements of the bases from their equilibrium positions along the direction of the hydrogen bonds that connect the two bases in a pair. Other approaches describe the rotational motions of bases in both polynucleotide chains of the DNA molecule. The Yakushevich model (Y-model) [4] belongs to this type of models and in a certain sense can be considered as an improved version of the Englander model, it is the continuous limit of a discrete model describing two chain of disks connected with one another by lengthwise and transverse springs. It consists in a system of two coupled nonlinear partial differential equations and at present nobody has succeeded in finding the general solution of this nonlinear system. However, after some assumptions the model can be reduced to the well known sine-Gordon equation or to the double sine-Gordon equation. In this context, kink soliton solutions have been found by Newton’s method and by Hereman’s method [5]. The sine-Gordon equation also appears in many branches of physics, where the known solutions have been applied, although to our knowledge some of them have not been related to the physics of DNA. New Classes of Stable Exact Solutions for a Nonlinear Rotational . .. 2003 A method for finding some families of particular solutions is based on the introduction of some specific Ansätze that yield to a reduction of the system to a a well known partial differential equation or to ordinary differential equations. Sometimes, the Ansätze are suggested by particular functional regimes of the DNA molecule, as when it is imposed to some external constrains. Fortunately, there exists a method for obtaining a diversity of Ansätze which is based on the invariance properties of the system under Lie groups of transformations [7,8]. Following this procedure we have been able to find new classes of particular solutions some of which are not traveling waves and we have analyzed some stability properties. This paper is organized as follows: Section 2 introduces the Y-model, and we perform a symmetry group analysis of this model obtaining the classical symmetries. In Section 3 we show that it is possible to obtain two classes of similarity solutions. One of them is constituted by traveling wave solutions, some of them can be obtained explicitly, but other ones are implicitly determined by integrals. The other class of solutions are not traveling waves solutions and we have carried out a qualitative study and analyzed some stability properties. Finally, a summary and some conclusions are presented in Section 4. 2. Symmetry properties of the Yakushevich model The Y-model is a continuous model describing the dynamics of the DNA molecule, it takes into account the rotational motions of bases in both polynucleotide chains of DNA. As stressed by Yakushevich, it should be seen as the first step in a hierarchy of more and more realistic models, but its study deserves attention just due to its simplicity, opening the possibility of discovering some essential features of DNA torsion dynamics without unessential complications due to more realistic modeling. The Y-model is defined by the following two coupled nonlinear partial differential equations: α(ϕ1)tt =β(ϕ1)zz +δ(2 sin(ϕ1)−sin(ϕ1+ϕ2)) , α(ϕ2)tt =γ(ϕ2)zz +δ(2 sin(ϕ2)−sin(ϕ1+ϕ2)) ,(1) where ϕ1(z, t)and ϕ2(z, t)represent the angles of rotation at the point z and time tof the polynucleotide chains 1 and 2, respectively, and α, β, δ, γ are constant parameters. System (1) is rather involved and only special classes of solutions are known. Most of the known solutions have been obtained after some assumptions are made, which are suggested by special regimes of the DNA dynamics, and by introducing some simple specific Ansätze. In most cases, only some specific classes of solutions of the reduced equations have been considered as applied to DNA. 2004 J. Ramírez et al. It is expected that the most probable conformational distortions in DNA resulting in an open state configuration are obtained through an evolution with ϕ1=−ϕ2=ϕor with ϕ1=ϕ2=ϕ. Both assumptions make possible the reduction of the two-field model (1) to a one-field equation. The first assumption gives rise to the well known sine-Gordon equation and has been considered by many authors (see [2]). The second assumption gives rise to a double sine-Gordon equation. Nevertheless, to our knowledge nobody has connected sufficiently the richness of the solutions of both equations with the physics of DNA, and in particular with this rotational model. If, for instance, we suppose that ϕ1=−ϕ2=ϕ, system (1) reduces to (γ−β)ϕzz = 0 , 2αϕtt −(β+γ)ϕzz −4δsin(ϕ)=0.(2) This equation describes the torsion dynamics when the two DNA strands move out of phase by π. If in addition we suppose β=γ, system (2) reduces to ϕtt =aϕzz +bsin(ϕ),(3) where a=β αand b= 2 δ α. Equation (3) is the sine-Gordon equation and, in the context of the Y-model, kink–(antikink) soliton solutions have been found following Newton’s Method and Hereman’s Method (see [2]). This equation also appears in relativistic field theory, condensed matter physics and in a number of other physical applications, as dislocations in crystals or the propagation of fluxons in Josephson junctions. In this section we investigate the symmetry properties of system (1) under Lie groups of transformations [7,8]. The knowledge of the invariance properties usually permits some reductions of the system and to obtain the so-called similarity solutions. Since, nowadays, this is a standard method in the study of ordinary and partial differential equations, we omit the intermediate details and directly give the relevant reductions that the symmetry groups generate. We denote by Z, T, U, V the infinitesimals corresponding to the variables z, t, ϕ1, ϕ2, respectively. The classical symmetries of system (1) can be obtained imposing the invariance conditions. The procedure leads to a complicated system of equations, from which we obtain the corresponding Lie group. Two cases appear, depending of the values of the parameters β and γ: 1. In the general case (β6=γ), we obtain a 2-parameter Lie group given by Z=c1, T =c2, U = 0 , V = 0 . New Classes of Stable Exact Solutions for a Nonlinear Rotational . . . 2005 The infinitesimal generators of this group are v1=∂ ∂z ,v2=∂ ∂t .(4) 2. In the special case of β=γwe obtain a 3-parameter Lie group given by Z=c3β t +c1, T =c3α z +c2, U = 0 , V = 0 .(5) The generators of this group are v1,v2and v3=β t ∂ ∂z +α z ∂ ∂t . The commutator table is v1v2v3 v10 0 αv2 v20 0 βv1 v3−αv2−βv10 (6) As we can observe, for both cases the model lacks nontrivial symmetries. Nevertheless, the known symmetries are useful for finding new classes of stable exact solutions. We dedicate the next section to obtain the similarity solutions. Some of them are also kink–(antikink) soliton solutions, which are obtainable either explicitly or implicitly, and other ones are not traveling waves. 3. Similarity solutions The known symmetries make it possible to obtain two classes of solutions, characterized by a different type of dependence with respect to the variables z, t. 3.1. Traveling wave solutions The above results inform that, for both the cases β6=γand β=γ, system (1) is invariant under the translations of v1and v2. This property leads us to look for traveling wave solutions. These solutions are of the form ϕ1=f1(z+λt), ϕ2=f2(z+λt).(7) Substituting, it can be seen that f1and f2must satisfy the following system of ordinary differential equations a1f00 1−2 sin(f1) + sin(f1+f2)=0, a2f00 2−2 sin(f2) + sin(f1+f2)=0,(8) where a1=α λ2−β δ;a2=α λ2−γ δand the derivatives are with respect to the similarity variable x=z+λt. 2006 J. Ramírez et al. For obtaining similarity solutions of system (8) we have to determine its classical symmetries. However, we have investigated that system (8) lacks nontrivial symmetries. Therefore, Lie symmetries are of no help to reduce system (8) and the solutions, we are able to find, are consequence of some specific Ansätze. Let us suppose that f2=f1+ 2kπ, being k∈Z, and a1=a2; in this case system (8) reduces to the equation a1f00 1−2 sin(f1) + sin(2f1) = 0 .(9) With the appropriate change of the dependent variable, this equation can be written in a rational form or in a non transcendent form. •Rational form: It can be achieved by making a change of the dependent variable in the form f1= 2parctan(h), where pis a non null integer number. For each different value of p, a different ordinary differential equation for his obtained. In the following, we consider the cases corresponding to p= 1,2,3. For any positive integer p, the corresponding possible values of the angle ϕ1are in the interval (−pπ, pπ). 1. Taking p= 1, the change of the dependent variable is f1= 2 arctan(h), equation (9) becomes a11 + h2h00 −2h2h2+a1h02= 0 .(10) We have found several families of solutions of this equation that determine the following solutions of (9) in terms of the Jacobi elliptic functions [9]: f1=±2 arctan r−1−m 2sn 2 m−1rm+1 a1 (x+c1), m!, f1=±2 arctan r−1−m 2msn−12 m−1rm+1 a1 (x+c1), m!, f1=±2 arctan s1−2m 2(m−1)cn 2r1−2m a1 (x+c1), m!, f1=±2 arctan s2−m 2(m−1)dn 2√m−2 m√a1 (x+c1), m!, f1=±2 arctan rm−2 2dn−12√m−2 m√a1 (x+c1), m!, f1=±2 arctan rm−1 2sd 2r1−2m a1 (x+c1), m!.(11) New Classes of Stable Exact Solutions for a Nonlinear Rotational . . . 2007 All the functions appearing in these solutions are bounded. As an example, we represent in Fig. 1 the profile and the graphic of the first traveling wave solution in (11), for the parameters a1=−1,m=c1= 0. -10 -5 5 10 -3 -2 -1 1 2 3 f x -4 0 4 -4 0 4 -2 0 2 t z Fig. 1. Representation of the profile and the graphic of the first solution in (11) that corresponds to m= 0,c1= 0 and a1=−1. The ranges of the dependent variables are x=z+ 2t∈[−5,5] (left), z∈[−5,5] and t∈[−5,4] (right). 2. Taking p= 2, the change of the dependent variable is f1= 4 arctan(h), equation (9) becomes a1h00 1 + h23−2a1hh021 + h22+16h3−1 + h2= 0 .(12) This equation admits a family of solutions that is implicitly given by x+c2=± h Zξ2+ 1 q8 a1ξ4+c1(1 + ξ2)4dξ , (13) where c1and c2are arbitrary constants. The existence domains of the functions h, implicitly determined by equation (13), are given by the convergence of the integral in the right-hand side of (13). A solution with domain (−∞,∞) requires the divergence of the integral. For c16= 0 the integrals are convergent, as improper integrals, at ±∞. If the radicand have only complex roots, the solution his defined in a bounded interval. If the radicand have only simple real roots ai, it changes the sign at aiand the integrand behaves as (ξ−ai)−1/2when ξ→a+ i; in this case the improper integral is convergent and 2008 J. Ramírez et al. the solution is only defined in a finite interval. Therefore, to obtain solutions defined for all x, the radicand must have real roots of multiplicity greater or equal than 2; this happens when the radicand and its derivative have common roots. Two cases appear: A. When c1= 0, the right member of (13) is an improper integral at h= 0 and h= +∞. In this case it is easy to calculate the integral and we get the solution f1=±4 arctan √2(x+c2)±pa1+ 2(x+c2)2 √a1!.(14) In Fig. 2 we show the profile and the graphic of the kink soliton solution corresponding to λ=−4,a1=1 and c2=0,i.e. f1= 4 arctan √2(z−4t)−p1 + 2(z−4t)2 A!.(15) -40 -20 20 40 -6 -5 -4 -3 -2 -1 x f -40 -20 0 20 40 -4 -2 0 2 4 -6 -4 -2 0 t z Fig. 2. Representation of the profile and the graphic of the solution (15), λ=−4, a1= 1 and c2= 0,x=z−4t∈[−50,50],z∈[−50,50] and t∈[−5,5]. B. When c1=−1 2a1the radicand in (13) has the double roots ξ=−1and ξ= +1. The corresponding integral is divergent at ξ=±1and it is convergent at ξ=±∞. The branch that corresponds to the interval (−1,1) is defined for x∈R. Therefore, the solution f1of (9) determined by his implicitly given by the equation x+c2=±1 √−2a1 tan(f1/4) Z 0 ξ2+ 1 pξ8+ 6ξ4+ 1dξ , (16) New Classes of Stable Exact Solutions for a Nonlinear Rotational . . . 2009 where a1<0and c2are arbitrary constants. The class of solutions given by (16) corresponds to a kink soliton solution of (9). In Fig. 3 we exhibit the profile of one of these solutions. -4 -2 0 2 4 -3 -2 -1 0 1 2 3 x f Fig. 3. Representation of the profile hof the solution f1that corresponds to (16), c2= 0,a1=−1/2,x∈[−5,5],f1∈(−π, π). 3. Finally, taking p= 3, the change of the dependent variable is f1= 6 arctan(h), equation (9) becomes −4h3h2−333h2−1+ 3a1h2+ 14h00(h2+ 1) −2hh02= 0 . (17) This equation admits solutions that can be given implicitly by x+c2=± h Z−3a1(ξ2+ 1)2 p−2a1(2ξ6−3ξ4+12ξ2+1)(3ξ2−1)2+ 9c1a2 1(1+ξ2)6dξ . (18) This case is similar to the former one. The radicand has double roots when c1= 0 and c1= 2/(9a1). For c1= 0 the double roots are ξ=±1/√3; in this case there is a solution that is defined for x∈R. When c1= 2/(9a1)the double roots are ξ= 0 and ξ=±√3; in this case there are two solutions that are defined for x∈R; these solutions correspond to the intervals (−√3,0) and (0,√3). •Non transcendent form: The change of the dependent variable f1= (2p+ 1) arctan(h), where pis an integer number, allows us to write equation (9) in a non transcendent form.