% Sorry, no figures in this TeX file !!!
%
% to appear in proceedings of St. Petersbourg DNA conference (held June 92)
\font\tafont = cmbx10 scaled\magstep2
\font\tbfont = cmbx10 scaled\magstep1
\font\sixrm = cmr6
\font\petit = cmr9
%\magnification=\magstep1
%\vsize=19.3752 truecm
%\hsize=11.833 truecm
%\hfuzz=2pt
%\tolerance=500
%\abovedisplayskip=0 mm plus6pt minus 4pt
%\belowdisplayskip=2 mm plus6pt minus 4pt
%\abovedisplayshortskip=0mm plus6pt minus 2pt
%\belowdisplayshortskip=2 mm plus4pt minus 4pt
%\predisplaypenalty=0
\clubpenalty=10000
\widowpenalty=10000
\newcount\notenumber
\def\clearnotenumber{\notenumber=0\relax}
\def\fnote#1{\advance\notenumber by 1
\footnote{$^{\the\notenumber}$}{\petit #1}}
\newcount\fignumber
\def\clearfignumber{\fignumber=0\relax}
\def\fgn{\advance\fignumber by 1
\the\fignumber }
\newcount\tabnumber
\def\cleartabnumber{\tabnumber=0\relax}
\def\tabn{\advance\tabnumber by 1
\the\tabnumber }
\newcount\titbnumber
\def\cleartitbnumber{\titbnumber=0\relax}
\def\tbn{\global\advance\titbnumber by 1
\the\titbnumber }
\def\titlea#1{{~ \vskip 1 truecm \tafont {\centerline {#1}}
} \vskip 1 truecm }
\def\titleb#1{ \bigskip {\tbfont #1 } \medskip}
\parindent=10pt
\parskip=0pt
\def\Z{{\bf Z}} %integers
\def\R{{\bf R}}%reals
\def\i{{(i)}}
\def\unm{{1 \over 2}}
\def\grad{\nabla}
\def\Lapl{\triangle}
\def\dd{{\rm d}} %"d dritto"; to be used in integrals
\def\und{{1 \over \delta}}
\def\pa{\partial}
\def\dint{\int_{- \infty}^{+ \infty} }
\def\arcsinh{{\mathop{\rm arcsinh}\nolimits}}
\def\W{{\cal W}}
\def\N{{\cal N}}
\def\L{{\cal L}}
\def\A{\AA} %Angstrom
\def\er{~{\rm eV/rad^2}}
\def\ea{~{\rm eV/\A^2}}
\def\ec{~{\rm erg/cm^2}}
\def\eb{~{\rm erg/rad^2}}
\def\ps{{\rm ~ ps}}
\def\a{\alpha}
\def\b{\beta}
\def\g{\gamma}
\def\d{\delta}
\def\eps{\varepsilon}
\def\k{\kappa}
\def\l{\lambda}
\def\r{\rho}
\def\s{\sigma}
\def\th{\theta}
\def\om{\omega}
\def\Om{\Omega}
\def\.#1{{\dot #1}}
\def\~#1{{\widetilde #1}}
\def\^#1{{\widehat #1}}
\def\=#1{{\bar #1}}
\def\[{\left[}
\def\]{\right]}
\def\({\left(}
\def\){\right)}
%\nopagenumbers
\titlea{Simple Models of DNA Dynamics}
\vfill
\centerline{Giuseppe Gaeta}
\centerline{\it Centre de Physique Theorique}
\centerline{\it Ecole Polytechnique, F-91128 Palaiseau (France)}
\centerline{\tt gaeta@orphee.polytechnique.fr}
\vskip 1 truecm
\titleb{Summary}
In this note I briefly discuss the possibility of considering meaningful
{\it simple} models describing DNA nonlinear dynamics; the focus is mainly on
base opening dynamics -
which is easily accessible to experimental observation - and by simple I mean
models which have only
one or two degrees of freedom per segment (i.e. nucleoside plus base) of the
double helix. After
discussing the feature such a model is required to have, I consider two models
(Y and PB models)
recently introduced, and show that despite their semplicity (relative, since
they are not solvable),
they succeed to grasp some essential feature of DNA dynamics, and exhibit good
agreement with
experimental observation. A more complete version of the discussion and results
presented here will
be given in the review paper [6].
%\medskip {\tafont I - Modelling DNA dynamics}
\titleb{\tbn. Structure of DNA}
The molecule of DNA is both extremely complicate and big; although there are
attempts of "ab initio"
molecular computations of such an intricate system, we will just set for an
attempt to identify the
(few) degrees of freedom relevant to the processes we want to study, which are
essentially the {\it transcription} of DNA sequence by (Y model) and the
{\it spontaneous denaturation} of DNA (PB model), and in general relevant at
body temperature, and
build an effective model for the evolution of these. Indeed, as we will see in
the following, such an
approach can be quite successfull.
The chemical structure of DNA is well known and will not be described here; for
informations on the structure and properties of DNA, the reader is referred to
[1]. We just recall that the
distance between successive sites of the double helix, i.e. the lenght of each
segment, is about
3.4 Angstrom, and the pitch of the helix is usually of 10 sites. In order to
form the DNA molecule,
two such helices couple, by hydrogen bonds linking bases at corresponding sites
on the two
helices; bases
linked by H bonds are said to form a pair.
It should be immediately remarked - and stressed - that atoms between which the
H
bonds exist are not at all isolated: they are part of a very complicate
molecule,
which moreover behaves in a way strongly dependent on the physico-chemical
environment it lives in: this is to say that (although it is tempting to do
this,
and sometimes we will actually do it for the sake of simplicity) one is not
justified in modelling the H bonds between bases as if they were between
isolated
atoms.
A similar remark holds also for the other interactions between atoms of the
molecule; anyway, as we will remark in a moment, the excitation energy of many
degrees of freedom is well above body temperature, so that it seems quite safe
to consider
these degrees of freedom as frozen. \eject
\titleb{\tbn. Interactions in DNA}
The chemical (covalent) bonds between atoms of the DNA molecule have actually,
as
just mentioned, for the greater part of cases excitation energies well above
the
temperature at which DNA lives. In these conditions, we should
consider the backbone as rigid, with no "small oscillations"
possible: although we are going to discuss classical mechanics models for the
degrees of freedom of interest, we should not forget that the molecule obeys
quantum mechanics, so that even in classical modelling it would be physically
{\it
wrong} to allow for motions of these "hard" degrees of freedom.
Let us look more closely at the structure of DNA, with special attention to
"soft"
degrees of freedom, i.e. to degrees of freedom excited at body temperature.
These
are of three kinds:
{\it a)} the H bonds between pairing bases;
{\it b)} the bases are essentially free to rotate around the N-C bond linking
them
to the nucleoside\fnote{From now on this will be "the" N-C bond, for ease
of notation.} (free with respect to covalent interactions; interactions of
other
nature - e.g. electrostatic - with the other atoms in the molecule give a
restoring torque, as we will see later);
{\it c)} the sugar ring in the nucleoside has nearly degenerate equilibrium
positions, separated by a potential barrer, and can tunnel between the
equilibrium positions on the two sides of the plane. This is known as
{\it sugar puckering}, and the excitation energy of such sugar puckering modes
is
similar to that of the H bonds. Sugar puckering affects the geometry of the
double
helix, and one of its effects is to take apart bases in a H bonded pair, i.e.
to
move radially (with respect to the axis of the double helix) the base attached
to
the sugar affected by puckering.
The identification of these degrees of freedom as the only ones which should
actually be affected by body-temperature dynamics is actually enough to
introduce
and justify the models we want to discuss. Indeed, due to the great complexity
of
the molecule, the models (and the parameters appearing in them) should be
regarded
as effective ones for the models of these selected degrees of freedom, i.e. we
should not try to deduce details of the models and value of the parameters from
first principles and the structure of the molecule, but just find a model which
describes coherently with experimental observations the motion of these degrees
of
freedom.
We will anyway mention some other details on interactions in DNA, before
determining
the models; this will permit to fix terminology.
Usually, the interactions between bases are divided into {\it pairing} and {\it
stacking} ones; the former are essentially due to H bonds, while the latter are
due to a number of effects in the molecule, and their net effect is to keep
stacked one over the other the bases on each helix. The stacking interactions
also
greatly influence the pairing.
The energy of H bond pairing also
greatly differs for the pairs (A-T) and (G-C): indeed, this is no surprise,
since in one case the
bases are couples by three H bonds, while in the other one by two H bonds.
The stacking interactions can be thought as responsible for the restoring
forces
trying to keep the bases in their equilibrium positions independently of the H
bonds, i.e. the forces opposing the "rotational" and "vibrational" movements of
points {\it b)} and {\it c)} above.
Indeed, one observes experimentally that once three consecutive pairs of bases
are
bridged by H bonds, stacking forces to continue the closing of the double
helix;
this also gives a brute extimate of 1/3 for the ratio between pairing effects
due to
stacking and direct pairing forces due to H bonds.
\titleb{\tbn. Local openings of the double helix}
The local openings of the double helix, i.e. the breakings of the H bonds
between
pairing bases, are accessible to experimental observation. This is usually done
through study of absorption of UV rays; see e.g. [1] for details. From our
point of view, the base pairs openings will be the central object of our
study.
It should be ramarked, in relation to the above discussion of accessible
degrees
of freedom, that the H bonds can be broken either by increasing the distance
between pairing bases (e.g. by the geometrical effect related to sugar
puckering,
as described above), either by rotation of the bases around the N-C bond
connecting them to the relative nucleoside. Indeed, H bonds are often strongly
directional, so that a change in direction of the bases, although comporting
only
a small variation in the relative distance, can be quite effective in breaking
the
H bonds.
As already remarked, the great complexity of the molecule suggests to keep the
above considerations at a very qualitative level, as the H bonds have no
apriori reason to behave like those between isolated atoms.
\titleb{\tbn. Models for DNA dynamics}
>From the above short discussion, it follows that our simple models for DNA
dynamics should have several features:
{\it a)} Consider only the two degrees of freedom per helix segment, i.e. per
base, considered above; from now on we will cal "vibrational" the dynamics
associated to sugar puckering and rectilinear (radial with respect to
the axis of the double helix) separation of the bases in a pair, and
"torsional"
the dynamics associated to rotations of the bases around the C-N bond linking
them to the sugar ring in the nucleoside.
{\it b)} Consider direct pairing and stacking forces.
{\it c)} Not depend on details of the forces between bases, but only on
qualitative features. On the other side, and always due to our lack of
understanding of such a complex molecule, the constants (coupling constants and
parameters) appearing in the model should be seen as effective ones, i.e. we
have
some freedom in their choice.
We will moreover choose to describe the dynamics we are interested in in terms
of
classical mechanics rather than in quantum terms (but, as explained before, the
restriction to the degrees of freedom we actually consider is justified by
quantum
mechanical considerations); this is admittedly a matter of preference, and
other
approaches would be possible, considering quantum mechanical models or even
looking at the double helix as an Ising system (the two states of the
"spin" at a site would correspond then to "closed" or "open" H bonds between
bases
at the site), see e.g. [2] and [3].
Although a model having these features would already be extremely simplified
(we
are in fact considering by now only two degrees of freedom per base), we want
to operate a further
simplification. Namely, we want to consider separately models for the
vibrational
dynamics alone and for the torsional dynamics alone (vibrational and torsional
in
the sense explained above), so to have only {\it one} degree of freedom per
base.
We will moreover introduce still another further severe simplification: we will
assume that all the bases are equal: i.e. they have the same mass and moment of
inertia, and all couplings of the same type are equal. (We will moreover assume
the
double helix to be of infinite lenght, as we suppose side effects to be
inessential).
Such models were introduced by Peyrard and Bishop [4] for the vibrational
dynamics, and by Yakushevich [5] for the torsional one; we will sometime refer
to
them as "planar" ones for reasons to be evident in the following. They will
also
be indicated, for short, as PB model and Y model respectively.
It must be stressed that while the reduction to two degrees of freedom per base
discussed above was well justified physically, the present splitting of the
dynamics in two noninteracting parts is admittedly due to an attempt to deal
with
a manageable problem rather than to physical considerations. We will see anyway
that, somehow rather surprisingly, these "one degree of freedom models" are
able
to capt relevant experimentally observed features of DNA dynamics; this will
give
a justification {\it a posteriori} of our splitting.
We will write down the models in term of kinetic ($T$) plus potential ($V$)
energy.
If $\phi_n^\i$ ($i=1,2$, $n \in \Z$) is the generalized coordinate describing
the
base at site $n$ on the $i$-th chain of the double helix - which we suppose of
infinite lenght - the Hamiltonian will be written as
$$ H = T (\.\phi ) + V (\phi ) \eqno(1) $$
where we assume
$$ T = \sum_{i,n} \unm \~m (\.\phi_n^\i )^2 \eqno(2) $$
with $\~m$ the mass or the moment of inertia of the bases, according to the
degree
of freedom taken in consideration.
As for the potential energy, we will decompose it into three parts: a pairing -
or
{\it transverse} - interaction energy $V_T$, a stacking interaction energy
$V_S$,
and a third kind of energy (helicoidal interaction energy) $V_H$ to be
discussed
later on and which for the moment will just be set to zero. So we will have
$$ V = V_T + V_S \eqno(3) $$
(models which satisfy (3), i.e. for which no "helicoidal" interaction is
present,
will be called {\it planar}), and we assume that
$$ V_T = \sum_n U_T (\phi_n^{(1)} , \phi_n^{(2)} ) ~~;~~ V_S = \sum_{i,n} U_S (
\phi_{n+1}^\i , \phi_n^\i ) \eqno(4) $$
We will use coordinates such that $\{ \phi_n^\i = 0 \}$ correspond to the
equilibrium stable configuration.
It should be stressed that, since stacking is considerably
stronger\fnote{It is maybe worth remarking that there is no contradiction
with the extimate of 1/3 for the ratio of pairing effects due to stacking to
direct pairing forces given above: here we consider the strenght of the
stacking
interaction themselves, and not just of their pairing effect.} than H
bonding, the nonlinear effects will appear essentially only in the direct
pairing
(transverse) interactions, i.e. the dynamic will be such to stretch the
stacking
interactions much less that the H bonds. These observations suggest to consider
the
harmonic potential approximation for stacking interactions, which we will
indeed
adopt. It should be stressed that, as it will be clear in the following, this
leads
to the appearance of the Laplacian once we pass to the continuum limit.
We pass now to describe two prototypical models of torsional and vibrational
DNA
dynamics, due respectively to Yakushevich and to Peyrard and Bishop. Earlier
proposals of similar models were formulated by several authors, but will not be
discussed here.
%{\tafont II - Planar models}
\titleb{\tbn. A model of DNA torsion dynamics}
In the model proposed by L. Yakushevich [5], the pairing interaction is
considered
in the approximation of harmonic potential, depending on the distance between
the
endpoints of the bases, i.e. the atoms bridged by H bonds. Although the force
is
linear in this distance, when we describe the position of the bases by means of
an
angle $\th$ of rotation around the N-C bond the force becomes a nonlinear
function
of such an angle.
The details of pairing interactions in such a model are shown in fig.\fgn.
Here the distance $\ell$ between points $p_1$ and $p_2$ satisfies
$$ \ell^2 = \[ d + 2r - r( \cos \th_1 + \cos \th_2 ) \]^2 + \[ r (\sin \th_1 +
\sin \th_2 ) \]^2 \eqno(1) $$
It is then convenient to pass to coordinates
$$ \phi_\pm = {\th_1 \pm \th_2 \over 2 } \eqno(2) $$
in terms of which, by standard trigonometric identities, the distance is
written as
$$ \ell^2 = d^2 + + 4 r^2 \( 1 + \cos^2 \phi_- - 2 \cos \phi_+ \cos \phi_- \)
\eqno(3) $$
Yakushevich does also suggest to consider the approximation $ d \to 0$.
In other words we are considering the stretchning of H bonds due to small
rotations
of the bases and extend this to the full range $[0, 2 \pi ]$ of angular
coordinates.
Notice that we are trascurating the directionality effects, as only distance is
taken into account.
The potential describing the pairing interaction will then be
$$ U_T = \unm \a [ 1 + \cos^2 \phi_- - 2 \cos \phi_+ \cos \phi_- ] \eqno(4) $$
(notice that this is $2 \pi$-periodic in both angular coordinates, as it
should);
this is sketched in fig.\fgn by level curves.
As for stacking interactions, we will also {\it a fortiori} keep the harmonic
approximation for them; if $\th_n^\i$ is the angular coordinate relative to the
$n$-th base along the $i$-th chain ($i=1,2$), the stacking interaction
potential
between bases at sites $n$ and $n+1$ will be
$$ U_S = \unm \b \[ ( \th_{n+1}^{(1)} - \th_n^{(1)} )^2 + (\th_{n+1}^{(2)} -
\th_n^{(2)} )^2 \] \eqno(5) $$
Passing again to the coordinates (2), i.e.
$$ \th_n^{(1)} = \phi_n^{(+)} + \phi_n^{(-)} ~~;~~ \th_n^{(2)} = \phi_n^{(+)} -
\phi_n^{(-)} \eqno(6) $$
the above reads
$$ V_S \unm \b \[ ( \phi_{n+1}^{(+)} - \phi_n^{(+)} )^2 + ( \phi_{n+1}^{(-)} -
\phi_n^{(-)} )^2 \] \eqno(7) $$
In the following, in order to avoid confusion and to lighten notation, we will
prefer to write $\psi$ for $\phi_+$ and $\chi$ for $\phi_-$; with this the
planar
Yakushevich model is defined by the Hamiltonian
$$ H = T + V_S + V_T \eqno(8) $$
$$ T = \sum_n \unm I \[ ( \.\th_n^{(1)} )^2 + ( \.\th_n^{(2)} )^2 \] = \sum_n
\unm
I \( \.\psi_n^2 + \.\chi_n^2 \) \eqno(9) $$
$$ V_S = \sum_n \unm \b \[ ( \psi_{n+1} - \psi_n )^2 + ( \chi_{n+1} - \chi_n
)^2
\] \eqno(10) $$
$$ V_T = \sum_n \unm \a \[ 1 + \cos^2 \chi_n - 2 \cos \psi_n \cos \chi_n \]
\eqno(11) $$
where we introduced the $\psi , \chi$ coordinates in the kinetic term $T$ as
well. The equilibrium (stable) configuration is represented by $\{ (\psi_n ,
\chi_n ) = (0,0) ~ \forall n \}$; $I$ is the moment of inertia of the bases
around the N-C bond.
The above hamiltonian gives as equations of motion
$$ \eqalign{
I {\ddot \psi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) - \a \sin \psi_n
\cos \chi_n \cr
I {\ddot \chi_n} =& \b (\chi_{n+1} - 2 \chi_n + \chi_{n-1} ) - \a \sin \chi_n (
\cos \psi_n - \cos \chi_n ) \cr } \eqno(12) $$
which for small oscillations reduce to
$$ \eqalign{
I {\ddot \psi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) - \a \psi_n \cr
I {\ddot \chi_n} =& \b (\chi_{n+1} - 2 \chi_n + \chi_{n-1} ) \cr } \eqno(13)
$$
\titleb{\tbn. The continuum limit of planar Yakushevich model}
We do now want to consider the continuum limit of this model: the arrays $\{
\psi_n (t) \}$, $\{ \chi_n (t) \}$ will then be substituted by fields $\psi
(x,t)
, \chi (x,t)$ which take values in the circle $S^1 = [0,2 \pi ]$ and such that
$$ \eqalign{ \psi ( n \d , t ) &= \psi_n (t) \cr \chi ( n \d , t ) &= \chi_n
(t)
\cr } \eqno(1) $$
where $\d$ is the unit distance between successive bases along the axis of the
double helix, $\d \simeq 3.4 \A$.
It should be stressed that in considering the continuum approximation, it makes
no
sense to consider wavelenght $\l$ shorter than the interbase distance, i.e.
than
$\d$.
We can now rewrite the various terms of the hamiltonian in the field
formulation;
we get
$$ \eqalignno{
T =& \und \int \unm I ( \psi_t^2 + \chi_t^2 ) \dd x &(2) \cr
V_T =& \und \int \unm \a [1 + \cos^2 \chi - 2 \cos \psi \cos \chi ] \dd x &(3)
\cr
V_S =& \und \int \unm \d^2 \b \[ ( \grad \psi )^2 + ( \grad \chi )^2 \] + O
(\d^4
) \dd x &(4) \cr } $$
where $\grad \equiv \pa / \pa x$ and the $(1/\d )$ factors are only for
normalization purposes, i.e. to keep the normalization of the field hamiltonian
equal to that of the discrete one.
The evolution equations for the fields $\psi , \chi $ issuing from this
hamiltonian (or, equivalently, the continuum limit of the equations of motion
for
$\psi_n , \chi_n$) are
$$ \eqalign{ I \psi_{tt} =& \b \d^2 \Lapl \psi - \a \sin \psi \cos \chi +
O(\d^4 )
\cr
I \chi_{tt} =& \b \d^2 \Lapl \chi - \a \sin \chi ( \cos \psi - \cos \chi ) +
O(\d^4 ) \cr } \eqno(5) $$
where $\Lapl$ is the spatial Laplacian, $\Lapl \equiv \pa^2 / \pa x^2 $
We stress that here $\d$ is a physical parameter, fixed by nature, so that we
cannot simply consider the limit $\d \to 0$ in order to get rid of the $0(\d^4
)$
terms. We will not discuss this point in detail here (see [6] for a
discussion), but
just state that if $\l_0$ is the shorter wavelenght appearing in the Fourier
decomposition of $\psi$,
we have to ask not only $\l_0 > \d$, but
$$ \l_0 >> 2 \d \eqno(6) $$
\titleb{\tbn. Solitons in planar Yakushevich model}
We want now to focus our attention on soliton solutions for the Yakushevich
model. The soliton solutions propagate without changing their form, i.e. we
have
$$ \eqalign{ \psi (x,t) =& \psi (x - vt) \equiv \psi (z) \cr
\chi (x,t) =& \chi (x - vt) \equiv \chi (z) \cr } \eqno(1) $$
This ansatz, once used into the field evolution equations, produces the {\it
equations for Yakushevich
solitons}:
$$ \eqalign{
(v^2 - \k ) \psi_{zz} =& - {\pa W \over \pa \psi} \equiv - \sin \psi \cos \chi
\cr
(v^2 - \k ) \chi_{zz} =& - {\pa W \over \pa \chi} \equiv - \sin \chi ( \cos
\psi
- \cos \chi ) \cr} \eqno(2) $$
Corresponding, the hamiltonian and lagrangian for Yakushevich solitons read
$$ \eqalignno{
H_S =& {I \over \d} \int \unm (v^2 + \k ) \[ (\grad \psi )^2 + (\grad \chi )^2
\]
+ W(\psi , \chi ) \dd x &(8) \cr
L_S =& {I \over \d} \int \unm (v^2 - \k ) \[ (\grad \psi )^2 + (\grad \chi )^2
\]
- W(\psi , \chi ) \dd x &(9) \cr } $$
where now $\grad \equiv \pa / \pa z$.
On physical grounds, we would require that for $z \to \pm \infty$ the fields go
to the equilibrium (zero-energy) configuration; due to the $2 \pi$ periodicity
of
$W$ in both $\psi$ and $\chi$, this means
$$ \eqalign{
& \psi_z ( \pm \infty ) = \chi_z ( \pm \infty ) = 0 \cr
& \psi ( \pm \infty ) = 2 \pi n_\pm ~~;~~ \chi (\pm \infty ) = 2 \pi m_\pm \cr
}
\eqno(10) $$
With no loss of generality we can assume
$$ n_- = m_- = 0 \eqno(11) $$
so that the finite energy solutions are classified by $(n_+ , m_+ ) \equiv
(n,m)$ alone, i.e. by two integers. This corresponds indeed to the topological
classification of solitons $\R \to S^1 \times S^1$, as $n$ and $m$ are {\it
rotation numbers} for the fields $\psi$ and $\chi$.
Solitons in the $(n,m)$ class will be topologically stable, and
we can look for them by minimizing the energy with prescribed boundary
conditions. We will not, anyway, use the topological theory of
solitons [7], and
keep the discussion at a simple level.
By interpreting $z$ as a time, and $\~\mu \equiv (v^2 - \k )$ as an effective
mass, (7) aquires the meaning of a mechanical (Newton) equation for a particle
of
mass $\~\mu$ moving in the two dimensional plane (with coordinates $\psi ,
\chi$) under the effects of the potential $W$. The boundary conditions given
above fix then the total energy to be $E=0$; it is then immediate to see that
nontrivial solutions can exist only if the effective mass is negative; this
also gives an upper bound
on the speed of solitons: $$ \eqalign{ & \~\mu = v^2 - \k < 0 \cr & |v| <
\sqrt{\k } = v_{{\rm max}}
\cr } \eqno(12) $$
which in turn implies
$$ \mu \equiv - \~\mu \le \k = \mu_{{\rm max}} \eqno(13) $$
>From now on we will use $\mu = - \~\mu = \k - v^2 $; notice that in terms of
this
effective mass, the effective potential will be $-W$. Equations (7) are now
$$ \eqalign{
\mu \psi_{zz} =& {\pa W \over \pa \psi } \cr
\mu \chi_{zz} =& {\pa W \over \pa \chi } \cr } \eqno(14) $$
Multiplying the first one by $\psi_z$ and the second one by $\chi_z$, we get
$$ \unm \mu { \dd ~ \over \dd z} \( \psi_z^2 + \chi_z^2 \) = { \dd ~ \over \dd
z}
W( \psi , \chi ) \eqno(15) $$
so that integrating - and using the boundary conditions (10), which sets the
integration constant to zero - we get
$$ \unm \mu \[ ( \grad \psi )^2 + ( \grad \chi )^2 \] = W (\psi , \chi )
\eqno(16) $$
This relation is satisfied by soliton solutions $( \=\psi , \=\chi )$; using it
in (8), and writing $\nu = \k + v^2$, we have for the energy of such a solution
$$ \eqalign{
E =& \unm { I \over \d} \int (\nu + \mu ) \[ ( \grad \=\psi )^2 + ( \grad
\=\chi
)^2 \] \dd z \cr
E =& {I \k \over \d } \int \[ ( \grad \=\psi )^2 + ( \grad \=\chi )^2 \] \dd z
\cr } \eqno(17) $$
\titleb{\tbn. Explicit soliton solutions}
In the simplest cases, i.e. for the $(1,0)$ and the $(0,1)$ soliton solutions
(see the discussion of boundary conditions in the previous section), it is
actually possible to give an analytic expression for the soliton solutions.
This is due to the fact that in this case only one of the fields is nonzero,
so
that we are reduced to a one-dimensional mechanical problem; a short look at
the
effective potential $-W$, either in its analytic form either in fig.2, suffices
to
get convinced that $(n,0)$ and $(0,n)$ soliton solutions with $|n| > 1$ require
{\it both} fields to be not identically zero\fnote{Obviously, this is
{\it a fortiori} true if both the soliton numbers $(n,m)$ are not zero.}, so
that
we can not reduce to a one dimensional motion.
We will not give here the details of the (standard) computation of these
explicit
solutions, but just their final form; the reader is referred to [6,8] for
details.
In the case of the (1,0) soliton, we have with $ A = \sqrt{2 \a / \mu }$ the
explicit solution
$$ \psi (z) = \arccos \[ { \sinh^2 (Ax) - 1 \over \sinh^2 (Ax) + 1
} \] \eqno(1) $$
For the (0,1) soliton, we have with $B = 2 \sqrt{\a / \mu}$,
$$ \chi (z) = \arccos \[ {B^2 z^2 -1 \over B^2 z^2 + 1 } \] \eqno(2) $$
It should be remarked that both (1) and (2) obey simple scaling relations under
changes of $\a , \mu$ and of the lenght scale for $z$; this is not a fortuitous
case, but indeed corresponds to general scaling properties of soliton solutions
for the {\it planar} Yakushevich model, as discussed elsewhere [6,8].
\titleb{\tbn. A model of DNA vibration dynamics}
We pass now to consider the "vibrational" dynamics of DNA, associated to sugar
puckering modes; as recalled before, this corresponds to radial motion of the
bases with respect to the axis of the double helix; the coordinate $x_n^\i$
describing the position of the base will then be its distance from the axis of
the
double helix, or better - in order to keep $\{ x_n^\i = 0 \}$ as the
equilibrium
configuration - this distance minus the distance corresponding to the
equilibrium
position.
Such a model was proposed by Peyrard and Bishop [4], who suggested to use as
pairing potential the Morse potential; namely, if $x_n^\i$ is the coordinate
defined above for the $n$-th base on the $i$-th chain, we have
$$ U_T = D \[ e^{- \a (x^{(1)} + x^{(2)} )/2} - 1 \]^2 \eqno(1) $$
or, passing again to the coordinates
$$ \phi_\pm = { x^{(1)} \pm x^{(2)} \over 2 } ~~;~~ \phi_+ \equiv \psi ~,~
\phi_-
\equiv \chi \eqno(2) $$
we rewrite this as
$$ U_T = D \[ e^{- \a \psi} - 1 \]^2 \eqno(3) $$
As for stacking interactions, we keep to the harmonic potential approximation,
which gives
$$ U_S = \unm \b \[ \( x_{n+1}^{(1)} - x_n^{(1)} \)^2 + \( x_{n+1}^{(2)} -
x_n^{(2)} \)^2 \] \eqno(4) $$
or, passing to the $\psi , \chi$ coordinates,
$$ U_S = \unm \b \[ \( \psi_{n+1} - \psi_n \)^2 + \( \chi_{n+1} - \chi_n \)^2
\]
\eqno(5) $$
Finally, the hamiltonian for this (Peyrad-Bishop, or PB) model will be
$$ H = T + V_S + V_T \eqno(6) $$
$$ \eqalignno{
T =& \sum_{i,n} \unm m \( x_n^\i \)^2 = \sum_n \unm m \( \.\psi_n^2 +
\.\chi_n^2
\) &(7) \cr
V_S =& \sum_n \unm \b \[ \( \psi_{n+1} - \psi_n \)^2 + \( \chi_{n+1} - \chi_n
\)^2
\] &(8) \cr
V_T =& \sum_n \unm D \[ e^{-\a \psi_n } -1 \]^2 &(9) \cr } $$
Here $m$ is the mass of the bases; the parameters $\a , \b$ are {\it not} the
same
as in the Yakushevich model.
In discussing small oscillations around the equilibrium configuration $\psi_n =
\chi_n =0$, we need the quadratic part of $V_T$; this is simply
$$ V_T^{(2)} = \sum_n \unm ( D \a^2 ) \psi_n^2 \eqno(10) $$
The above hamiltonian gives as equations of motion
$$ \eqalign{
m {\ddot \psi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) - \a D \[ e^{- \a
\psi_n } - 1 \] \cr
m {\ddot \chi_n} =& \b (\chi_{n+1} - 2 \chi_n + \chi_{n-1} ) \cr } \eqno(11)
$$
which for small oscillations reduce to
$$ \eqalign{
m {\ddot \psi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) - \a^2 D \psi_n
\cr
m {\ddot \chi_n} =& \b (\chi_{n+1} - 2 \chi_n + \chi_{n-1} ) \cr } \eqno(12)
$$
\titleb{\tbn. The continuum limit of planar Peyrard-Bishop model}
In order to pass to the continuum limit of this model, we proceed similarly to
the
case of the Yakushevich model; i.e., the arrays $\psi_n (t) , \chi_n (t)$ are
substituted by fields $\psi (x,t) , \chi (x,t)$ such that
$$ \eqalign{ \psi (n \d ,t ) =& \psi_n (t) \cr
\chi ( n \d ,t) =& \chi_n (t) \cr} \eqno(1) $$
with $\d$ the distance between successive bases along the axis of the double
helix, $\d \simeq 3.4 \A$. Notice that $\chi$ takes value on the whole real
line
$\R$, while for $\psi$ the physical range is $\psi \ge -2 d_0$, where $d_0$ is
the
distance of the bases from the axis of the double helix at equilibrium.
Actually,
for both $\psi$ and $\chi$ it makes no sense to consider too big absolute
values;
th erelevant point here for later discussion is that, contrary to the
Yakushevich
case, the topology of the target space is trivial (this means in
particular that we have no topological solitons).
With the same remarks as in the Yakushevich case, we can write the hamiltonian
in
the field formulation as
$$ H = T + V_T + V_S \eqno(2) $$
$$ \eqalignno{
T =& \und \int \unm m \( \psi_t^2 + \chi_t^2 \) \dd x &(3) \cr
V_T =& \und \int \unm D \( e^{- \a \psi} - 1 \)^2 \dd x &(4) \cr
V_S =& \und \int \unm \d^2 \b \[ ( \grad \psi )^2 + ( \grad \chi )^2 \] \dd x
&(5)
\cr } $$
and the corresponding field evolution equations are
$$ \eqalign{
m \psi_{tt} =& \b \d^2 \Lapl \psi + \a D \( e^{-\a \psi} - 1 \) \cr
m \chi_{tt} =& \b \d^2 \Lapl \chi \cr } \eqno(6) $$
We recall that these are valid if and only if the solutions have wavelenghts
large
compared to $\d$.
%{\tafont III - Helicoidal models}
\titleb{\tbn. Helical geometry and helicoidal DNA models}
The models discussed up to now have been at several reprises qualified as {\it
planar}; this is due to the fact that they can be schematized as in fig.\fgn,
which
also means that they don't take into account the peculiar geometry of the
double
helix.
It turns out that taking this into account, even in a "minimal" way, produces
qualitative changes in the dynamics, quite relevant also in relation with
"observable" properties of the model (especially in the case of the PB model).
The modification we consider here was first proposed in the context of the
Yakushevich model, and then applied also to the Peyrard- Bishop one [9,10]; the
models modified in this way will be qualified as {\it helicoidal}. In
the following, we will write hY for "helicoidal Yakushevich", and hPB for
"helicoidal Peyrard- Bishop".
In a double helix (with pitch of the helix $2h$ in base units) it happens that
bases which are half-wind of the helix apart, i.e. $h$ sites away, come to be
spatially near, so that they can interact; this fact is shown schematically in
fig.\fgn. It should be stressed that such an interaction is actually present in
real
DNA, where it is realized by bridging water filaments between the concerned
bases.
>From the point of view of our models, and referring to the notation previously
employed,
this means that the hamiltonian will now be of the form
$$ H = T + V_T + V_S + V_H \eqno(1) $$
where $V_H$ represents these new or - as they will be called - {\it helicoidal}
interactions. These will appear through terms of the form, referring again to
section * for notation,
$$ U_H = k \( \phi_n^{(1)} , \phi_{n+h}^{(2)} \) + k \( \phi_n^{(2)} ,
\phi_{n+h}^{(1)} \) \eqno(2) $$
It should be remarked that helical geometry poses some constraint on the choice
of
coordinates for helicoidal models (the present author was actually responsible
for
introducing a wrong choice of coordinates). The choice of coordinates made
above for
both models we consider took this problem into account, and is therefore
correct.
As mentioned above, the helicoidal interaction is due in real DNA to bridged
water
filaments; given that the distance between bases coupled by these is $5 h
\simeq
17.0 \A$, from the geometry of water molecules and the size of O-H-O bonds one
can
estimate these filaments to be
made up of $\sim$ 8 ${\rm H}_2 {\rm O}$ molecules;
stretching the filament amounts to stretching $\sim$ 9 H bonds; this also means
that
the corresponding elastic constant should be $ \sim 1/9$ of that assumed for
single H bonds. Due to
this, these bonds will be stretched much less than the H bonds
pairing bases at corresponding sites on different helices (i.e. than the H
bonds
entering in $V_T$). We are therefore perfectly justified in considering the
harmonic approximation for $V_H$, which we will indeed do in the following.
As it is the case for all the interactions in our simple models, also the $V_H$
introduced here should be seen as an effective one; therefore one should be
most
careful in extending the above discussion, which justified the harmonic
approximation, to quantitative matters, i.e. to an extimation for the ratio
of the coupling constants relative to $V_T$ and $V_H$ (this should be $\sim 9$
by
the above argument). Indeed, e.g; in the context of the hPB model Dauxois gave
quite different estimates in order to fit experimantal findings [10].
We will leave the coupling constant for $V_H$, as already done for the other
constants appearing in the model, as an effecive constant to be adjusted. We
als
stress that, given the relatively long range of this effective interaction, it
is
potentially even more sensitive than the other ones to the effect of the
complicate
interactions in the molecule and with its environment.
\titleb{\tbn. Introduction of helicoidal terms in the Peyrard- Bishop model}
We will begin by shortly considering the effect of introducing "helicoidal"
terms
in the PB model, as was first done by Dauxois [10].
We add to the PB hamiltonian a term $V_H$ which, coherently with
the discussion of the previous section, is simply an harmonic potential:
$$ U_H = \unm \g \[ \( x_{n+h}^{(1)} - x_n^{(2)} \)^2 +\( x_{n+h}^{(2)} -
x_n^{(1)}
\)^2 \] \eqno(1) $$
Upon introduction of the $\psi, \chi$ coordinates, so that
$$ x_n^{(1)} = \psi_n + \chi_n ~~;~~ x_n^{(2)} = \psi_n - \chi_n \eqno(2) $$
the interaction energy reads, with straightforward algebra,
$$ U_H = \unm \g \[ \( \psi_{n+h} - \psi_n \)^2 + \( \chi_{n+h} + \chi_n \)^2
\]
\eqno(3) $$
Notice the different signs in the $\psi$ and $\chi$ parts of this.
The total hamiltonian of the helicoidal Peyrard- Bishop (hPB) model reads
therefore
$$ H = T + V_S + V_T + V_H \eqno(4) $$
with $T, V_S , V_T $ as in the PB model and
$$ V_H = \sum_n \unm \g \[ \( \psi_{n+h} - \psi_n \)^2 + \( \chi_{n+h} +
\chi_n
\)^2 \] \eqno(5) $$
with $h$ the half-pitch of the helix in base units (i.e. $h=5$ for ordinary
DNA).
The hPB hamiltonian gives as equations of motion
$$ \eqalign{
m {\ddot \psi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) + \g ( \psi_{n+h}
- 2
\psi_n + \psi_{n-h} ) - \a D \[ e^{- \a \psi_n } - 1 \] \cr
m {\ddot \chi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) - \g (
\psi_{n+h} - 2
\psi_n + \psi_{n-h} ) \cr } \eqno(6) $$
which for small oscillations reduce to
$$ \eqalign{
m {\ddot \psi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) + \g ( \psi_{n+h}
- 2
\psi_n + \psi_{n-h} ) - \a^2 D \psi_n \cr
m {\ddot \chi_n} =& \b (\chi_{n+1} - 2 \chi_n + \chi_{n-1} ) - \g ( \psi_{n+h}
- 2
\psi_n + \psi_{n-h} ) \cr } \eqno(7) $$
Let us now consider the continuum limit of this model. Clearly the helicoidal
interaction leads to the appearance of nonlocal (in $x$, the coordinate along
the
axis of the double helix) terms; the field evolution equations will be,
understanding that $\psi \equiv \psi (x,t)$, $\chi \equiv \chi (x,t) $ unless
otherways indicated, and writing $w = h \d$,
$$ \eqalign{
m \psi_{tt} =& \b \d^2 \Lapl \psi + \g \W_+ (\psi ) + \a D \( e^{-\a \psi} - 1
\) \cr
m \chi_{tt} =& \b \d^2 \Lapl \chi - \g \W_- (\chi ) \cr } \eqno(8) $$
with $\W_\pm$ defined as
$$\W_\pm (\phi ) = \phi (x+w , t) \mp 2 \phi (x,t) + \phi (x-w,t) \eqno(9) $$
\titleb{\tbn. Dispersion relations for the helicoidal Peyrard- Bishop model}
We set
$$\k = \b \d^2 / m ~~,~~ \s = D \a^2 / m ~~,~~ \r = \g / m ~;$$
with these the equations for small oscillations in the hPB model are
$$\eqalign{ \psi_{tt} =& \k \Lapl \psi + \r \W_+ (\psi ) - \s \psi \cr
\chi_{tt} =& \k \Lapl \chi - \r \W_- (\chi ) \cr} \eqno(1) $$
Passing to Fourier components, i.e. to
$$ \eqalign{
\psi (x,t) =& F_{q,\om} e^{i qx} e^{i \om t} \equiv F_{q,\om} | q, \om > \cr
\chi (x,t) =& G_{q,\om } e^{i qx} e^{i \om t} \equiv G_{q,\om } | q, \om > \cr}
\eqno(2) $$
we have to take care of the action of the nonlocal operators $\W_\pm$ on $| q,
\om
>$: this is
$$ \W_\pm \( | q, \om > \) = \( e^{iqw} \mp 2 + e^{-iqw} \) | q, \om > = 2 \[
\cos
(qw) \mp 1 \] | q, \om > \eqno(3) $$
With this, the dispersion relations are immediately seen to be
$$ \eqalign{
\om_\psi^2 =& \k q^2 + 2 \r \[1 - \cos (qw) \] + \s \cr
\om_\chi^2 =& \k q^2 + 2 \r \[1 + \cos (qw) \] \cr }
\eqno(4)$$
It is maybe preferable to use standard trigonometric formulas and reexpress the
above as
$$ \eqalign{
\om_\psi^2 =& \k q^2 + 4 \r \sin^2 \({qw \over 2 } \) + \s \cr
\om_\chi^2 =& \k q^2 + 4 \r \cos^2 \({qw \over 2 } \) \cr } \eqno(5)$$
The dispersion relations (4),(5) are plotted in fig.\fgn (with values of the
parameters following the estimates given by [10]). It is immediately clear, by
comparing with the dispersion relations of the planar PB model, that the
introduction of helicoidal terms has drastically changed the (small
oscillations)
dynamics.
Notice that in the PB and hPB models the $\psi$ and $\chi$ fields do exactly
decouple, so that we are actually in presence of separate dynamics; we stress
that this is true not only for small oscillations, but at the full nonlinear
level as
well.
Let us briefly comment on these dispersion relations by remarking some features
of
theirs.
{\it a)} Now there is no mode with $\om (q) =0 $ for $q \to 0$; this means that
we
have a threshold for creation of excitations of given (small) amplitude.
{\it b)} The crossing of the two branches\fnote{As already mentioned, in
the PB model the two branches are non interacting, so that the crossing points
have
no physical meaning; the situation will be different in the context of the Y
model.},
and their relative ordering for $q=0$, are controlled by the sign of the
quantity $$
s = 1 - {\a \over 2 \g } $$
{\it c) } the $\psi$ branch has absolute minimum at $q=0$, while the $\chi$
branch
has an absolute minimum at $q_0 \not= 0$ ($q_0 \simeq 4 \pi / w$, which
corresponds
to a wavelenght $\l_0 \simeq w/2$). This latter fact is an obvious consequence
of
the introduction of the helicoidal interaction and of the signs in the $\chi$
part
of $V_H$: $\chi$ waves with $\chi_{n+h} = - \chi_n$ cost less energy.
\titleb{\tbn. Nondispersive waves and breather modes in the helicoidal
Peyrard-
Bishop model}
The waves corresponding to extremal points of the dispersion relations have
zero
group velocity; i.e., they are {\it standing waves} or, in other words,
correspond
to {\it breather modes}: there are oscillatory local openings of the double
helix
which remain localized without spreading or travelling.
The PB model is aimed at describing the denaturation of DNA; this process
arises
indeed by local opening of "bubbles" in the double helix, i.e. opening of
consecutive pairs of bases.
The experiments tell us that the size of the bubbles is of a few base units,
and
the time frequencies involved are of the order of picoseconds;
it is hoped that the breather modes of the hPB model can indeed model this
process.
We discuss here numerical findings\fnote{More
complete and refined numerical simulations of this model are being performed by
Peyrard and collaborators; these also concern estimates of parameters.} based
on parameter estimates given in [10]. It turns out
that the hPB model is rather well compatible with experimental findings.
Notice also that, due to the virial theorem, the energy of waves of amplitude
$A$
and obeying the dispersion relation $\om = \om (q)$ is just
$$ E = A^2 m \om^2 (q) \eqno(1) $$
Let us now write for ease of notation
$$ \xi = { qw \over 2} \eqno(2) $$
and rewrite the dispersion relatons, in terms of the physical coupling
parameters,
as
$$ \eqalign{
\Om_\psi (\xi ) \equiv & m \om^2_\psi ( q(\xi )) = 4 \g \[ B \xi^2 + C + \sin^2
(\xi
) \] \cr
\Om_\chi (\xi ) \equiv & m \om^2_\chi ( q(\xi )) = 4 \g \[ B \xi^2 + \cos^2
(\xi )
\] \cr } \eqno(3) $$
where the parameters $B,C$ are given by
$$ B = {\b \over h^2 \g } ~~;~~ C = {\a D \over 4 \g } \eqno(4) $$
The extremal points $\xi = \xi_0 $ of these are obtained by solving $\pa \Om /
\pa
\xi = 0$, namely
$$ 2 B \xi \mp \sin (2 \xi ) = 0 \eqno(5) $$
with the minus sign for the $\psi$ branch and the plus one for the $\chi$
branch.
Notice that necessarily $\xi_0 \le 1/B = h^2 \g /\b$.
The wavelenght corresponding to $\xi_0$ is simply given by
$$ \l (\xi_0 ) = { 2 \pi \over q (\xi_0 ) } = { \pi h \over \xi_0 } \d \eqno(6)
$$
(Notice that the above bound on $\xi_0$ ensures then $\l_0 \ge [( \pi \b ) / (h
\g )
] \d $; with the value of $\b / \g$ considered in this section, this also
ensures
$\l_0 >> \d$, as required for the validity of our continuum
approximation).
The time frequency $\om$ and period $T$ corresponding to $\xi_0$ are
$$ \eqalignno{
\om ( \xi_0 ) =& \sqrt{{\Om ( \xi_0 ) \over m }} & (7) \cr
T ( \xi_0 ) =& {2 \pi \over \om ( \xi_0 ) } & (8) \cr } $$
>From (1), and expecting each field mode to have an average energy $ = \unm k
T$,
we have ($T \simeq 37 ~^0 C$) the expected amplitude of oscillations
corresponding
to $\xi_0$:
$$ A ( \xi_0 ) = \sqrt{ { k T /2 \over \Om ( \xi_0 ) }} \eqno(9) $$
As for the parameters $\a , \b , \g , m$, we stick to the estimates given by
Dauxois [10]; adapting the notation to our present one, these are given in
table I (we report $\a^2 D$ rather than $\a$ since the former is the coefficent
of the linear term in the equations of motion).
The corresponding numerical results are reported, for the $\psi$ and $\chi$
branches, respectively in tables II and III. As anticipated, the orders of
magnitude that we obtain in this way for $\l$ and $T$ are in good agreement
with
experimental observations (notice that opening of the double helix correspond
to the $\psi$ mode).
$$ \matrix{
\a & \simeq & 0.8 \ea & \simeq & 1.3 \cdot 10^4 \ec \cr
\b & \simeq & 1.5 \ea & \simeq & 2.4 \cdot 10^4 \ec \cr
\g & \simeq & 0.5 \ea & \simeq & 8.0 \cdot 10^3 \ec \cr
m & \simeq & 300 ~ {\rm m.u.} & \simeq & 5 \cdot 10^{-22} ~ {\rm g} \cr} $$
\medskip
\centerline{\petit Table I - Parameter choice for hPB model}
\bigskip
$$\matrix{
\xi & 1.4 & 3.7 & 4.0 & 0. \cr
\Omega (\xi_0 ) & 8.5 ~ 10^{3} & 1.5 ~10^{4} & 2.2 ~10^{4} & 3.2 ~10^{4} \cr
\om (q_0 ) & 4.1 ~10^{12} & 5.5 ~10^{12} & 6.6 ~10^{12} & 8.0 ~10^{12} \cr
T (q_0 ) & 1.5 \ps & 1.1 \ps & 0.9 \ps & 0.8 \ps \cr
A (q_0 ) & 0.16 \A & 0.12 \A & 0.10 \A & 0.08 \A \cr
\lambda (q_0 ) & 11.2 ~\d & 4.2 ~\d & 3.9 ~\d & \infty \cr} $$
\medskip
\centerline{\petit Table II - Physical quantities ($\chi $ branch, hPB model)}
\bigskip
$$ \matrix{
\xi & 1.8 & 2.8 & 0 \cr
\Omega (q_0 ) & 6.8 ~10^4 & 5.9 ~10^4 & 2.6 ~10^4 \cr
\om (q_0 ) & 1.2 ~10^{13} & 1.1 ~10^{13} & 7.2 ~10^{12} \cr
T (q_0 ) & 0.5 \ps & 0.6 \ps & 0.9 \ps \cr
A (q_0 ) & 5.6 ~10^{-2} \A & 6.0 ~10^{-2} \A & 9.0 ~10^{-2} \A \cr
\lambda (q_0 ) & 8.7 ~\d & 5.6 ~\d & \infty \cr}$$
\medskip
\centerline{\petit Table III - Physical quantities ($\psi$ branch, hPB model)}
\bigskip
\titleb{\tbn. Introduction of helicoidal terms in the Yakushevich model}
Let us now consider the effects of introducing helicoidal interactions in the
Yakushevich model; here again we add to the Yakushevich hamiltonian
an harmonic potential modelling the helicoidal interaction; i.e. we have
$$ U_H = \unm \g \[ \( \th_{n+h}^{(1)} - \th_n^{(2)} \)^2 +\( \th_{n+h}^{(2)} -
\th_n^{(1)} \)^2 \] \eqno(1) $$
Here again we pass to normal coordinates $\psi , \chi$, with
$$ \th_n^{(1)} = \psi_n + \chi_n ~~;~~ \th_n^{(2)} = \psi_n - \chi_n \eqno(2)
$$
in terms of which (1) reads
$$ U_H = \unm \g \[ \( \psi_{n+h} - \psi_n \)^2 + \( \chi_{n+h} + \chi_n \)^2
\]
\eqno(3) $$
Notice here as well the different signs with which $\psi$ and $\chi$ enter in
this interaction.
The total hamiltonian of the helicoidal Yakushevich (or hY) model reads
therefore
$$ H = T + V_S + V_T + V_H \eqno(4) $$
where $T, V_S , V_T $ are as in the Y model and
$$ V_H = \sum_n \unm \g \[ \( \psi_{n+h} - \psi_n \)^2 + \( \chi_{n+h} +
\chi_n
\)^2 \] \eqno(5) $$
with $h$ the half-pitch of the helix in base units ($h=5$).
This hamiltonian gives as equations of motion
$$ \eqalign{
I {\ddot \psi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) + \g ( \psi_{n+h}
- 2
\psi_n + \psi_{n-h} ) - \a \sin \psi_n \cos \chi_n \cr
I {\ddot \chi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) - \g (
\psi_{n+h} - 2
\psi_n + \psi_{n-h} ) - \a \sin \chi_n ( \cos \psi_n - \cos \chi_n ) \cr }
\eqno(6) $$
and small oscillations are described by the equations
$$ \eqalign{
I {\ddot \psi_n} =& \b ( \psi_{n+1} - 2 \psi_n + \psi_{n-1} ) + \g ( \psi_{n+h}
- 2
\psi_n + \psi_{n-h} ) - \a \psi_n \cr
m {\ddot \chi_n} =& \b (\chi_{n+1} - 2 \chi_n + \chi_{n-1} ) - \g ( \psi_{n+h}
- 2
\psi_n + \psi_{n-h} ) \cr } \eqno(7) $$
In order to consider the continuum limit of these, we introduce again the
nonlocal
(linear) operators $\W_\pm$,
$$\W_\pm (\phi ) = \phi (x+w , t) \mp 2 \phi (x,t) + \phi (x-w,t) \eqno(8) $$
and write the full equations as
$$\eqalign{ \psi_{tt} =& \k \Lapl \psi + \r \W_+ (\psi ) - \s \sin \psi \cos
\chi \cr
\chi_{tt} =& \k \Lapl \chi - \r \W_- (\chi ) - \s \sin \chi ( \cos \psi - \cos
\chi
) \cr} \eqno(9) $$
where we have used
$$\k \equiv \b \d^2 / I ~~,~~\r \equiv \g / I ~~;~~ \s = \a / I ~~;~~ w = h \d
\eqno(10) $$
with $\d$ the distance between consecutive bases along the axis of the double
helix
($\d \simeq 3.4 \A$).
As for the small oscillation equations, with the same notations these are given
simply by
$$\eqalign{ \psi_{tt} =& \k \Lapl \psi + \r \W_+ (\psi ) - \s \psi \cr
\chi_{tt} =& \k \Lapl \chi - \r \W_- (\chi ) \cr} \eqno(11) $$
Notice that these are formally analogous to the one obtained for the hPB model,
so
that we will not discuss the dispersion relations for the hY model, but
just refer to the corresponding discussion of the hPB model.
We stress that in the case of the Y and hY models the fields $\psi , \chi$ are
decoupled at the level of small oscillations (so that we can still call these
"normal coordinates"), but are instead coupled through nonlinear terms,
contrary to
what happens in PB and hPB models. Since the $\psi$ and $\chi$ fields are
now interacting through nonlinear terms, the "crossing points" (i.e. the points
at which
the two branches meet) do have physical meaning, corresponding to resonances
between the two fields [6,13].
\titleb{\tbn. Solitons in the helicoidal Yakushevich model}
As already recalled, the Yakushevich and helicoidal Yakushevich models have
been
primarily developed to describe travelling solitons, which are supposed to play
a
central role in the transcription of DNA [11].
Therefore, we will concentrate on solitons, and not consider
standing nondispersive waves in the hY model. These are waves corresponding to
the
extremal points of the dispersion relations, and their physical interpretation
is
again as breather modes giving local pulsating and nonpropagating openings of
the
double helix, i.e. breaking of the H bonds pairing bases at same site on
different
helices. The numerical computations, presented elsewhere [12], show that indeed
experimental data on spontaneous denaturation of DNA are best fitted by the hPB
model.
The full field equations for the hY model were obtained to be
$$\eqalign{
\psi_{tt} =& \k \Lapl \psi + \r \W_+ (\psi ) - \s \sin \psi \cos \chi \cr
\chi_{tt} =& \k \Lapl \chi - \r \W_- (\chi ) - \s \sin \chi ( \cos \psi - \cos
\chi
) \cr} \eqno(1) $$
to which correspondhamiltonian $H$ and lagrangian $L$ given by
$$ \eqalign{
L =& {I \over \d } \int \unm \( \psi_t^2 + \chi_t^2 \) - \unm \k \[ ( \grad
\psi )^2
+ ( \grad \chi )^2 \] - W (\psi , \chi ) \dd x \cr
H =& {I \over \d } \int \unm \( \psi_t^2 + \chi_t^2 \) + \unm \k \[ ( \grad
\psi )^2
+ ( \grad \chi )^2 \] + W (\psi , \chi ) \dd x \cr } \eqno(2) $$
where now $W( \psi , \chi )$ is made up of a local and a nonlocal ($\N$) part:
$$ \eqalign{ & W (\psi , \chi ) = \unm \s \[ 1 + \cos^2 \chi - 2 \cos \psi \cos
\chi
\] + \N (\psi , \chi ) \cr
& \N ( \psi , \chi ) = \unm \r \[ \[ \psi (x+w) - \psi (x) \]^2 + \[ \chi
(x+w)
+ \chi (x) \]^2 \] \cr } \eqno(3) $$
The appearance of nonlocal terms in (1),(2) is the main new qualitative feature
encountered in the helicoidal version of the model. Physically, this correspond
to
the introduction of a second lenght scale in the model: beside the lattice
spaceing $\d$, we now also have introduced the (half) pitch of the helix, by
the
scale $w=h \d$.
We want now to concentrate on solitons, i.e. on solutions of the form
$$ \eqalign{ \psi (x,t) =& \psi (x - vt) \equiv \psi (z) \cr
\chi (x,t) =& \chi (x - vt) \equiv \chi (z) \cr } \eqno(4) $$
Notice that the discussion of boundary conditions conducted in the planar case
does
still apply here, as well as the related classification of soliton solutions by
topological (or soliton) integer numbers $(n,m)$.
Inserting the ansatz (4) into (2), we get the hamiltonian and the lagrangian
for
helicoidal Yakushevich solitons, where $\grad \equiv \pa / \pa z$:
$$ \eqalign{
H =& {I \over \d} \int \unm (v^2 + \k ) \[ (\grad \psi )^2 + (\grad \chi )^2 \]
+ W(\psi , \chi ) \dd z \cr
L =& {I \over \d} \int \unm (v^2 - \k ) \[ (\grad \psi )^2 + (\grad \chi )^2 \]
- W(\psi , \chi ) \dd z \cr } \eqno(5) $$
We write again $\~\mu = v^2 - \k = - \mu$, and the bounds on $v,\mu $ given in
the
planar case also apply here.
The field evolution equations are now\fnote{With more precise notation, one
should have $\d W / \d \psi$, $\d W / \d \chi$, as these are now functional
derivatives.}
$$ \eqalign{
\mu \psi_{zz} =& {\pa W \over \pa \psi } \cr
\mu \chi_{zz} =& {\pa W \over \pa \chi } \cr } \eqno(6) $$
By the same procedure as in the planar case, these give
$$ \unm \mu \[ ( \grad \psi )^2 + ( \grad \chi )^2 \] = W (\psi , \chi )
\eqno(7) $$
which again shows, once used in (5), that the energy of soliton solutions
$\=\psi ,
\=\chi$ is given by
$$ E = {I \k \over \d } \int \[ ( \grad \=\psi )^2 + ( \grad \=\chi )^2 \] \dd
z
\eqno(8) $$
\titleb{\tbn. Exact solution for the helicoidal (1,0) soliton}
The nonlocal terms make it quite difficult to obtain
explicit soliton solutions for the hY model; one can nevertheless get some
scaling
relations obeyed by these solutions,
which are also of use in order to obtain the explicit form of the simplest of
them.
The idea is to scale down the lenght $w = h \d$ to $\d$, which produces an
"effective planar lagrangian", and then apply the continuous approximation. The
limitation to this approach come from the special sign in the nonlocal
interaction of $\chi$ field, which prevents from obtaining an effective
lagrangian
of the Yakushevich type, and would introduce terms of the type $\chi \grad
\chi$
(see the discussion in [6]).
The limitations to work with an "effective planar lagrangian" for the rescaled
hY
model do not apply if we just consider the (1,0) solitons, due to $\chi \equiv
0$.
In this case we can obtain the explicit soliton solution just by the one
obtained
for the planar Y model, by scaling arguments.
Indeed, by scaling arguments [6,8], one proofs that by redefining once again an
effective mass
$$ \mu_H \equiv \mu - \r \d^2 \equiv (\k - \r ) \d^2 - v^2 \eqno(1) $$
the equations are formally analogous to the equation encountered in the planar
case;
i.e. we get
$$ \unm \mu_H ( \grad \psi )^2 = \a (1 - \cos \psi ) \eqno(2) $$
This can be integrated as shown earlier; notice that in these terms the
introduction of helicoidal interactions amounts to a rescaling of the constant
$A$
appearing in the explicit solution, which is now
$$ A = \sqrt{ { 2 \a \over \mu_H }} = \sqrt{{2 \a \over \mu - \r \d^2 }}
\eqno(3) $$
The dependence of (1,0) soliton energy on speed is given by
$$ E ( \l , \g ) = { A' \over B' - \g^2 } E ( \d , 0 ) \eqno(4) $$
where
$$ A' = 1 + { \g \over \b } ~~;~~B' = 1 + {\g \over h^2 \b } \eqno(5) $$
(see [6,8] for details). This means that the scaling of energy with speed is
qualitatively analogous to the one obtained for the planar Y model; similar
considerations apply to the relation of size versus speed of soliton
solutions.
It should be stressed that the (1,0) solitons correspond to those whose role in
DNA
transcription was first suggested by Englander [11], so that the discussion
of the present section is of particular relevance, although concerning only one
specific soliton solution.
\titleb{\tbn. Numerical results for solitons in the helicoidal Yakushevich
model}
For solitons other than the (1,0) ones, we are not able to determine explicit
analytical solutions, nor exact scaling laws at fixed $\l$; we must therefore
revert to numerical computations.
We fix the parameters $\a , \b , \g , I$ appearing in the model, and let $\mu$
vary between $0$ and
$\mu_{{\rm max}}$ (i.e. we consider different values of $\g$; it should be
recalled that physically
$\g \simeq 10^{-8}$). For any given value of $\mu$, one can determine
numerically, on a
sufficently large discrete lattice, the configuration $\=\psi , \=\chi$ which
makes
extremal - among those with prescribed boundary conditions - the Yakushevich
soliton
lagrangian (see eq.(2) of section 22).
Let $S(\mu )$ be this extremal value of the lagrangian,
$$ S(\mu ) = - { I \over \d} \int \unm \mu \[ ( \grad \=\psi )^2 + ( \grad
\=\chi
)^2 \] + W ( \=\psi , \=\chi ) \dd z \eqno(1) $$
The corresponding energy, by comparing the soliton lagrangian and hamiltonian
and
recalling that for soliton solutions
$$ \unm \mu \[ ( \grad \=\psi )^2 + ( \grad \=\chi )^2 \] = W ( \=\psi , \=\chi
)
\eqno(2) $$
is given simply by
$$ H ( \mu ) = {\mu_{{\rm max}} \over \mu } S (\mu ) \eqno(3) $$
In terms of $\g$, this reads
$$ H (\g ) = {1 \over 1 - \g^2 } S \( \mu ( \g ) \) \equiv {1 \over 1 - \g^2 }
S \(
(1 - \g^2 ) \mu_{{\rm max}} \) \eqno(4) $$
The minimization of the lagrangian is performed by a Montecarlo algorithm; the
values of the parameters are given in table IV and discussed in [8]; the
results of the numerical
computations are given in fig.\fgn for the (1,1) soliton (as recalled above,
the $\psi$ field is
supposed to matter for the solitons we are looking for: we want therefore a
nonzero $\psi$ soliton
number); here we plot $H(\g )$ in terms of $\g$, together with the best fit by
functions of
the form
$$ H(\g ) = { A \over B - \g^2 } \eqno(5) $$
This shows good qualitative agreement with the situation seen above for the
(1,0)
soliton, which is in turn qualitatively equivalent to the one for planar Y
solitons.
\titleb{Conclusions}
We have discussed the relevant degrees of freedom in DNA dynamics at body
temperature;
we have then considered in detail two simple models compatible with our general
discussion.
These give behaviours which are in quite good agreement with experimental
findings; in
particular the hPB model correctly predicts the size and time frequency order
of magnitude
for spontaneous openings of base pairs, while the Y and hY models have soliton
solutions as those
supposed to play a central role in DNA transcription, the size of these
matching the
experimental observations.
\titleb{Acknowledgements}
I would like to thank the organizers of the conference for inviting me in St.
Petersbourg and for giving me this opportunity to present my views on DNA
models. My interest in DNA nonlinear dynamics was raised by papers of, and kept
awake by correspondence with, prof. L. Yakushevich, which I would like to
thank here. I would also like to thank C. Reiss and M. Peyrard for the
discussions we had together and for permit to use part of a preliminary version
of our joint work [6] here.
\vfill \eject
\titleb{References} \parskip = 5 pt \parindent=10pt
{
\item{[1]} W. Saenger, {\it "Principles of Nucleic Acid Structure"}, Springer,
New York, 1984; second corrected printing 1988
\item{[2]} S. Yomosa, {\it Phys. Rev. A} {\bf 27} (1983), 2120; {\it Phys.
Rev. A} {\bf 30} (1984), 474
\item{[3]} M.Ya. Azbel, {\it J. Phys. A} {\bf 12} (1979), L29
\item{[4]} M. Peyrard and A.R. Bishop, {\it Phys. Rev. Lett.} {\bf 62} (1989),
2755
\item{[5]} L.V. Yakushevich, {\it Phys. Lett. A} {\bf 136} (1989), 413
\item{[6]} G. Gaeta, C. Reiss and M. Peyrard, "Nonlinear DNA dynamics"; in
preparation.
\item{[7]} A. Dubrovin, S.P. Novikov and A. Fomenko, {\it "Geometrie
Contemporaine I, II \& III"}, Mir,
Moscow, 1982 \& 1987; {\it "Modern Geometry I \& II"}, Springer, New York,
1984
\item{[8]} G. Gaeta, "Solitons in planar and helicoidal Yakushevich model";
preprint C.P.Th. Ecole Polytechnique, 1992
\item{[9]} G. Gaeta, {\it Phys. Lett. A} {\bf 143} (1990), 227
\item{[10]} T. Dauxois, {\it Phys. Lett. A} {\bf 159} (1991), 390
\item{[11]} S.W. Englander {\it et al.}, {\it Proc. Natl. Acad. Sci. USA} {\bf
77} (1980), 7222
\item{[12]} G. Gaeta, "An amended version of helicoidal models for DNA
nonlinear dynamics"; preprint C.P.Th. Ecole Polytechnique, 1992
\item{[13]} G. Gaeta, "Crossing points and resonance in DNA models"; in
preparation
}\bye