Influences of structure disorders and temperatures of systems on the bio-energy transport in the protein molecules (II)

Abstract

The influence of molecular structure disorders and physiological temperature on the states and properties of solitons as transporters of bio-energy are numerically studied through the fourth-order Runge-Kutta method and a new theory based on my paper [Front. Phys. China, 2007, 2(4): 469]. The structure disorders include fluctuations in the characteristic parameters of the spring constant, dipole-dipole interaction constant and exciton-phonon coupling constant, as well as the chain-chain interaction coefficient among the three channels and ground state energy resulting from the disorder distributions of masses of amino acid residues and impurities. In this paper, we investigate the behaviors and states of solitons in a single protein molecular chain, and in α-Helix protein molecules with three channels. In the former we prove first that the new solitons can move without dispersion, retaining its shape, velocity and energy in a uniform and periodic protein molecule. In this case of structure disorder, the fluctuations of the spring constant, dipole-dipole interaction constant and exciton-phonon coupling constant, as well as the ground state energy and the disorder distributions of masses of amino acid residues of the proteins influence the states and properties of motion of solitons. However, they are still quite stable and are very robust against these structure disorders, even in the presence of larger disorders in the sequence of masses, spring constants and coupling constants. Still, the solitons may disperse or be destroyed when the disorder distribution of the masses and fluctuations of structure parameters are quite great. If the effect of thermal perturbation of the environment on the soliton in nonuniform proteins is considered again, it is still thermally stable at the biological temperature of 300 K, and at the longer time period of 300 ps and larger spacing of 400 amino acids. The new soliton is also thermally stable in the case of motion over a long time period of 300 ps in the region of 300−320 K under the influence of the above structure disorders. However, the soliton disperses in the case of a higher temperature of 325 K and in larger structure disorders. Thus, we determine that the soliton’s lifetime and critical temperature are 300 ps and 300−320 K, respectively. These results are also consistent with analytical data obtained via quantum perturbed theory. In α-helix protein molecules with three channels, results obtained show that these structure disorders and quantum fluctuations can change the states and features of solitons, decrease their amplitudes, energies and velocities, but they still cannot destroy the solitons, which can still transport steadily along the molecular chains while retaining energy and momentum when the quantum fluctuations are small, such as in structure disorders and quantum fluctuations of , \text{Δ}W=\pm 8\%\overline{W}, \text{Δ}J=\pm 1\%\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=\pm 3\%({\overline{\chi }}_1+{\overline{\chi }}_2) and \text{Δ}L=\pm 1\%\overline{L} and Δε₀ = ε|β_n|, ε=0.1 meV, |β_n| < 0.5. Therefore, the solitons in the improved model are quite robust against these disorder effects. However, the solitons may be dispersed or disrupted in cases of very large structure disorders. When the influence of temperature on solitons is considered, we find that the new solitons can transport steadily over 333 amino acid residues in the case of motion over a long time period of 120 ps, and can retain their shapes and energies to travel forward along protein molecules after mutual collision of the solitons at the biological temperature of 300 K. Therefore, the soliton is also very robust against thermal perturbation of the α-helix protein molecules at 300 K. However, the soliton disperses in cases of higher temperatures at 325 K and in larger structure disorders. Thus, their critical temperature is about 320 K. When the effects of structure disorder and temperature are considered simultaneously, the soliton has high thermal stability and can transport for a long time along the protein molecular chains while retaining its amplitude, energy and velocity, even though the fluctuations of the structure parameters and temperature of the medium increase continually. However, the soliton disperses in the larger fluctuations of , \text{Δ}({\chi }_1+{\chi }_2)=\pm 2\%({\overline{\chi }}_1+{\overline{\chi }}_2), \text{Δ}J=\pm 1.3\%\overline{J}, \text{Δ}W=\pm 6\%\overline{W}, \text{Δ}L=\pm 1.5\%\overline{L} and Δε₀ = ε |β_n|, ε=0.82 meV, \left|{\beta }_n\right|⩽0.5 at T =300 K, and at temperatures higher than 315 K when the fluctuations are , \text{Δ}({\chi }_1+{\chi }_2)=\pm 1\%({\overline{\chi }}_1+{\overline{\chi }}_2), \text{Δ}J=\pm 0.7\%\overline{J}, \text{Δ}W=\pm 7\%\overline{W}, \text{Δ}L=\pm 0.8\%\overline{L} and Δε₀ = ε |β_n|, ε=0.4 meV, \left|{\beta }_n\right|⩽0.5. This means that the critical temperature of the soliton is only 315 K in this condition. In a word, we can conclude from the above investigations that the soliton in the improved model is very robust against the structure disorders and thermal perturbation of proteins at the biological temperature of 300 K in α-helix protein molecules, and is a possible bio-energy transport carrier; the improved model is a possible candidate for the mechanism of this transport.

Graphical abstract

Keywords

protein molecule / biological energy / soliton stability / structure disorder / temperature / exciton / Runge-Kutta way / numerical simulation

1 Introduction

We know that bio-energy transport is a fundamental process in living systems, and a lot of biological phenomena, for example, muscle contraction, DNA replication, neural information transfer along cell membranes, and the work of sodium and calcium pumps, are all associated with this process, where the energy is released by hydrolysis of adenosin triphosphate (ATP). Therefore, to study bio-energy transport in protein molecules is a very important problem in biology. However, understanding the mechanism of energy transport is a long-standing problem that remains of great interest up to now. One can assume that the energy is stored as vibrational energy in the C=O stretching mode (amide-I) of α-helix protein molecules with three channels, as shown in Fig. 1. Thus, Davydov first proposed a theory of bio-energy transport in proteins in the 1970s [2−4]. Following Davydov’s idea [2−4], the mechanism of energy transport can be described as follows: The stretching vibration of C=O in the proteins results in the deformation of amino acid molecules; a coupling interaction between the amide-I vibrational quantum (exciton) and a deformation of the amino acid lattice occurs; through the nonlinear coupling interaction, the exciton is then self-trapped as a soliton [2−4]. It, together with the deformed lattice, can move over macroscopic distances along the protein molecular chains, retaining its wave shape, energy and momentum, as well as the other properties of a quasi-particle. This is just the Davydov model for bio-energy transport which he first proposed in α-helix protein molecules in Fig. 1 in the 1970s [1]. Davydov’s idea yields a compelling picture of the mechanism of bio-energy transport in protein molecules and consequently has been the subject of a large number of work [2−16]. Problems related to the Davydov model, including the foundation and accuracy of his theory, the quantum and classical properties, and the thermal stability and lifetime of the Davydov soliton, have been extensively studied by many scientists [5−35]. Further, considerable controversy has arisen in recent years over whether the Davydov soliton is sufficiently stable in the region of biological temperatures to provide a viable explanation for bio-energy transport. Many numerical simulations [13−34] have been based essentially on classical motion and are subject to the criticism that they are likely to yield equations with unreliable estimates for the stability of the soliton since the dynamics of the soliton are not determined by the Schr¨odinger equation [6]. For the thermal equilibrium properties of the Davydov soliton there is quantum Monte Carlo simulation [31, 32]. In the simulation, correlation characteristics of soliton-like quasi-particles occur only at low temperatures (about T <10 K) for widely accepted parameter values. This is consistent at a qualitative level with the result of Cottingham et al. [33, 34]. The latter is a straightforward quantum-mechanical perturbation calculation, in which the lifetime of the Davydov soliton obtained is too small (about 10⁻¹²−10⁻¹³ s) to be useful in biological processes. This shows clearly that the Davydov soliton is not a true wave function of biological systems. It is then necessary to reform Davydov’s wave function. Scientists have thought that the soliton with a multiquantum state (n > 2), for example, Brown et al.’s coherent state [13−17], and Kerr et al.’s [29, 30] and Schweitzer et al.’s [33, 34] multiquantum state, and Cruzeiro–Hansson’s [22−24] and F¨orner’s [25−28] two quantum state, and so on, would be thermally stable in the region of biological temperatures, and could provide a realistic mechanism for the bio-energy transport in the protein molecules. However, the assumption of the standard coherent state is unsuitable or impossible for protein molecules because the number of particles in this state is innumerable and one cannot retain conservation of the number of particles of the system. The assumption of a multiquantum state (n > 2) along with a coherent state also entails that the energy released in ATP hydrolysis (about 0.43 eV) can excite only two quanta of amide-I vibration. On the other hand, Cruzeiro-Hansson [9] thought that F¨orner’s two-quantum state in the semiclassical case was not exact. We prove that the Cruzeiro-Hansson’s ansatz contains four excitons (quanta) [17], instead of two excitons. Obviously, it is not possible at all to create the four excitons from the energy released in ATP hydrolysis alone.

Based on the work of Cruzeiro-Hansson and F¨orner, and so on, we improve and develop [36−76] the Davydov model by simultaneously changing the Hamiltonian and the wave function of the system. We add a new coupling interaction between the acoustic phonon and amide-I vibrational modes into the original Davydov’s Hamiltonian, and replace the one-quantum (exciton) state in Davydov’s wave function with a quasi-coherent two-quantum state. Thus, the equation of motion and the properties of the soliton excited in the improved model are completely different from that in the Davydov model. I suppose that this model could resolve the controversy on the thermal stability and lifetime of the soliton excited in protein molecules. In our theory, the wave function and Hamiltonian of the protein molecules with one-channel is represented by [36−76]

\begin{array}{}\text{(1)}& \begin{array}{ll}\hfill & |\Phi (t)\rangle =|a(t)\rangle |\beta (t)\rangle \hfill \\ \hfill & =\frac{1}{\lambda }\left[1+{\displaystyle \sum _n}{a}_n(t)B_n^{+}+\frac{1}{2!}{\left({\displaystyle \sum _n}{a}_n(t)B_n^{+}\right)}^2\right]|0{\rangle }_{ex}\hfill \\ \hfill & \cdot \exp \left\{-\frac{\text{i}}{\hslash }{\displaystyle \sum _n}\left[{\beta }_n(t)P_n-{\pi }_n(t)u_n\right]\right\}|0{\rangle }_{ph}\hfill \end{array}\end{array}

and

\begin{array}{}\text{(2)}& \begin{array}{ll}\text{H}\hfill & =H_{\text{ex}}+H_{\text{ph}}+H_{\text{int}}\hfill \\ \hfill & ={\displaystyle \sum _n}[{\epsilon }_0{B}_n^{+}{B}_n-J(B_n^{+}{B}_{n+1}+B_n{B}_{n+1}^{+})]\hfill \\ \hfill & +{\displaystyle \sum _n}\left[\frac{P_n^2}{2M}+\frac{1}{2}W\cdot {(u_n-u_{n-1})}^2\right]\hfill \\ \hfill & +{\displaystyle \sum _n}[{\chi }_1(u_{n+1}-u_{n-1})B_n^{+}{B}_n\hfill \\ \hfill & +{\chi }_2(u_{n+1}-u_n)(B_{n+1}^{+}{B}_n+B_n^{+}{B}_{n+1})]\hfill \end{array}\end{array}

where B_n^{+} and B_n are the Boson creation and annihilation operators for the exciton, respectively, and {|0\rangle }_{\text{ex}} and |0{\rangle }_{\text{ph}} are the ground states of the exciton and phonon, respectively; u_n and P_n are the displacement and momentum operators of the amino acid residue in site n, respectively; the a_n(t)\mathrm{,}{\beta }_n(t)=\langle \Phi \left|u_n\right|\Phi (t)\rangle and {\pi }_n(t)=\langle \Phi (t)\left|P_n\right|\Phi (t)\rangle are three sets of unknown functions, and λ is a normalization constant, where {\epsilon }_0=\hslash {\omega }_0=1665{\text{cm}}^{-1}=0.2035\text{eV} is the excitation energy of an isolated amide-I oscillator or energy of the exciton (the C = 0 stretching mode). The present nonlinear coupling constants are χ₁ and χ₂, which represent the modulations of the one-site energy and resonant (or dipole-dipole) interaction energy for the excitons caused by the displacement of the amino acid residue, respectively. M is the mass of an amino acid residue, W is the elasticity constant of the protein, and J is the dipole-dipole interaction energy between neighbouring amino acids, while |0{\rangle }_{\text{ph}} and {|0\rangle }_{\text{ex}} are vacuum states of the phonon and exciton, respectively. H_ex here describes Boson-type Frenkel excitons excited by the energy released from ATP hydrolysis in the protein molecules, H_ph describes a harmonic feature of amino acid residues, H_int represents the interaction between the two modes of motion. Usually for all parameters in Eqs. (1) and (2) site-independent mean values are used. The average value of the dipole-dipole coupling between neighboring amide-I oscillators is \overline{J}\text{=0}\text{.967}\text{meV}. The average spring constant of the hydrogen bonds is taken usually to be \overline{W}=13\text{ N/m}. The average mass \overline{M} is taken as that of myosine (\overline{M}=114m_{\text{p}}, m_p is a proton mass). For {\overline{\chi }}_1 the experimental value is 62 pN, {\overline{\chi }}_2=10-15\text{pN} for the protein molecules with one channel.

The Hamiltonian and wave function shown in Eqs. (1) and (2) are different from Davydov’s Haniltonian and wave function. We have added a new interaction term, \begin{array}{l}{\displaystyle \sum _n}{\chi }_2(u_{n+1}-u_n)(B_{n+1}^{+}{B}_n+B_n^{+}{B}_{n+1})\hfill \end{array}, into the original Davydov’s Hamiltonian. Thus, the Hamiltonian now has better correspondence for the interactions and can also represent features of the mutual correlations of the collective excitations and motions of quasi-particles in the protein molecules. The present wave function of the exciton in Eq. (1) is not an excitation state of a single-particle, but a coherent state, or more accurately, a quasi-coherent state. It retains only three terms of the expansion of a standard coherent state, which is mathematically justified in the case of small a_n(t) (i. e., |a_n(t)|\ll 1), and it can be viewed as an effective truncation of a standard coherent state. Therefore, we call |a(t)\rangle a a quasi-coherent state. However, it is not an eigenstate of the number operator, \widehat{N}={\text{Σ}}_n{B}_n^{+}{B}_n, but a coherent superposition of the excitonic state with two quanta and the ground state of the exciton. However, in this state the number of quanta is determinate, instead of innumerable. To find out how many excitons this state contains, we have to compute the expectation value of the number operator \widehat{N} in this state and sum over the sites. The average number of excitons for this state is

\begin{array}{}\text{(2)}& N=\langle a(t)\left|\widehat{N}\right|a(t)\rangle ={\displaystyle \sum _n\langle a(t)|B_n^{+}{B}_n|a(t)\rangle =2}\end{array}

Therefore, it really contains two excitons. Thus, we see that the improved model is completely different from the Davydov model, and the results obtained from this model are fully different from that of the Davydov model. The distinctions between the two models are shown in Table 1 [36−55]. From Table 1 we know that the improved model repulses and avoids the shortcoming of the Davydov model [2−4]; the new soliton in the model is thermally stable at the biological temperature of 300 K and has a long enough lifetime, thus it can play an important role in biological processes, i.e., the new soliton is a very good carrier for bio-energy transport in protein molecules of living systems.

However, the above results were obtained analytically, in which the protein molecules are single chains and thought to be periodic systems; all physical parameters of the protein molecules were used in terms of their average values, and some approximation methods, containing long-wave approximation, continuum approximation, and so on, were used in the calculations. In practice, the biological proteins consist of 20 different amino acid residues with molecular weights between 75 m_p(glycine) and 204 m_p (tryptophan), which correspond to variations in mass between , where \overline{M}=114m_{\text{p}} is the average mass of an amino acid residue and m_p is the proton mass. Therefore, the protein molecules are not periodic but disordered or nonuniform in structure. This structural aperiodicity necessarily results in fluctuations of the spring constant, the dipole-dipole interaction, the exciton-phonon coupling constant, the diagonal disorder and the chain-chain interaction among the molecular channels. In α-Helix protein molecules with three channels, there are also chain-chain interactions among the three channels, which carry a dispersive effect and further influence the states of the new solitons. Thus, the states of new solitons will be changed because Careri et al.’s investigations [77, 78] appear to indicate that even relatively small amounts of disorder in an amorphous film of acetanilide (ACN), a protein-like crystal (i.e., molecular structure of acetanilide crystal is quite analogous to an α-helix protein), is enough to destroy the spectral signature of a “soliton”. In such a case, it is very necessary to study the properties of the soliton in the bio-energy transport process in disordered protein molecules. In this paper we will study the influence of structure disorders and temperatures of systems on the features of solitons excited in a single protein chain as well as in α-helix protein molecules with three channels (shown in Fig. 1) via numerical simulation and the fourth-order Runge-Kutta method [79, 80] in which no approximation methods are used. We will see that the soliton is still very robust against such structural disorders of the protein molecules. In Section 2 we first investigate the states and behaviors of the solitons in a single protein molecular chain. In Section 3 we further study the properties and states of the solitons in the α-Helix protein molecules with three channels. In Section 4 we state the conclusions of this paper.

2 The states and behaviors of the solitons in proteins with single channels

2.1 Calculation method

Utilizing Eqs. (1), (2) and from the following Schr¨odinger equation and Heisenberg equation:

\begin{array}{}\text{(3)}& \text{i}\hslash \frac{\partial }{\partial t}|\Phi (t)\rangle =H|\Phi |(t)\rangle \end{array}

and

\begin{array}{}\text{(4)}& \text{i}\hslash \frac{\partial }{\partial t}\langle |\Phi (t)|u_n|\Phi (t)\rangle =\langle \Phi (t)|[u_n\mathrm{,}H]|\Phi (t)\rangle \end{array}

\begin{array}{}\text{(5)}& \text{i}\hslash \frac{\partial }{\partial t}\langle \Phi (t)|P_n|\Phi (t)=\langle \Phi (t)|[P_n\mathrm{,}H]|\Phi (t)\rangle \end{array}

we can obtain

\begin{array}{}\text{(6)}& \begin{array}{ll}\text{i}\hslash {\dot{a}}_n(t)=\hfill & {\epsilon }_0{a}_n(t)-J[a_{n+1}(t)+a_{n-1}(t)]\hfill \\ \hfill & +{\chi }_1[q_{n+1}(t)-q_{n-1}(t)]a_n(t)+{\chi }_2[q_{n+1}(t)\hfill \\ \hfill & -q_n(t)][a_{n+1}(t)+a_{n-1}(t)]+\frac{5}{2}\{w(t)\hfill \\ \hfill & -\frac{1}{2}{\displaystyle \sum _m}{q}_m(t){\pi }_m(t)-{\dot{\pi }}_m(t){\dot{q}}_m(t)]\}{a}_n(t)\hfill \end{array}\end{array}

\begin{array}{}\text{(7)}& \begin{array}{ll}M{\ddot{q}}_n(t)=\hfill & W[q_{n+1}(t)-2q_n(t)+q_{n-1}(t)]\hfill \\ \hfill & +2{\chi }_1[|a_{n+1}{|}^2-|a_{n-1}{|}^2]\hfill \\ \hfill & +2{\chi }_2\{a_n^{*}(t)[a_{n+1}(t)-a_{n-1}(t)]\hfill \\ \hfill & +a_n(t)[a_{n+1}^{*}(t)-a_{n-1}^{*}(t)]\}\hfill \end{array}\end{array}

using transformation: a_n(t) → a_n exp (iε₀t/h) we can eliminate the term ε₀a_n(t) in Eq. (6). Again making a transformation: a_n(t) = a_n(t)r_n + ia(t)i_n, then Eqs. (6) and (7) become

\begin{array}{}\text{(8)}& \begin{array}{ll}\hslash \dot{a}{r}_n=\hfill & -J(a{i}_{n+1}+a{i}_{n-1})+{\chi }_1(q_{n+1}-q_{n-1})a{i}_n\hfill \\ \hfill & +{\chi }_2(q_{n+1}-q_n)(a{i}_{n+1}+a{i}_{n-1})\hfill \end{array}\end{array}

\begin{array}{}\text{(9)}& \begin{array}{ll}-\hslash \dot{a}{i}_n=\hfill & -J(a{r}_{n+1}+a{r}_{n-1})+{\chi }_1(q_{n+1}-q_{n-1})a{r}_n\hfill \\ \hfill & +{\chi }_2(q_{n+1}-q_n)(a{r}_{n+1}+a{r}_{n-1})\hfill \end{array}\end{array}

\begin{array}{}\text{(10)}& {\dot{q}}_n=y_n/M\end{array}

\begin{array}{}\text{(11)}& \begin{array}{ll}{\dot{y}}_n=\hfill & W(q_{n+1}-2q_n+q_{n-1})\hfill \\ \hfill & +2{\chi }_1(a{r}_{n+1}^2+a{i}_{n+1}^2-a{r}_{n-1}^2-a{i}_{n-1}^2)\hfill \\ \hfill & +4{\chi }_2[a{r}_n(a{r}_{n+1}-a{r}_{n-1})\hfill \\ \hfill & +a{i}_n(a{i}_{n+1}-a{i}_{n-1})]\hfill \end{array}\end{array}

\begin{array}{}\text{(12)}& {|a_n|}^2={|a{r}_n|}^2+{|a{i}_n|}^2\end{array}

where ar_n and ai_n are real and imaginary parts of a_n. The solutions of Eqs. (8)−(11) are listed in the Appendix.

The above equations can determine states and behaviors of the new soliton. There are four equations for one peptide group [79, 80]. Therefore, for the protein molecules constructed by N amino acids there are 4N associated equations. When the fourth-order Runge-Kutta method is used to numerically calculate the solutions of the above equations we should discretize them. Thus the n is replaced by j and let the time be denoted by n, the step length of the space variable is denoted by h in the above equations. The system of units, eV for energy, Å for length, and ps for time, are proven to be suitable for the numerical solutions of Eqs. (6) and (7). For a numerical simulation obtained by the fourth-order Runge-Kutta method [18], we require that the following conditions must be satisfied: the total energy E=\langle \Phi (t)\left|H\right|\Phi (t)\rangle =\text{constant} (up to 0.0012 %); a possible imaginary part of the energy which can occur due to numerical inaccuracies is zero to an accuracy of 0.001 feV;and the norm is conserved up to 0.3pp (parts per million). An initial excitation is required in this calculation, and it is chosen as, a_n(0) = A sech [(n−n₀) (χ₁ + χ₂)²/(4JW)](where A is the normalization constant) at the size n, for the applied lattice, q_n(0) = π_n(0)=0. The molecular chain is fixed, N is chosen to be N = 50, and a time step size of 0.0195 is used in the simulations. Total numerical simulation is performed through data parallel algorithms and MATLAB language.

2.2 Influence of structure disorder on the soliton

2.2.1 Numerical results for free motion in the uniform periodic chains

In utilizing the above average values for the parameter, \overline{M}, \overline{W}, {\overline{\epsilon }}_0, \overline{J}, {\overline{\chi }}_1 and {\overline{\chi }}_2, we calculate numerically the solution of Eqs. (8)−(11) by using the fourth-order Runge-Kutta method [79, 80] in uniform and periodic proteins, where |a_j|² is the probability of the soliton occurring at the jth amino acid molecule. Thus, we can plot the state of the soliton in time-place. The result is shown in Fig. 2. This figure shows that the amplitude of the solution can retain constancy. In Figs. 3 and 4 we show the propagation behavior of the solution for a long time period of 300 ps, and the collision property of two solitons, respectively. From the figures we see that the solution is very stable while in motion for a long time period. Therefore, Eqs. (9) and (10) have exactly soliton solutions in uniform and periodic proteins.

2.2.2 Influence of mass disorders and the spring constant on the new soliton

As mentioned above, the biological protein consists of 20 different amino acid molecules with weights between 75 m_p (glycine) and 204 m_p (tryptophane), which correspond to weight variations between 0.67\overline{M} and 1.80\overline{M}. Therefore it is not a periodic, but an aperiodic system, and there is a disordered distribution in mass sequences of the amino acids. In order to study the influence of a random series of masses in the proteins on the new soliton, here we introduce a parameter α_k, representing the size of the disorder distribution, which is a random-number generator with equal probability within a prescribed interval and can denote the mass at each point in the molecular chain, i.e. M_k={\alpha }_k\overline{M}. Numerical simulation shows that when the α_k is small, for example, a natural \begin{array}{l}0.67\overline{M}⩽M_n⩽1.80\overline{M}\hfill \end{array}, the state of the new soliton is not influenced. Up to the α_k intervals of, for example, 0.67⩽{\alpha }_k⩽300, the stability of the new soliton may still remain, but in the case of large intervals such as 0.67⩽{\alpha }_k⩽700, the vibrational energy is dispersed, as shown in Fig. 5. The interval of 0.67⩽{\alpha }_k⩽300, over which the motion of the new soliton is unperturbed, is evidently larger than the above natural interval of masses of the amino acids. Thus, the new soliton is very robust against mass disorder in the proteins.

However, in the above calculation we do not consider the changes in structural features of proteins arising from the disorder of masses of amino acids, which will in turn result in changes of the structural parameters, for example, the spring constant W, dipole-dipole interaction constant J, coupling constant (χ₁ + χ₂), and the ground state energy ε₀ in Eqs. (4) and (5). According to Forner’s method [18, 19] the changes in parameters are represented by their fluctuations of average increments, \text{Δ}W=W-\overline{W}, \text{Δ}J=J-\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=({\chi }_1+{\chi }_2)-({\overline{\chi }}_1+{\overline{\chi }}_2) and Δε₀ = ε₀−ε₀, respectively, where W, J, (χ₁ + χ₂) and ε₀ are the values of the parameters in the protein molecules without the structural disorders. However, for the variation in ground state energy arising from imported impurities or from the mass disorder, we usually use the random number generator,|βn|, to designate its random feature, i.e., the Δε₀ is denoted by Δε₀ = ε |β_n|. In the following we study the collective effects of these fluctuations on the new soliton along this idea of research.

In Fig. 6 we show the changes in stability of a new soliton, with increasing fluctuations of the spring constant W in the case of a mass interval of 0.67⩽{\alpha }_k⩽2. We see that up to a random variation of \pm 40\%\overline{W}, the dynamics of the soliton have not changed. For \pm 50\%\overline{W}, the soliton disperses somewhat, but its velocity is only somewhat diminished, when compared to the earlier \overline{W} case \overline{W}. Finally, for \pm 70\%\overline{W}, the soliton disperses and its propagation is irregular, as shown in Fig. 6.

2.2.3 Influence of disorder of mass sequence and fluctuations of interaction constant J

The soliton in the improved model is, in general, more sensitive to the variations in the dipole-dipole interaction constant J caused by the structural disorder when compared with the other parameters. The simulation shows that for a variation in J alone, the soliton is stable up to \pm 9\%\overline{J}, and it disperses at \text{Δ}J=\pm 15\%\overline{J}. If we simultaneously consider the collective effects of disorder of mass sequence and fluctuation of J on the new soliton, then its state is obviously changed. In Fig. 7 we show the collective effect of the fluctuations \text{Δ}J=\pm 5\%\overline{J}, =\pm 10\%\overline{J}, =\pm 15\%\overline{J}, =\pm 20\%\overline{J}, on the soliton in the case of a mass interval of 0.67⩽{\alpha }_k⩽2. From these figures we see clearly that the new soliton is stable at \text{Δ}J⩽9\%\overline{J}, but it disperses at \text{Δ}J=\pm 15\%\overline{J}, and disperses significantly at \text{Δ}J=\pm 20\%\overline{J}.

2.2.4 Variations of the states of the new solitons resulting from the mass disorder and variations of coupling constant (χ₁ + χ₂)

The numerical calculation shows that, arising from the disorder of structure, if the coupling constant (χ₁ + χ₂) alone is changed, (χ₁ + χ₂) can be varied to ±25% ({\overline{\chi }}_1+{\overline{\chi }}_2), and in this case the new soliton does not disperse. However, for a fluctuation together with natural mass variation, the stability of the new soliton will be changed. In Fig. 8 we illustrate the changes in the states of new solitons with increasing fluctuations of (χ₁ + χ₂) at 0.67\overline{M}⩽M⩽2\overline{M}. We see from this figure that only at are the new solitons stable, but they obviously disperse at \text{Δ}({\chi }_1+{\chi }_2)=35\%({\overline{\chi }}_1+{\overline{\chi }}_2).

2.2.5 Synthetic influence of mass disorders and fluctuations of the spring constant and interaction constant

We also study the collective effect of fluctuations of the spring constant and dipole-dipole interaction on the stability of a new soliton in the case of a mass interval of 0.67⩽{\alpha }_k⩽2. In Fig. 9 we show the changes in stability of the new soliton with an increasing spring constant at \text{Δ}J=\pm 10\%\overline{J}. From this figure we see that up to a random variation of \pm 30\%\overline{W}, the dynamics of the new soliton can be maintained. For \pm 40\%\overline{W}, the soliton disperses somewhat, and its velocity is somewhat diminished. Finally, for \pm 50\%\overline{W}, the soliton disperses and its propagation becomes irregular. In Fig. 10 we plot the evolution of new solitons with an increasing spring constant in the case of a mass interval of 0.67⩽{\alpha }_k⩽2 and \text{Δ}J=\pm 15\%\overline{J}. We know from this figure that the stability of the new soliton is changed when compared with Fig. 9. The soliton is only stable at a random variation of \pm 20\%\overline{W}, disperses at \pm 30\%\overline{W}, and is destroyed at \pm 50\%\overline{W}. Obviously, this is due to the increasing fluctuations of the dipole-dipole interactions in the proteins.

2.2.6 Combined influence of the mass disorder and changes in the coupling constant and spring constant on the new solitons

In Fig. 11 we show variations in the properties of the new solitons arising from mass disorder and fluctuations of coupling and spring constants. From this figure we find that the soliton is only stable at fluctuations of \text{Δ}W=\pm 30\%\overline{W} and \text{Δ}({\chi }_1+{\chi }_2)=20\%({\overline{\chi }}_1+{\overline{\chi }}_2) in the case of a mass interval of 0.67⩽{\alpha }_k⩽2. With an increasing spring constant the soliton disperses gradually. Up to \text{Δ}({\chi }_1+{\chi }_2)=25\%({\overline{\chi }}_1+{\overline{\chi }}_2), we find that the new soliton is still stable, at 0.67⩽{\alpha }_k⩽2 and \text{Δ}W=\pm 10\%\overline{W}, but it disperses at 0.67⩽{\alpha }_k⩽2 and \text{Δ}W=\pm 40\%\overline{W}, and is disrupted at 0.67⩽{\alpha }_k⩽2 and \text{Δ}W=\pm 70\%\overline{W} and \begin{array}{l}\text{Δ}({\chi }_1+{\chi }_2)=20\%({\overline{\chi }}_1+{\overline{\chi }}_2)\hfill \end{array}. Meanwhile, its propagation is irregular at 0.67⩽{\alpha }_k⩽2 and \text{Δ}W=\pm 50\%\overline{W}. These results are shown in Fig. 12.

2.2.7 Collective effects of mass disorder and fluctuations of the spring constant, interaction constant and coupling constant on the new soliton

It is quite necessary to collect the combined effect of random variations of the above physical parameters resulting from structural disorders on the properties of the new soliton. The changes of states of solitons with increasing interaction constants in the case of 0.67⩽{\alpha }_k⩽2 and \text{Δ}({\chi }_1+{\chi }_2)=10\%({\overline{\chi }}_1+{\overline{\chi }}_2) and \text{Δ}W=\pm 20\%\overline{W} are plotted in Fig. 13. We see from this figure that the soliton is stable at \text{Δ}J=\pm 2.5\%\overline{J}, begins to disperse at \text{Δ}J=\pm 10\%\overline{J}, obviously disperses at \text{Δ}J=\pm 20\%\overline{J}, and is destroyed at \text{Δ}J=\pm 25\%\overline{J}. We show the states of the solitons for different ΔJ and ΔW in the case of 0.67⩽{\alpha }_k⩽2 and \text{Δ}J=\pm 10\%\overline{J} in Fig. 14. We find no change in the dynamics of the new soliton at \text{Δ}W=\pm 25\%\overline{W} and \text{Δ}J=\pm 10\%\overline{J}, but it begins dispersing at ΔW = \pm 30\%\overline{W}, \text{Δ}J=\pm 15\%\overline{J}, and disperses considerably at \text{Δ}W=\pm 25\%\overline{W}, \text{Δ}J=\pm 15\%\overline{J}.

2.2.8 Combined effects of mass disorder and diagonal disorder, and fluctuations of the spring constant, interaction constant and the coupling constant on the new soliton

We are more interested in the collective effect of simultaneous random-variations of the above five parameters resulting from structural disorders on the new soliton. In general, the new soliton is very sensitive to the diagonal disorder, which is the change in ground state energy, Δε₀, caused by different side groups of amino acids and local geometric distortions due to the impurities imported. We found that for an isolated impurity in the middle of the chain, which causes the change in the energy to be Δε₀ = εδ_n, the soliton can pass the impurity only if ε <1 meV. In other cases it is re-flected or dispersed. In the case of a random sequence, \text{Δ}{\epsilon }_0=\epsilon \left|{\beta }_n\right|\mathrm{,}\left|{\beta }_n\right|⩽0.5, β_n is a random parameter, and only when ε < 1 meV can the soliton pass the impurity to propagate along the chains. For higher values of ε the soliton is dispersed. When the diagonal disorder occurs together with fluctuations of the other four parameters, the state of the new soliton is changed obviously. The states of solitons with an increasing spring constant in the case of \begin{array}{l}0.67⩽{\alpha }_k⩽\mathrm{2,}\text{Δ}({\chi }_1+{\chi }_2)=5\%({\overline{\chi }}_1+{\overline{\chi }}_2)\hfill \end{array}, \text{Δ}J=\pm 5\%\overline{J} and meV are given in Fig. 15. We see from this figure that when the five parameters are all randomly varied, a maximal disorder that would occur still with the soliton in motion, is \text{Δ}W=\pm 10\%\overline{W}, \text{Δ}J=\pm 5\%\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=\pm 5\%({\overline{\chi }}_1+{\overline{\chi }}_2), 0.67\overline{M}⩽M⩽2\overline{M} and \text{Δ}{\epsilon }_0=\epsilon \left|{\beta }_n\right|\mathrm{,}\left|{\beta }_n\right|⩽\mathrm{0.5,}\epsilon =1 meV. In other cases the soliton is dispersed or reflected by the impurity. This phenomenon is observed in Fig. 16. In this figure we exhibit clearly the changes in state of the new soliton due to simultaneous variations in the above five parameters. However, the soliton is reflected by the impurity at \text{Δ}W=\pm 20\%\overline{W}, \text{Δ}J=\pm 5\%\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=\pm 7\%({\overline{\chi }}_1+{\overline{\chi }}_2), 0.67\overline{M}⩽M⩽2\overline{M}, ε =1 meV, and \text{Δ}{\epsilon }_0=\epsilon \left|{\beta }_n\right|, \left|{\beta }_n\right|⩽0.5, and it is disrupted at \text{Δ}W=\pm 20\%\overline{W}, \text{Δ}J=\pm 10\%\overline{J}, \begin{array}{ll}\text{Δ}({\chi }_1+{\chi }_2)\hfill & =\pm 9\%({\overline{\chi }}_1+{\overline{\chi }}_2)\hfill \end{array}, and \text{Δ}{\epsilon }_0=\epsilon \left|{\beta }_n\right|\mathrm{,}\left|{\beta }_n\right|⩽\mathrm{0.5,}\epsilon =1\text{meV}.

In summary, we have looked at the stability of new solitons under the influence of various structural disorders in the proteins via the fourth-order Runge-Kutta method which can identify and check the stability of state of a new soliton transporting bio-energy. The above results definitely show that the new solitons in the improved model are very robust against these structural disorders. However, the Davydov soliton with one quantum is not so, and smaller structural disorders in the proteins will destroy its stability [18, 19]. Forner’s results [18, 19], obtained by same method, showed that in the case of a mass disorder of \begin{array}{l}0.67\overline{M}⩽M_k⩽1.80\overline{M}\hfill \end{array}, or with fluctuations of \pm 30\%\overline{W}, or \pm 5\%\overline{J}, or \pm 20\%\overline{{\chi }_1}, or ε < 0.5 meV in Δε_0n = εδ_n, no change in the Davydov soliton dynamics were found, but for the other conditions greater than these values, the Davydov soliton was either dispersed or destroyed. Obviously, the structural disorders needed for the Davydov soliton to be stable are smaller than that of the new soliton in the improved model.

What is this, or in saying, what are the reasons why the new soliton in the improved model is comparably more stable than Davydov’s? Very clearly, this is due to the fact that the improved model with a quasi-coherent two-quanta state differs from the Davydov model with only one quantum state. Although Eqs. (4) and (5) can become the dynamic equations in the Davydov model when χ₂ =0 and \sqrt{2}{a}_n are replaced by A_n, (in such a case, the normalization condition of the Davydov wavefunction then \begin{array}{rrrr}\hfill \text{becomes}& \hfill {\displaystyle {\sum }_n{\left|A_n(t)\right|}^2}& \hfill =& \hfill 2{\displaystyle {\sum }_n{\left|a_n(t)\right|}^2=2}\end{array}. This again shows clearly that the new wavefunction in the improved model contains exactly two quanta, instead of three quanta [21, 22], when compared with the Davydov wave function containing one quantum), with the nonlinear coupling energy, G_p, and a binding energy, E_BP, determination of the features of the new soliton are greatly increased and larger by about three and twenty times than Davydov’s soliton due to simultaneous changes in the Hamiltonian and wave function. In other word, a new interaction between the acoustic vibration of the amino acid and amide-I is added in the Hamiltonian, wave function with the quasi-coherent two-quantum state are used for the system in the improved model, thus the new soliton becomes very stability.

2.3 Infuence of temperature on a new soliton in a single chain

2.3.1 The simulation method and results at 300 K in case of no structural disorders

We now study the influences of thermal perturbation at biological temperatures on the behavior of the new soliton.

Since biological proteins are always in an environment (or heat bath) with a biological temperature of 300 K, we should firstly decide whether or not the thermal motion of the lattice still permits soliton motion in such a case. In accordance with thermodynamic theory, from a C=O vibrational energy of ε₀ = 3.28 × 10⁻²⁰N·m and a temperature of 300 K, one finds that the Boltzmann factor is 3×10⁻⁴ in this case. This means that only three out of 10 000 excitons are excited in the thermal equilibrium at 300 K. Then one can safely assume that the temperature of environment affects primarily soliton motion via only the lattice of the amino acid. Thus, prior to soliton commencement the system is in equilibrium with the temperature of environment and the lattice is in thermal motion (which can be described as a linear combination of its normal modes), while the exciton system is in its ground state. With soliton commencement a non-equilibrium state is created, the state of the lattice changes, and the behavior of the soliton is influenced by the thermal perturbation via the lattice. This amounts to assuming that the time the soliton needs to travel through the protein is small compared to the time the temperature of environment needs to re-establish equilibrium with the system, since the soliton velocity is high. However, how do we represent the effect of the heat bath on the lattice? Here we adopt Lomdahl and Kerr’s method [29−31]. Lomdahl, and Kerr, and Lawrence and co-workers [29−31] found that the Davydov soliton is destroyed at 300 K by using its corresponding dynamic equations, thus the validity and availability of the Davydov theory was questioned in their results. In Lomdahl and Kerr’s numerical calculation methods [29, 30], the decay term, M\Gamma {\dot{q}}_n, arising from the solution medium of proteins and random noise term, F_n(t), resulting from the temperature of environment, T, were added into the displacement equation of the amino acid molecules. Utilizing this idea, Eq. (7) becomes

\begin{array}{}\text{(13)}& \begin{array}{ll}M{\ddot{q}}_n(t)=\hfill & W[q_{n+1}(t)-2q_n(t)+q_{n-1}(t)]\hfill \\ \hfill & +2{\chi }_1[|a_{n+1}{|}^2-|a_{n-1}{|}^2]\hfill \\ \hfill & +2{\chi }_2\{a_n^{*}(t)[a_{n+1}(t)-a_{n-1}(t)]\hfill \\ \hfill & +a_n(t)[a_{n+1}^{*}(t)-a_{n-1}^{*}(t)]\}\hfill \\ \hfill & -M\Gamma {\dot{q}}_n+F_n(t)\hfill \end{array}\end{array}

where Γ is a dissipation coefficient for the vibration of amino acids, which is about 10⁸s⁻¹ for the proteins. Thus M\Gamma {\dot{q}}_n and F_n(t) should also be added into Eq. (11). Then we can find the solutions of Eqs. (8)−(13), a(t)r_n and a(t)i_n, and the corresponding |a_n(t)|² by using the above method. In such a case we must give an explicit representation of F_n(t).

In accordance with statistical physics, the thermal noise term F_n(t) is related with the temperature of the systems, and the average value of its correlation function can be represented by

\begin{array}{}\text{(13)}& \langle F(x\mathrm{,}t)F(\mathrm{0,0})\rangle =2M\Gamma K_{\text{B}}T\delta (x)\delta (t)/{\tau }^{\prime }\end{array}

where τ′ is the damping constant. We assume that the deviation of the random noise satisfies the normal distribution with criterion deviation and has zero expectation value, thus it can be expressed by N(F_n)=\frac{1}{\sqrt{2\text{π}\sigma }}\text{exp}[-F_n^2/(2\sigma )], where σ = 2M K_BT Γ /τ′, τ′ is the time constant, and Γ is measured by an inverse number of the time constant of environment with biological temperature. Thus, F_n can be denoted by F_n(t)=\sqrt{\sigma }{\displaystyle \sum _{r=1}^L\left[X_{nr}(t)-\frac{1}{2}\right]}; here the random number X_nr (t) is in the region of (0⩽X_{nr}⩽1). We here assume that L = 12, then the deviation of [X_nr (t)−1/2] is 1/12, the domain of the random noise force is just |F_n(t)|⩽6\sqrt{\sigma }. Thus, F_n(t) can be represented by \sqrt{\sigma }. Therefore the temperature of the system is considered in this calculation.

From Eqs. (6) and (13) we can study the influence of the dissipative and random noise forces on the states of the new soliton by the above fourth-order Runge-Kutta method. When the soliton commences (t=0), i.e., the initial condition of a_n(0)=\text{Asech}\left[(n-n_0)({\overline{\chi }}_1+{\overline{\chi }}_2{)}^2/(4\overline{JW})\right] is added in one end of the protein molecule. Although the lattice energy fluctuations associated with the environment are larger by roughly three orders of magnitude than the local lattice energies associated with the soliton motion, we can still see that the soliton moves through the chain completely undisturbed at 300 K. This result is shown in Fig. 17, where above average values for the physical parameters in Eqs. (6) and (13) are used. We know from this figure that despite the large lattice-energy fluctuations occurring due to the environment with biological temperature, the nonlinear coupling between the lattice and exciton is still able to stabilize the soliton, or in other words, the thermal perturbation at 300 K cannot destroy the new soliton, which is consistent with the analytical results [59−71]. This stable property of the new soliton can still be retained in cases of a long time period of 300 ps at 300 K or at a higher temperature at 310 K, as shown in Fig. 18 and Fig. 19, respectively. These results show clearly that the new soliton in the improved model is thermally stable in the region between 300 K to 310 K. The lifetime of the new soliton is at least about 300 ps at the biological temperature of 300 K. What does this mean? As mentioned above, the characteristic unit of time for the model is \begin{array}{l}{\tau }_0=r_0/v_0=(\overline{M}/\overline{W}{)}^{1/2}\approx 0.98\times {10}^{-13}\text{s}\hfill \end{array}, and L/r₀ = 100, τ/τ₀ > 500 is a reasonable criterion for the soliton to be a possible mechanism for energy transport in proteins, where τ is the lifetime of the soliton. τ = 300 ps means τ/τ₀ > 1000. Thus, the new soliton in this improved model is very thermally stable at 300 K, i.e., the new soliton has a long enough lifetime enabling it to play an important role in biological processes. This result agrees with analytic data in Table 1 [60−67]. At the high temperatures of 320 K and 325 K, the new soliton disperses as shown in Figs. 20 and 21, respectively. Thus we estimate that the critical temperature of the new soliton is about 320 K.

However, the Davydov soliton with one quantum does not have the above features. For a comparison between them, we give the propagation behaviors of the Davydov soliton at T = 30 K, 40 K and 300 K, respectively, obtained by Forner [25−28] in Figs. 22−24 with the same values of the physical parameters. From Figs. 22−24 we see that the Davydov soliton is only thermally stable at T⩽40\text{K}, but is destroyed completely at 300 K. Therefore, the Davydov soliton is not thermally stable at the biological temperature of 300 K; its critical temperature is only 40 K. Thus, the Davydov soliton is considerably different from the new soliton with two quanta in the improved model.

2.3.2 The effects of structure disorders on the soliton at biological temperatures

However, the influence of nonuniformity or disorder of structure in protein molecules on the new soliton has not been considered in the above calculations. As mentioned above, the nonuniformity of structure arising from the disorder in the mass sequence of amino acid molecules is always existent in the proteins. Therefore, it is quite necessary to further study the influence of disorders in the mass sequence and fluctuations of the force constant W, the dipole-dipole interaction constant J, the coupling constants, χ₁ and χ₂, between the exciton and phonon, and the ground state energy ε₀ arising from the structure nonuniformity, on the new soliton in the region of the biological temperatures 300−310 K. In such a case we must introduce the random number generators, α_k and |β_n|, to designate the random features of the mass sequences and the ground state energy in the nonuniform proteins as mentioned above, respectively. The states of the new solitons in the case of simultaneous variations of the above five physical parameters at biological temperatures arising from structure nonuniformity are numerically simulated through Eqs. (6) and (13) and the above method. The behavior of the new soliton is shown in Fig. 25, when the disorder of the mass sequence is in the region of , or 0.67 < α_k < 2, where, α_k is a random-number generator with M_k={\alpha }_k\overline{M} equal probability with a prescribed interval, and the fluctuations of the coupling constant (χ₁ +χ₂), dipole-dipole interaction constant J, force constant W and ground state energy ε₀ are about \text{Δ}({\chi }_1+{\chi }_2)=\pm 5\%\overline{({\chi }_1+{\chi }_2)}, \text{Δ}J=\pm 5\%\overline{J}, \text{Δ}W=\pm 10\%\overline{W}, Δε₀ = ε |β_n|, ε=0.4 meV, \left|\beta \right|⩽0.5, for T =300 K, T =310 K, T =315 K and T =320 K, respectively. From these figures we can see clearly that the new soliton is still thermally stable at T < 320 K, but the new soliton begins to disperse at T = 320 K. Thus, the critical temperature of the new soliton is at least 315 K when structure nonuniformity exists. Therefore, we can conclude that the new soliton is very robust against thermal perturbation and structure nonuniformity among protein molecules at biological temperatures.

Obviously, the above thermal stability of the new soliton is due to the fact that the new interaction between the acoustic vibration of the amino acid and amide-I was added into the Hamiltonian of the system, and the wavefunction of the quasi-coherent two-quantum state for the excitons was used simultaneously in the new model. However, if we adopt only the wavefunction of the quasi-coherent two-quantum state or the excitation state of a single particle with two quanta, for the excitons in Eq. (1) and letting χ₂=0, i.e., the new interaction term in Eq. (2) is not considered, the new soliton is not thermally stable at 300 K, begins to disperse at 280 K, and is only stable at 270 K using the above method, as is shown in Fig. 26. On the contrary, if we let {\chi }_2\ne 0 in Eq. (2), but only for two terms of the wavefunction of the exciton in Eq. (1), or the excitation state of a single particle with one quantum is adopted, then the state of the new soliton is only thermally stable at 150 K, begins to disperse at 170 K, but is destroyed at 180 K, as is shown in Fig. 27. These results show clearly and again that the Hamiltonian and wavefunction in Eqs. (1) and (2) in the improved model are very correct and relevant for the protein molecules; these can result in complete thermally stability of the new soliton at biological temperatures. Very clearly, this phenomenon is due to the new improved model we used as mentioned above, in which a new interaction term is added into the Hamiltonian and the wavefunction of the quasi-coherent two-quantum state is used simultaneously. In such a case, the nonlinear coupling energy G_p and binding energy, E_BP, determining the features of the new soliton, are greatly increased relative to those in Davydov model. In the improved model G_{\text{p}}=2G_{\text{D}}\left[1+2\left(\frac{{\chi }_2}{{\chi }_1}\right)+{\left(\frac{{\chi }_2}{{\chi }_1}\right)}^2\right] and E_{\text{BP}}=16E_{\text{BD}}\left[1+4\left(\frac{{\chi }_2}{{\chi }_1}\right)+6{\left(\frac{{\chi }_2}{{\chi }_1}\right)}^2+4{\left(\frac{{\chi }_2}{{\chi }_1}\right)}^3+{\left(\frac{{\chi }_2}{{\chi }_1}\right)}^4\right], respectively, where G_{\text{D}}=4{\chi }_1^2/[(1-s^2)W] and E_{\text{BD}}=-{\chi }_1^4/(3JW) are corresponding values of the Davydov soliton. Thus, G_p and E_BP are larger than G_D and E_BD in Davydov’s soliton. by about three and twenty times, respectively. Then the stability of the new soliton is greatly enhanced in the improved model. We concretely calculated values using the above parameters for the α-helix protein molecules [1−3, 55−67] showing that the binding energy of the new soliton E_BP = 7.8 × 10⁻²¹ J is also great by about two times the thermal energy K_BT = 4.14×10⁻²¹J at 300 K. Then, the binding energy can completely suppress the destroying effect of thermal perturbation on the new soliton. Therefore, the new soliton with a quasi-coherent two-quantum state is very robust against thermal fluctuation, and is thermally stable at biological temperatures. However, the binding energy of the Davydov soliton is only E_BD = 0.188 × 10⁻²¹J by same parameter values, which is smaller than the above E_BP, and the thermal energy K_BT at 300 K is less by about forty times. Thus, it is easily destroyed by thermal perturbation. Hence, the Davydov soliton is not thermally stable at 300 K, as verified by Forner and Lomdahl et al.’s numerical simulations [24−30] and Cottingham et al.’s [33, 34] analytic results.

From the above results we see that the new solitons in the improved model are quite stable up to 320 K. However, how do the interaction parameters change if the stability is given for temperatures smaller than 100 K? On the other hand, if the critical temperature of the soliton is 330 K, then how do we choose the interaction parameters? We further study these problems. If these interaction parameters are chosen as \overline{J}=9.4\times {10}^{-4}\text{eV}, \overline{W}=19\text{N/m}, {\overline{\chi }}_1=5.6\times {10}^{-11}\text{N} and {\overline{\chi }}_2=5\times {10}^{-12}\text{N}, {\epsilon }_0={\overline{\epsilon }}_{0}0.205\text{ev}, \overline{M}=1.25\times {10}^{-25} kg, then the new soliton is stable at 90 K, but is dispersed at 95 K, as shown in Fig. 28. If these interaction parameters are chosen to be \overline{J}=9.98\times {10}^{-4}\text{eV}, \overline{W} = 13 N/m, {\overline{\chi }}_1=6.2\times {10}^{-11}\text{N} and {\overline{\chi }}_2=20\times {10}^{-12}\text{N}, {\epsilon }_0={\overline{\epsilon }}_0=0.207\text{eV}, \overline{M}=1.57\times {10}^{-25}\text{kg}, then the new soliton is stable at 330 K, but is dispersed at 335 K. The results are shown in Fig. 29. Therefore in such a case, the critical temperature of the new soliton is about 335 K (or 62℃). Thus we can conclude from these results that if the values of the coupling constants, χ₂ and χ₂, are reduced but the value of the elasticity constant is increased, then the soliton has a lower critical temperature. On the contrary, the soliton has a higher critical temperature if the values of the coupling constants, χ₂ and χ₂, are increased and the value of the elasticity constant is reduced. However, here we should point out that the protein molecules are already degenerated at high temperatures, thus their biofunctions disappear. Therefore, discussing these states of the new solitons in conditions of high temperatures of 330−400 K have only a theoretical meaning, but no practical and biological effects.

3 The states and properties of solitons in α-helix protein molecules with three channels

3.1 Numerical simulation method

For the α-helix protein molecules with three channels, the Hamiltonian and the wave function in Eqs. (1) and (2) are respectively replaced by [60−76, 81]:

\begin{array}{}\text{(14)}& \begin{array}{ll}H\hfill & =H_{\text{ex}}+H_{\text{ph}}+H_{\text{inx}}\hfill \\ \hfill & ={\displaystyle \sum _n}[{\epsilon }_0{B}_{n\alpha }^{+}{B}_{n+1\alpha }-J(B_{n\alpha }^{+}{B}_{n+1\alpha }+B_{n\alpha }{B}_{n+1\alpha }^{+})]\hfill \\ \hfill & +{\displaystyle \sum _n}\left[\frac{P_{n\alpha }^2}{2M}+\frac{1}{2}W{(q_{n\alpha }-q_{n-1\alpha })}^2\right]\hfill \\ \hfill & +{\displaystyle \sum _n}[{\chi }_1(q_{n+1\alpha }-q_{n-1\alpha })B_{n\alpha }^{+}{B}_{n\alpha }\hfill \\ \hfill & +{\chi }_2(q_{n+1\alpha }-q_{n-1\alpha })(B_{n+1\alpha }^{+}{B}_{n\alpha }+B_{n\alpha }^{+}{B}_{n+1\alpha })\hfill \\ \hfill & +L(B_{n\alpha }^{+}{B}_{n+1\alpha }+B_{n\alpha }^{+}{B}_{n-1\alpha })]\hfill \end{array}\end{array}

\begin{array}{}\text{(15)}& \begin{array}{ll}\hfill & |\Phi (t)\rangle =|\Phi (t)\rangle |\beta (t)\rangle \hfill \\ \hfill & =\frac{1}{\lambda }\left[1+{\displaystyle \sum _{n\alpha }}{a}_{n\alpha }(t)B_{n\alpha }^{+}+\frac{1}{2!}{\left({\displaystyle \sum _{n\alpha }}{a}_{n\alpha }(t)B_{n\alpha }^{+}\right)}^2\right]|0{\rangle }_{\text{ex}}\hfill \\ \hfill & \cdot \exp \left\{-\frac{\text{i}}{\hslash }{\displaystyle \sum _n}[q_{n\alpha }(t)P_{n\alpha }-{\pi }_{n\alpha }(t)u_{n\alpha }]\right\}|0{\rangle }_{\text{ph}}\hfill \end{array}\end{array}

where subscript α=1,2,3 denote the numbers of the three channels, and L is the coefficient of the chain-chain interactions among the three channels in the protein molecules. From the time-dependent Schr¨odinger equation

\begin{array}{}\text{(16)}& H|\Phi \rangle =\text{i}\hslash \frac{\partial }{\partial t}|\Phi \rangle \end{array}

with the Hamiltonian Eq. (3), the above time-dependent wave function Eq. (4), and

\begin{array}{}\text{(17)}& \text{i}\hslash \frac{\partial }{\partial t}\langle |\Phi (t)|u_n|\Phi (t)\rangle =\langle \Phi (t)|[u_n\mathrm{,}H]|\Phi (t)\rangle \end{array}

\begin{array}{}\text{(18)}& \text{i}\hslash \frac{\partial }{\partial t}\langle \Phi (t)|P_n|\Phi (t)\rangle =\langle \Phi (t)|[P_n\mathrm{,}H]|\Phi (t)\rangle \end{array}

And further considering the neighboring interactions among the three channels, we can find that

\begin{array}{}\text{(19)}& \begin{array}{ll}\text{i}\hslash {\dot{a}}_{n\alpha }(t)=\hfill & {\epsilon }_0{a}_{n\alpha }(t)-J[a_{n+1\alpha }(t)+a_{n-1\alpha }(t)]\hfill \\ \hfill & +{\chi }_1[q_{n+1\alpha }(t)-q_{n-1\alpha }(t)]a_{n\alpha }(t)\hfill \\ \hfill & +{\chi }_2[q_{n+1\alpha }(t)–q_{n-1\alpha }(t)][a_{n+1\alpha }(t)+a_{n-1\alpha }(t)]\hfill \\ \hfill & +\frac{5}{2}\{w(t)-\frac{1}{2}{\displaystyle \sum _m}[q_{m\alpha }(t){\pi }_{m\alpha }(t)\hfill \\ \hfill & -{\dot{\pi }}_{m\alpha }(t){\dot{q}}_{m\alpha }(t)]\}{a}_{n\alpha }(t)\hfill \\ \hfill & +L[a_{n\alpha +1}(t)+a_{n\alpha -1}(t)]\hfill \end{array}\end{array}

\begin{array}{}\text{(20)}& \begin{array}{ll}M{\ddot{q}}_{n\alpha }=\hfill & W[q_{n+1\alpha }(t)-2q_{n\alpha }(t)+q_{n-1\alpha }(t)]\hfill \\ \hfill & +2{\chi }_1[|a_{n+1\alpha }(t{)|}^2-|a_{n-1\alpha }(t{)|}^2]\hfill \\ \hfill & +2{\chi }_2\{a_{n\alpha }^{*}(t)[a_{n+1\alpha }(t)-a_{n-1\alpha }(t)]\hfill \\ \hfill & +a_{n\alpha }(t)[a_{n+1\alpha }^{*}(t)-a_{n-1\alpha }^{*}(t)]\}\hfill \end{array}\end{array}

We can eliminate the term containing ε₀ in Eq. (8) by the following transformation:

\begin{array}{}\text{(21)}& {\phi }_{n\alpha }(t)=a_{n\alpha }(t)\text{exp}(-\text{i}{\epsilon }_{0}t/\hslash )\end{array}

Because a_n(t) in Eqs. (8) and (9) is a complex function, we can thus make the transformation

\begin{array}{}\text{(21)}& a_{n\alpha }(t)=a{r}_{n\alpha }(t)+\text{i}a{i}_{n\alpha }(t)\end{array}

with

\begin{array}{}\text{(22)}& {|a_n|}^2={|a{r}_n|}^2+{|a{i}_n|}^2\end{array}

Thus, Eqs. (8) and (9) change into

\begin{array}{}\text{(23)}& \begin{array}{ll}\hslash \dot{a}{r}_{n\alpha }=\hfill & -J(a{i}_{n+1\alpha }+a{i}_{n-1\alpha })+{\chi }_1(q_{n+1\alpha }-q_{n-1\alpha })a{i}_{n\alpha }\hfill \\ \hfill & +{\chi }_2(q_{n+1}-q_{n-1\alpha })(a{i}_{n+1\alpha }+a{i}_{n-1\alpha })\hfill \\ \hfill & +L[a{i}_{n\alpha +1}(t)+a{i}_{n\alpha -1}(t)]\hfill \end{array}\end{array}

\begin{array}{}\text{(24)}& \begin{array}{ll}-\hslash \dot{a}{i}_{n\alpha }=\hfill & -J(a{r}_{n+1\alpha }+a{r}_{n-1\alpha })\hfill \\ \hfill & +{\chi }_1(q_{n+1\alpha }-q_{n-1\alpha })a{r}_{n\alpha }\hfill \\ \hfill & +{\chi }_2(q_{n+1\alpha }-q_{n-1\alpha })(a{r}_{n+1\alpha }+a{r}_{n-1\alpha })\hfill \\ \hfill & +L[a{i}_{n\alpha +1}(t)+a{i}_{n\alpha -1}(t)]\hfill \end{array}\end{array}

\begin{array}{}\text{(25)}& {\dot{q}}_{n\alpha }=\frac{y_{n\alpha }}{M}\end{array}

\begin{array}{}\text{(26)}& \begin{array}{ll}{\dot{y}}_{n\alpha }=\hfill & W[q_{n+1\alpha }-2q_{n\alpha }+q_{n-1\alpha }]\hfill \\ \hfill & +2{\chi }_1[a{r}_{n+1\alpha }^2+a{i}_{n+1\alpha }^2-a{i}_{n-1\alpha }^2-a{i}_{n-1\alpha }^2]\hfill \\ \hfill & +4{\chi }_2[a{r}_{n\alpha }(a{r}_{n+1\alpha }-a{r}_{n-1\alpha })\hfill \\ \hfill & +a{i}_{n\alpha }(a{i}_{n+1\alpha }-a{i}_{n-1\alpha })]\hfill \end{array}\end{array}

where ar_n and ai_n are real and imaginary parts of a_n(t) Eqs. (23)−(26) are equations we use in simulation calculations. The dynamic properties of a peptide group or amino acid molecule are described by the above four-equations [79, 80]. Thus, the protein molecules consisting of N amino acid molecules should associatively solve 4N equations. From these equations we can find their solutions, and ar_n and ai_n, by numerical simulation with the fourth-order Runge-Kutta method and by using the following initial condition at a point n₀α: aφ_nα(t) = A sech[(n−n₀)(χ₁ + χ₂)²/(4JW)] for the α-helix protein molecules with three channels, where α=1,2,3, and A is the normalization factor. As such, we can determine the solutions of Eqs. (19) and (20). However, in the simulation calculations, the above boundary conditions must be satisfied, i.e., (1) the energy of the soliton must remain constant up to 0.0012 %; (2) the probability of the soliton must be normalized at any time; (3) the energy of the soliton is real, its imaginary part must approach zero up to an accuracy of 0.001 feV. In accordance with these three criteria and utilizing the above equations and the above initial conditions we can calculate the evolution of time and space for the probability by using MATLAB language and data-parallel programming, where the time step size is chosen as 0.01 ps. In this calculation the values of the physical parameters we use are as follows: The mass \overline{M}=5.73\times {10}^{-25}\text{kg=114}\times \text{3 amu} (atomic mass units), 114 amu is the mass of myosine, \overline{W}\text{=39 N/m}, {\overline{\epsilon }}_0=0.2035\text{eV}, \overline{J}=9.68\times {10}^{-4}\text{eV}, {\overline{\chi }}_1=6.2\times {10}^{-11}\text{N}, {\overline{\chi }}_2=(10-18)\times {10}^{-12}\text{N} and \overline{L}=1.5\text{meV} for the α-helix protein molecules with three channels [5, 82]. We numerically calculate their solutions related to time and the probability of the soliton occurring at the nth amino acid molecule, |a_j|².

3.2 The properties of solitons in uniform protein molecules

When the above initial condition is imported from the end of the molecular chain, the numerical solutions of Eqs. (23)−(26), obtained by using the fourth-order Runge-Kutta method [79, 80], and the above average values of the parameters for the α-helix protein molecules with three channels in the improved model, are shown in Fig. 30. In Fig. 30(a) we show the behaviors of motion of the solution, with the initial condition of φ_nα(t) =A sech[(n−n₀)(χ₁ + χ₂)²/(4JW)], where α=1,2,3 are simultaneously motivated on the first ends of the three channels. From this figure we see that this solution can retain a clock shape while moving over a long distance in the range of the spacing of 400 amino acid residues at a time of 40 ps without dispersion along the molecular chains, i.e., this solution is a soliton. Therefore, Eqs. (19) and (20) have exactly soliton solutions with clock shapes. This is similar to the analytic results obtained for proteins with single channels in the continuum approximation in this model in which the dynamical equation is a standard nonlinear Schrodinger equation [61−70]. In Fig. 30(b) and (c) we plot the features of motion of the solutions, with the initial conditions a_nα(t = 0) = 0, where α=1,2, φ_nα(t) = A sech[(n−n₀)(χ₁ +χ₂)²/(4JW)] and φ_nα(t) = A sech[(n−n₀)(χ₁ + χ₂)²/(4JW)], where α=1,2, and a_n3(t = 0) = 0 are used, respectively. These initial conditions denote that the first ends of one channel and two channels are motivated, and that the other two channels and the single channel are not linked, respectively. We see from Fig. 30(b) that two new waves with small amplitudes are generated, except for one soliton occurring on the channel linked by the above initial conditions. Obviously, the two new waves are still excited from the above initial conditions through the interactions among the three channels. Although the two excitations are small, they can move over long distances along the two chains while keeping their amplitudes. Therefore, they are still solitons with small amplitudes and clock shapes. However, in Fig. 30(c) we see a strange phenomenon in that the amplitudes of solitons generated in the two motivated chains are small, and in an unmotivated chain the soliton formed from the superposition of waves induced by the other two chains is greater. From this study we already know that the solitons formed have similar and higher energy when the initial conditions are simultaneously motivated on the first ends of the three channels, corresponding to practical cases, but the soliton feature of the solutions is substandard in the other two cases and do not occur in practical cases.

In order to confirm further the soliton feature of the solutions of Eqs. (19) and (20), we further study the collision property of solitons with clock shapes, set up from opposite ends of the channels for α-helix protein molecules with three channels; the result is shown in Fig. 31(a), where the above initial conditions simultaneously motivate the opposite ends of the three channels as mentioned above. From this figure we can see clearly that the initial two solitons with clock shapes separating 100 amino acid spacings in each channel collide with each other at about 17 ps. After this collision, the two solitons in each channel go through each other without scattering and retain their clock shapes to propagate toward and separately along the three chains, satisfying the rules on collision of macroscopic particles. Thus, we can also judge from this result that the solutions of Eqs. (19) and (20) in the α-helix protein are exact solitons. In Fig. 31(b) and (c) we plot the collision features of the two solitons generated from opposite ends of the channels with the above initial conditions motivating the opposite ends of one channel and two channels, and where they are not linked with the opposite ends of the other two channels and the single channel, respectively. In these two cases the soliton feature of the solutions after collision is substandard, especially where the initial condition motivates only the opposite ends of the single channel, as shown in Fig. 30(b). Thus, the solutions of Eqs. (19) and (20) have better soliton feature, where the opposite ends of the three channels are motivated simultaneously by the above initial conditions in the α-helix protein molecules. Hence, we conclude from Figs. 30 and 31 that the soliton excited in the α-helix protein molecules has higher stability in the case of simultaneous motivation of the initial condition on the three channels, and then work only needs to be done in case of the following.

In the above simulation we study only the behavior of the soliton under the condition of a short time period of 40 ps, which clearly exhibits the transport feature of this soliton. However, what is its behavior in cases of longer time periods and larger spacings? We thus study further the behaviors of long time periods for the solutions of Eqs. (23)−(26) in α-helix proteins with three channels. In Fig. 32 we show the results of soliton solutions of Eqs. (19) and (20) obtained at a time of 120 ps and at 300 amino acid spacings, where the above initial conditions are simultaneously linked on the first ends of the three channels. We can see clearly from Fig. 32 that in such cases, the soliton still retains its amplitude and clock shape while moving. This result shows that the lifetime of the soliton is at least 120 ps. What is the meaning of a lifetime of 120 ps? In the light of above calculation method we can find out that the lifetime of the soliton, τ = 120 ps, corresponds to τ /τ₀ > 700, which thus is larger that the reasonable criterion for the soliton of τ/τ₀ > 500. This means that the soliton in the improved model is a possible carrier of bio-energy transport in proteins. This conclusion also agrees with analytic results in Table 1 [60−76, 81]. This shows that our analytic results and the improved model are both correct [60−76, 81].

3.3 The properties of solitons in disordered α-helix protein molecules

In the investigation of the influence of quantum fluctuations of physical parameters arising from structure disorder, we should first introduce the random number generators, α_k and |β_n|, to designate the disorder features of the mass sequences or to create random sequences of the masses of amino acids and ground state energy as mentioned above. we should use M_k={\alpha }_k\overline{M} to designate the mass of amino acid residue at each point in a random series of masses within a prescribed interval, and use Δε₀ = ε−ε₀ = ε |β_n| to represent the random features of the ground state energy caused by different amino acid side groups and corresponding local geometric distortions due to the impurities imported and changes of side radicals in the disordered protein. For an isolated impurity at the nth site in the molecular chains, the change in ground state energy is denoted by Δε_0n = εδ_n. The quantum fluctuations of the spring constant, the exciton-phonon coupling constant, the chain-chain interaction coefficient and the dipole-dipole interaction arising from the structure disorder are still represented by \text{Δ}W=W-\overline{W}, \begin{array}{l}\text{Δ}({\chi }_1+{\chi }_2)=[({\chi }_1+{\chi }_2)-({\overline{\chi }}_1+{\overline{\chi }}_2)]\hfill \end{array}, \text{Δ}L=L-\overline{L} and \text{Δ}J=J-\overline{J}, respectively. In the following, we will use these representations to study the influence of these fluctuations and structure disorders on the solitons excited in the α-helix proteins using the fourth order Runge-Kutta method [79, 80].

3.3.1 The individual influences of different fluctuations and disorders on the motion of solitons

In accordance with the above method we study the in-fluences of these nonuniformities, whether of the smaller or larger intervals for α_k, such as, 0.67⩽{\alpha }_k⩽100 or 0.67⩽{\alpha }_k⩽300, on the soliton. This result shows that this nonuniformity does not significantly affect the stability of the new soliton and the vibrational energy is not dispersed, as shown in Fig. 33., where other parameters take mean values. The interval over which the new soliton moves unperturbed is evidently larger than the variation in the masses of the natural amino acids (0.67⩽{\alpha }_k⩽1.80). Therefore, the soliton in Pang’s model is very robust against mass nonuniformities in α-helix protein molecules.

Numerical simulation results show that we do not find any changes in the dynamics of the soliton up to a random fluctuation of \pm 25\%\overline{W}. The results of our calculations are shown in Fig. 34, in the cases where fluctuations are \pm 15\%\overline{W} and \pm 25\%\overline{W}, respectively. For \pm 30\%\overline{W}, the soliton velocity is only somewhat diminished when compared with the case of \overline{W}. Finally, for ± 40% \overline{W}, the soliton disperses completely and propagation is irregular.

The soliton is very sensitive to variations in the dipole-dipole interaction J in the α-helix protein molecules. The soliton is stable up to 5% as shown in Fig. 35, where we denote the features of the solitons under changes of \pm 3\%\overline{J} and \text{Δ}J=\pm 5\%\overline{J}, but the soliton disperses at \text{Δ}J=\pm 7\%\overline{J}.

In Fig. 36 we exhibit the changes of phase of the soliton under the influence of fluctuations of \text{Δ}({\chi }_1+{\chi }_2)=\pm 6\%({\overline{\chi }}_1+{\overline{\chi }}_2) and \pm 9\%({\overline{\chi }}_1+{\overline{\chi }}_2). These figures show that the new solitons are stable when (χ₁ + χ₂) varies up to \pm 9\%({\overline{\chi }}_1+{\overline{\chi }}_2). Meanwhile, the solitons disperse and split gradually into some small waves when (χ₁ + χ₂) are larger than this value.

With \left|{\beta }_n\right|⩽0.5 for ε=0.57 meV and 1.26 meV, we show the behaviors of the soliton in Fig. 37. These figures show that in the case of a random sequence, only if ε <1.3 meV and \left|{\beta }_n\right|⩽0.5 can the soliton pass through the chain and is stable. For higher values of ε, the soliton is reflected or dispersed.

In Fig. 38 we show the variations in features of the solitons due to the fluctuations of \text{Δ}L=\pm 5\%\overline{L} and \pm 6\%\overline{L}. These figures show that the soliton is stable only if . From the simulation we know that the smaller the fluctuation of the chain-chain interaction, the higher the stability of the soliton. When ΔL = 0 and L = 0, the stability of the soliton is the highest.

3.3.2 Associated effects of the fluctuations of six structural parameters on the solitons

We now study the associated influences of the change of mass of amino acids and the fluctuations in dipole-dipole interaction, the exciton-phonon coupling constant, the spring constant, the diagonal disorder, and the chain-chain interaction arising from the structure disorder of the α-helix protein molecules, on the behaviors of the soliton. In Fig. 39 we show the changes of property of the solitons under the influence of different fluctuations of structure parameters, in which Fig. 39(a) denotes the states of the soliton with fluctuations of 0.67 <α_k < 2, \text{Δ}W=\pm 5\%\overline{W}, \text{Δ}J=\pm 1\%\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=\pm 1\%({\overline{\chi }}_1+{\overline{\chi }}_2) and \text{Δ}L=\pm 1\%\overline{L}; Fig. 39(b) results when 0.67 < α_k < 2, \text{Δ}W=\pm 10\%\overline{W}, \text{Δ}J=\pm 2\%\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=\pm 3\%({\overline{\chi }}_1+{\overline{\chi }}_2) and \text{Δ}L=\pm 2\%\overline{L}; Fig. 39(c) results when 0.67 < α_k < 2, \text{Δ}W=\pm 10\%\overline{W}, \text{Δ}J=\pm 2\%\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=\pm 5\%({\overline{\chi }}_1+{\overline{\chi }}_2) and \text{Δ}L=\pm 2\%\overline{L};Fig. 39(d) results when 0.67 < α_k < 2, \text{Δ}W=\pm 10\%\overline{W}, \text{Δ}J=\pm 2\%\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=\pm 6\%({\overline{\chi }}_1+{\overline{\chi }}_2) and \text{Δ}L=\pm 2\%\overline{L}. From these figures we see clearly that the solitons shown in Fig. 39(a), (b) and (c) are very stable under the actions of different structure disorders, but the soliton shown in Fig. 39(d) is already unstable due to the influence of the structure disorder its amplitude changes, thus it is already dispersed. This manifests that the solitons excited in α-helix protein molecules are not stable in cases of large structure disorders. In comparing these results in Fig. 39 with those shown in Fig. 30(a) or Figs. 33−38, we find that the collective effects of these structure disorders and quantum fluctuations change the states and features of solitons, making their amplitudes, energies and velocities decrease, although such influences cannot destroy the solitons, which can still transport steadily along the molecular chains while retaining energy and momentum when the quantum fluctuations are small, i.e., the solitons are quite robust against these disorder effects. However, the solitons may be dispersed or disrupted in cases of very large structure disorders.

3.3.3 The collective effects of various nonuniformities on the motion of the soliton

We finally study the collective effects of the change in mass of amino acids and the fluctuations of dipole-dipole interaction, the exciton-phonon coupling constant, the spring constant, the diagonal disorder and the chain-chain interactions arising from the structure disorder of α-helix protein molecules, on the properties of the soliton. The results obtained are denoted in Fig. 40, in which Fig. 40(a) represents the features of the soliton formed under fluctuation conditions of 0.67 < α_k < 2, \text{Δ}W=\pm 3\%\overline{W}, \text{Δ}J=\pm 1\%\overline{J}, \begin{array}{l}\text{Δ}({\chi }_1+{\chi }_2)=\pm 1\%({\overline{\chi }}_1+{\overline{\chi }}_2)\hfill \end{array}, \text{Δ}L=\pm 1\%\overline{L} and Δε₀ = ε|β_n|, ε=0.35 meV, |β_n| < 0.5;Fig. 40(b) results when 0.67 < α_k < 2, \text{Δ}W=\pm 6\%\overline{W}, \text{Δ}J=\pm 1\%\overline{J}, \begin{array}{l}\text{Δ}({\chi }_1+{\chi }_2)=\pm 2\%({\overline{\chi }}_1+{\overline{\chi }}_2)\hfill \end{array}, \text{Δ}L=\pm 1\%\overline{L} and Δε₀ = ε|β_n|, ε=0.1 meV, |β_n| < 0.5; Fig. 40(c) results when 0.67 < α_k < 2, \text{Δ}W=\pm 8\%\overline{W}, \text{Δ}J=\pm 1\%\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=\pm 3\%({\overline{\chi }}_1+{\overline{\chi }}_2) and \text{Δ}L=\pm 1\%\overline{L} and Δε₀ = ε|β_n|, ε=0.1 meV, |β_n| < 0.5. We see from this figure that in these conditions the nature of the soliton during bio-energy transport in the α-helix protein molecules can still be maintained, but the soliton begins to disperse when these structure disorders are larger than the values 0.67 < α_k < 2, \text{Δ}W=\pm 8\%\overline{W}, \text{Δ}J=\pm 1\%\overline{J}, \begin{array}{ll}\text{Δ}({\chi }_1+{\chi }_2)\hfill & =\pm 3\%({\overline{\chi }}_1+{\overline{\chi }}_2)\hfill \end{array} and \text{Δ}L=\pm 1\%\overline{L} and Δε₀ = ε|β_n|, ε=0.1 meV, |β_n| < 0.5. Therefore, we can conclude that the new soliton in Pang’s model is very robust against various structure disorders of the α-helix protein molecules.

3.4 The states of the solitons at biological temperatures

However, as mentioned above, the α-helix protein molecules in living systems always work at the biological temperature of 300 K, therefore, we must study the behaviour of the soliton during transport at 300 K, and we should add the effect of temperature on the soliton to the above equations. In the light of Lomdahl and Kerr’s method [29, 30] we should add the decay term M\Gamma {\dot{q}}_n and random noise term, F_n(t), resulting from the temperature into the displacement equation of the amino acid molecules, Eq. (20). Thus, the latter can now be represented by

\begin{array}{}\text{(27)}& \begin{array}{ll}M{\ddot{q}}_n(t)=\hfill & W[q_{nt1}(t)-2q_n(t)+q_{n-1}(t)]\hfill \\ \hfill & +2{\chi }_1[|a_{n+1}{|}^2-|a_{n-1}{|}^2]\hfill \\ \hfill & +2{\chi }_2\{a_n^{*}(t)[a_{n+1}(t)-a_{n-1}(t)]\hfill \\ \hfill & +a_n(t)[a_{n+1}^{*}(t)-a_{n-1}^{*}(t)]\}\hfill \\ \hfill & -M\Gamma {\dot{q}}_n+F_n(t)\hfill \end{array}\end{array}

where Γ is the dissipation coefficient of vibration of the amino acids. In light of the above method we can give F_n(t) as \sqrt{\sigma }, where σ = 2M K_BT Γ/τ. Thus, we can now find the soliton solution of the equations of motion in Eqs. (19) and (27) with decay effect and random noise force using the above method and the fourth-order Runge-Kutta method [79, 80]. The results at 300 K are shown in Fig. 41 for the α-helix protein molecules with three channels, with the above initial conditions simultaneously linked on the first ends of the three channels. From this figure we can see that the new soliton in the improved model can move along the three channels at a constant speed and amplitude without dispersion. So the soliton is still thermally stable at the biological temperature of 300 K. In Fig. 42 we also show the results of soliton motion at a long time period of 120 ps and large spacings of 1000 sites (i.e., 333 amino acid residues contained in each channel) at 300 K for the α-helix protein molecules, when the above initial conditions are simultaneously linked on the first ends of the three channels. We see from this figure that the solitons are undisturbed in such conditions, and really move over a long time period and through large spacings along the protein molecular chains while retaining their amplitudes and velocities at bio-temperatures. In Fig. 43 we plot the collision behaviors of the solitons with clock shapes, set up from opposite ends of the channels in the α-helix protein molecules, when the above initial conditions are simultaneously linked on the opposite ends of the three channels. From this figure we see clearly that the initial two solitons with clock shapes, separating 100 amino acid spacings in each channel, collide with each other at about 16 ps. After the collision, the two solitons in each channel go through each other still retaining their clock shapes and propagating toward and separately along the three chains. These results show clearly that although there are large lattice fluctuations in the protein molecules due to the influence of temperature, the nonlinear coupling interaction between the amino acids and excitons is still able to stabilize the soliton, therefore this soliton is very robust against thermal perturbation of the environment. In this case the lifetime of the new soliton is also, at least, 120 ps. This means that the new soliton could play an important role in biological processes.

We now study the changes of states of the solitons with increasing temperature. In Fig. 44 we exhibit the transport properties of the soliton in the α-helix proteins at different temperatures of 295 K, 305 K, 310 K, 315 K, 320 K and 325 K, respectively. This figure shows that the amplitudes of the soliton decrease with increasing temperature, and it begins to disperse at 320 K. This means that the soliton must expend a part of its energy to retain and suppress the increasing destructive effect of thermal perturbation arising from the rise in temperature of the system, thus its amplitude or energy does decrease. We can estimate that the critical temperature of the soliton is about 320 K in this case. At the same time, we also find that the transport velocity of the soliton also decreases with the increase in temperature of the system. In Table 2 we give concrete data of these velocities at different temperatures. In the meantime, the temperature-dependence of the soliton velocity is shown in Fig. 45. Here, the decrease in soliton velocity with increasing temperature is also obvious. Evidently, this is due to the enhancement of disordered thermal motion of the medium resulting from the rise in temperature, which increases the resistance of motion of the soliton. Thus, the velocity of the soliton necessarily decreases.

The above results of our investigation manifest clearly that the nonlinear coupling interaction between the amino acids and excitons is still able to stabilize the soliton, although it undergoes destructive influences due to increases of chain-chain interaction and thermal perturbation. Therefore, the soliton is very robust against the influence of chain-chain interaction and thermal perturbation of the environment.

However, the structure nonuniformity of the protein molecules is not considered in the above calculation. Thus, we should further study the influence of structure nonuniformity of proteins on the new soliton at the biological temperature of 300 K by using the fourth-order Runge-Kutta method [18] and Eqs. (19) and (27). In Fig. 46, we plot the states of the new soliton in the nonuniform α-helix proteins at T = 300 K, when the disorder of the mass sequence is in the region of , or 0.67 < α_k < 2, and the fluctuations of (χ₁ + χ₂), J, W and ε₀ are about \text{Δ}({\chi }_1+{\chi }_2)=\pm 4\%\overline{({\chi }_1+{\chi }_2)}, \text{Δ}J=\pm 2\%\overline{J}, \text{Δ}W=\pm 4\%\overline{W},Δε₀ = ε |β_n|, ε=0.5 meV, \left|{\beta }_n\right|⩽0.5. From these figures we see clearly that the new soliton is undisturbed and still thermally stable at 300 K when these structure nonuniformity occurs in the proteins. Therefore, we can conclude that the new soliton is robust against thermal perturbation and structure nonuniformity of the protein molecules at biological temperatures. Thus, the new soliton in the improved model can be a carrier for the bio-energy transport in protein molecules.

4 Conclusions

In this paper the influences of the changes of the itself features of proteins, i.e., molecular structure disorders and environment conditions, i.e., physiological temperature and damping effect of medium on the states and properties of the soliton in transporting bio-energy have been numerically studied using the fourth-order Runge-Kutta method and an improved theory, where the structure disorders include the fluctuations of the characteristic parameters of the spring constant, the dipole-dipole interaction constant and the exciton-phonon coupling constant, as well as the chain-chain interaction coefficient among the three channels and ground state energy resulting from the disorder distribution of masses of amino acid residues and impurities. Therefore, a different variety of influenced factors on the proteins are considered and any approximate methods and calculations also have been not used in these studies. Thus our investigations are complete, concrete and realistic, then the results obtained are available and credible. In this paper, we investigated the behaviors and states of the solitons in a single protein molecular chain and in α-Helix protein molecules with three channels. In the former, we proved first that the new solitons can move without dispersion, retaining their shape, velocity and energy in uniform and periodic protein molecules. In the case of structure disorder, the fluctuations of the spring constant, the dipole-dipole interaction constant, the exciton-phonon coupling constant, the ground state energy and the disorder distributions of masses of amino acid residues of the proteins influence the states and properties of motion of the solitons, but the new solitons are quite stable and very robust against these structure disorders. Even with larger disorders in the sequence of masses and stranger fluctuations of spring constants and coupling constants, the solitons are still not dispersed or destroyed. However, when the disorder distribution of the masses and the fluctuations are quite great, the solitons disperse. If the effect of thermal perturbation of the environment on the soliton in nonuniform proteins is considered further, then the soliton is still thermally stable at the biological temperature of 300 K even with a longer time period of 300 ps and a larger spacing of 400 amino acids, which is consistent with the analytic results obtained by quantum perturbed theory. Meanwhile, the new soliton is also thermally stable in cases of motion over a long time period of 300 ps in the region of 300−320 K under the influence of structure disorders. However, the soliton disperses in cases of a higher temperature of 325 K and in larger structure disorders. Thus, we have determined that the soliton’s lifetime and critical temperature are at least 300 ps and 320 K, respectively. These results are also consistent with analytical data. Furthermore, the thermal stability of the new soliton can still be retained at the biological temperature of T = 315 K in nonuniform proteins, where the disorder of the mass sequence is , or 0.67 < α_k < 2, and the fluctuations of the J, χ₁ + χ₂, W and ε₀ are \text{Δ}J=\pm 5\%\overline{J}, \text{Δ}({\chi }_1+{\chi }_2)=\pm 5\%\overline{({\chi }_1+{\chi }_2)}, \text{Δ}W=\pm 10\%\overline{W} and Δε₀ = ε |β_n|, ε=0.4 meV, |β_n| < 1, respectively. Thus, we see clearly from these results that the new soliton in the improved model is very robust against thermal perturbation and structure nonuniformity of protein molecules.

In the α-Helix protein molecules with three channels, the results obtained show that these structure disorders and the fluctuations of the spring constant, the dipole-dipole interaction constant and the exciton-phonon coupling constant, as well as the chain-chain interaction coefficient among the three channels and ground state energy resulting from the disorder distributions of masses of amino acid residues, side radicals and impurities, can change the states and features of solitons, and make their amplitudes, energy and velocities decrease, but they cannot destroy the solitons, i.e., the solitons can still transport steadily along the molecular chains, retaining energy and momentum when the quantum fluctuations are small, as for example, with the structure disorders and quantum fluctuations of 0.67 < α_k < 2, \text{Δ}W=\pm 8\%\overline{W}, \text{Δ}J=\pm 1\%\overline{J}, \begin{array}{l}\text{Δ}({\chi }_1+{\chi }_2)=\pm 3\%({\overline{\chi }}_1+{\overline{\chi }}_2)\hfill \end{array} and \text{Δ}L=\pm 1\%\overline{L} and Δε₀ = ε|β_n|, ε=0.1 meV, |β_n| < 0.5. Therefore, the solitons in the improved model are quite robust against these disorder effects. However, the solitons may be dispersed or disrupted in cases of very large structure disorders. When the influence of temperature on the soliton is again considered, we find that the new soliton can transport steadily over 333 amino acid residues in each channel in cases of motion over a long time period of 120 ps, and still retain their shapes and energies to travel forward along protein molecules after mutual collision of the solitons at the biological temperature of 300 K. Therefore, the soliton is very robust against thermal perturbation of the α-helix protein molecules at 300 K. However, the soliton disperses in the case of a higher temperature of 325 K and in larger structure disorders. Thus, the critical temperature of the soliton is about 320 K. When the effects of structure disorder and temperature are considered simultaneously, the soliton also has higher thermal stability and can transport for a long time along the protein molecular chains, retaining its amplitude, energy and velocity even though the fluctuations of the structural parameters and the temperature of the medium increase continually; however, it disperses at larger fluctuations of , \text{Δ}({\chi }_1+{\chi }_2)=\pm 2\%({\overline{\chi }}_1+{\overline{\chi }}_2), \text{Δ}J=\pm 1.3\%\overline{J}, \text{Δ}W=\pm 6\%\overline{W}, \text{Δ}L=\pm 1.5\%\overline{L} and Δε₀ = ε |β_n|, ε=0.82 meV, \left|{\beta }_n\right|⩽0.5 at T =300 K, and at temperatures higher than 315 K when the fluctuations are , \text{Δ}({\chi }_1+{\chi }_2)=\pm 1\%({\overline{\chi }}_1+{\overline{\chi }}_2), \text{Δ}J=\pm 0.7\%\overline{J}, \text{Δ}W=\pm 7\%\overline{W}, \text{Δ}L=\pm 0.8\%\overline{L} and Δε₀ = ε |β_n|, ε=0.4 meV, \left|{\beta }_n\right|⩽0.5. Thus we find that the critical temperature of the soliton is only 315 K under these conditions.

From these investigations we can gain the following conclusions:

(1) The chain-chain interaction among these channels can decrease the stability of the soliton excited in the α-helix protein molecules due to its dispersive effect, but the stability of the soliton formed in the case of simultaneous motivation of three channels by initial conditions in the α-helix proteins is better than that in other initial conditions, and is lower than that in a protein of one channel.

(2) The soliton in the new model can steadily transport over 333 amino acid residues in cases of motion for a long time period of 120 ps, and can retain its shape and energy to travel forward after mutual collision with other solitons at the biological temperature of 300 K in α-helix protein molecules. Therefore, the soliton is very robust against thermal perturbation of the system at 300 K.

(3) The features and states of the soliton excited in α-helix proteins change with increasing temperature of the system. In general, the amplitude and velocity of the soliton decrease with increases in temperature of the medium. These effects are due to the enhancement of thermal motion disorder of the medium arising from the rise in temperature, which increases the resistance of motion of the soliton.

(4) We find that a critical temperature of the soliton occurs in α-helix protein molecules using the numerical simulation method, and it is about 315 K under influences of structure nonuniformity and thermal perturbation. These are a few new results which have not been obtained in α-helix protein molecules until now.

From these investigations we can conclude that the new solitons in Pang’s model are very robust against the influence of structure disorders and thermal perturbation of proteins at the biological temperature of 300 K in α-helix protein molecules, and are possible carriers of bio-energy transport; the improved model is a possible candidate for the mechanism of this transport, which can thus play an important role in studies of the bio-energy transport in α-helix protein molecules at biological temperatures.

References

[1] X. F. Pang, Front. Phys. China, 2007, 2(4): 469

[2] A. S. Davydov, J. Theor. Biol., 1973, 38: 559

[3] A. S. Davydov, Phys. Scr., 1979, 20: 387

[4] A. S. Davydov, Solitons in Molecular Systems, Dordrecht: Reidel Publishing Comp., 1985; 2nd Ed., 1991

[5] S. Gheorghiu-Svirschevski, Phys. Rev. E, 2001, 64: 051907

[6] Z. Ivic, S. Zekiric, and Z. Przuly, Phys. Lett. A, 2002, 306:144

[7] D. Henning, Phys. Rev. E, 2000, 61: 4550

[8] D. Henning, Euro. Phys. J. B, 2001, 24: 377

[9] P. L. Christiansen and A. C. Scott, Davydov’s soliton revisited, New York: Plenum Press, 1990

[10] A. C. Scott, Phys. Rev. A, 1982, 26: 578

[11] A. C. Scott, Phys. Rep., 1992, 217: 1

[12] A. C. Scott, Physica D, 1990, 51: 333

[13] D. W. Brown, Phys. Rev. A, 1988, 37: 5010

[14] D. W. Brown, and Z. Ivic, Phys. Rev. B, 1989, 40: 9876

[15] Z. Ivic and D. W. Brown, Phys. Rev. Lett., 1989, 63: 426

[16] D. W. Brown, B. J. West, and K. Lindenberg, Phys. Rev. A, 1986, 33: 4104

[17] D. W. Brown, K. Lindenberg, and B. J. West, Phys. Rev. B, 1987, 35: 6169

[18] X. F. Pang, Chin. J. Biochem. Biophys., 1986, 18: 1

[19] X. F. Pang, Chin. J. Atom. Mol. Phys., 1986, 6: 1986

[20] X. F. Pang, Chin. Appl. Math., 1986, 10: 276

[21] L. Cruzeiro, J. Halding, P. L. Christiansen, O. Skovgard, and A. C. Scott, Phys. Rev. A, 1985, 37: 703

[22] L. Cruzeiro-Hansson, Phys. Rev. A, 1992, 45: 4111

[23] L. Cruzeiro-Hansson, Phys. Rev. Lett., 1994, 73: 2927

[24] L. Cruzeiro-Hansson, V. M. Kenker, and A. C. Scott, Phys. Lett. A, 1994, 190: 59

[25] W. Förner, Phys. Rev. A, 1991, 44: 2694

[26] W. Förner, J. Phys.: Condens. Matt., 1991, 3: 1915, 3235

[27] W. Förner, J. Phys.: Condens. Matt., 1992, 4: 4333

[28] W. Förner, J. Phys.: Condens. Matt., 1993, 5: 823

[29] P. S. Lomdahl and W. C. Kerr, Phys. Rev. Lett., 1985, 55:1235

[30] W. C. Kerr, and P. S. Lomdahl, Phys. Rev. B, 1989, 35:3629

[31] X. Wang, D. W. Brown, and K. Lindenberg, Phys. Rev. A, 1998, 37: 3357

[32] X. Wang, D. W. Brown, and K. Lindenberg, Phys. Rev. Lett., 1989, 62: 1796

[33] J. P. Cottingham and J. W. Schweitzer, Phys. Rev. Lett., 1989, 62: 1792

[34] J. W. Schweitzer, Phys. Rev. A, 1992, 45: 8914

[35] H. Bolterauer and M. Opper, Z. Phys. B, 1991, 82: 95

[36] X. F. Pang, J. Phys.: Condens. Matt., 1990, 2: 9541

[37] X. F. Pang, Phys. Rev. E, 1994, 49: 4747

[38] X. F. Pang, Euro. Phys. J. B, 1999, 10: 415

[39] X. F. Pang, Physica D, 2001, 154: 138

[40] X. F. Pang, Chin. Phys. Lett., 1993, 10: 381

[41] X. F. Pang, Chin. Phys. Lett., 1993, 10: 437

[42] X. F. Pang, Chin. Phys. Lett., 1993, 10: 573

[43] X. F. Pang, Chinese Physics, 2000, 9: 108

[44] X. F. Pang, Commun. Theor. Phys., 2001, 35: 763

[45] X. F. Pang, The Theory of Nonlinear Quantum Mechanics, Chongqing: Chongqing Press, 1994

[46] X. F. Pang, Chinese Science Bulletin, 1993, 38: 1572

[47] X. F. Pang, Chinese Science Bulletin, 1993, 38: 1665

[48] X. F. Pang, Acta Phys. Slovaca., 1998, 48: 99

[49] X. F. Pang, Chin. J. Infrared Mill. Waves, 1993, 12: 377

[50] X. F. Pang, Chin. J. Infrared Mill. Waves, 1997, 16: 64

[51] X. F. Pang, Chin. J. Infrared Mill. Waves, 1997, 16: 232

[52] X. F. Pang, Acta. Phys. Sinica, 1993, 42: 1841

[53] X. F. Pang, Acta. Phys. Sinica, 1997, 46: 625

[54] X. F. Pang, Chin. J. Atom. Mol. Phys., 1987, 5: 383

[55] X. F. Pang, Chin. J. Atom. Mol. Phys., 1995, 12: 411

[56] X. F. Pang, Chin. J. Atom. Mol. Phys., 1996, 13: 70

[57] X. F. Pang, Chin. J. Atom. Mol. Phys., 1997, 14: 232

[58] X. F. Pang, Int. J. Infrared. Mill. Waves, 2001, 22: 277

[59] X. F. Pang, Int. J. Infrared. Mill. Waves, 2002, 23: 375

[60] X. F. Pang, Phys. Rev. E, 2000, 62: 6989

[61] X. F. Pang, Euro. Phys. J. B, 2001, 19: 297

[62] X. F. Pang, Commun. Theor. Phys., 2001, 35: 323

[63] X. F. Pang, Commun. Theor. Phys., 2002, 37: 715

[64] X. F. Pang, Commun. Theor. Phys., 2004, 41: 470

[65] X. F. Pang, Commun. Theor. Phys., 2004, 43: 367

[66] X. F. Pang, Int. J. Infrared Mill. Waves, 2001, 22: 291

[67] X. F. Pang, J. Phys. Chem. Solids, 2001, 62: 793

[68] X. F. Pang, Chin. J. BioMed. Engineering, 1999, 8: 39

[69] X. F. Pang, Chin. J. BioMed. Engineering, 2001, 10: 613

[70] X. F. Pang, Phys. Stat. Sol.(b), 2002, 226: 317

[71] X. F. Pang and Y. P. Feng, Quantum Mechanics in Nonlinear Systems, New Jersey: World Science Publishing Co., 2005: 471

[72] X. F. Pang, Soliton Physics, Chengdu: Sichuan Scien. Tech. Press, 2003: 673

[73] X. F. Pang, H. W. Zhang, J. F. Yu, and Y. H. Luo, Inter. J. Mod. Phys. B, 2005, 19: 4677

[74] X. F. Pang, H. W. Zhang, and Y. H. Luo, J. Phys.: Condens. Matt., 2006, 18: 613

[75] X. F. Pang, H. W. Zhang, J. F. Yu, and Y. H. Luo, Int. J. Mod. Phys. B, 2006, 20: 3027

[76] X. F. Pang, J. F. Yu, and Y. H. Luo, Int. J. Mod. Phys. B, 2007, 21: 13

[77] G. Careri, U. Buontempo, F. Galluzzi, A. C. Scott, E. Gratton, and E. Shydmsunder, Phys. Rev. B, 1984, 30: 4689

[78] J. C. Eibeck, P. S. Lomdahl, and A. C. Scott, Phys. Rev. B, 1984, 30: 4073

[79] J. Stiefel, Einfuhrung in die Numerische Mathematik, Stuttgart: Teubner Verlag, 1965

[80] K. E. Atkinson, An Introdution to Numerical Analysis, 2nd Ed., New York: Wiley, 1989

[81] X. F. Pang and M. J. Liu, Commun. Theor. Phys., 2007, 48: 369

Rights & permissions

Higher Education Press and Springer-Verlag 2008