H. Numerical calculation of the FDC

The calculation of the FDC has to be performed numerically due to the length of the DNA sequence. The computation of the DNA energy (Eq. 3.17) is a simple sum of terms that extends over all the base pairs of the molecular construct. So this large summation must be done numerically. Apart from that, the exponential terms that enter the partition function have a wide range of orders of magnitude and they must be treated correctly.

The partition function is a function of the distance (i.e., the total extension) $ x_\mathrm {tot}$, which is a variable that has to be discretized in order to calculate the value of $ Z$ at each position. The distance is divided into $ m$ equidistant points $ x_i$ $ i=0,\dots,m$ separated by a distance $ \Delta x$ (see Fig. H.1). A value of $ \Delta x=5$ nm is enough for our calculations. The details of the calculated FDC are missed for higher values of $ \Delta x$ and lower values of it do not improve the calculation.

Figure H.1: Discretization of the distance.
\includegraphics[width=8cm]{figs/appendix3/discretization.eps}

Now, for each value of $ x_i$ we have to find the minimum of Eq. 3.22 with respect to $ n$. In order to do this, we have to solve Eq. 3.23 for all values of $ n$ according to $ x_i=x_\mathrm{tot}(f,n)$, compute the energies and get the value $ n_i^*$ that minimizes the energy at fixed $ x_i$. Equation 3.23 is a trascendental equation and it can be solved numerically by using the Newton's method. This method gives the equilibrium force $ f$ of the system after few interations (no more than 10 to obtain a numerical estimation of $ f$ with a relative error smaller than 10$ ^{-5}$) of the following expression

$\displaystyle f_\mathrm{new}=f_\mathrm{old}-\frac{x_\mathrm{tot}(f,n)-x_i}{x'_\mathrm{tot}(f,n)}$ (H.1)

where $ f_\mathrm{new}$ and $ f_\mathrm{old}$ are the previous and the posterior values obtained from the iteration, respectively; $ x_\mathrm{tot}(f,n)$ is Eq. 3.23 evaluated at $ f$ and $ n$; and $ x'_\mathrm{tot}(f,n)=\frac{\partial}{\partial f} x_\mathrm{tot}(f,n)$ is its derivative with respect to $ f$. An initial guess of $ f=0$ ensures the convergence of the solution because the first and second derivatives of function $ x_\mathrm {tot}$ have positive values for all values of $ f>0$, which fulfills the Newton's method conditions (see Fig. H.2).

Figure H.2: Newton's method. Red curve shows the approximate elastic response of the system for $ n=4000$ and trap stiffness $ k=80$ pN$ \cdot \mu $m$ ^{-1}$ according to $ x_\mathrm {tot}(f,n)-x_i$. The value of $ x_i=4000$ determines the origin ordinate. The first iteration of the Newton's method starts at ( $ f,x_\mathrm {tot}-x_i$)=($ 0,-4000$) and after some steps it converges to the solution ( $ f,x_\mathrm {tot}-x_i$)=($ 14,0$). Since the first and second derivatives of $ x_\mathrm {tot}(f,n)-x_i$ are monotonically decreasing for any value of $ x_i>0$, $ f>0$ and $ n$, the shape of the function is very similar to the one depicted here. Therefore, the evolution of the iterations has a similar pattern and the Newton's method always converges when starting from $ f=0$.
\includegraphics[width=11cm]{figs/appendix3/newton.eps}

Once we have the combination of values ($ x_i,n_i^*$) we can now calculate the minimum energy $ G_\mathrm{min}(x_i)=G_\mathrm{tot}(x_i,n_i^*)$ according to Eq. 3.22. The discretized Eq. 3.27 can be written as

$\displaystyle f(x_i)$ $\displaystyle =$ $\displaystyle -k_B T \, \frac{\ln Z(x_{i+1})-\ln Z(x_i)}{x_{i+1}-x_i}$  
  $\displaystyle =$ $\displaystyle -k_B T \, \frac{\ln Z(x_{i+1})-\ln Z(x_i)}{\Delta x}$  
  $\displaystyle =$ $\displaystyle -\frac{k_B T}{\Delta x} \ln \frac{Z(x_{i+1})}{Z(x_i)}$ (H.2)

where the partition functions have to be calculated according to

$\displaystyle Z(x_i)=\sum_{n=0}^{N} \exp \left(-\frac{G_\mathrm{tot}(x_i,n)}{k_BT}\right)$ (H.3)

The problem with this equation is that the values of $ G_\mathrm{tot}(x_i,n)$ become larger as the value of $ x_i$ increases. When these energies are introduced in the exponential, the value of $ Z(x_i)$ can be orders of magnitude smaller than 1. Therefore the quotient $ \frac{Z(x_{i+1})}{Z(x_i)}$ of Eq. H.2 is numerically less accurate for large values of $ x_i$. The result is a loss of accuracy in the calculation of the FDC as $ x_i$ increases.

This problem can be fixed by taking advantage of the calculated value of $ G_\mathrm{min}(x_i)$. The idea consists in rewriting the partition function with all the energies referred to this state. So we have,

$\displaystyle Z(x_i)$ $\displaystyle =$ $\displaystyle \sum_{n=0}^{N} \exp \left(-\frac{G_\mathrm{tot}(x_i,n)}{k_BT}\right)$  
  $\displaystyle =$ $\displaystyle \sum_{n=0}^{N} \exp \left(-\frac{G_\mathrm{tot}(x_i,n)-G_\mathrm{min}(x_i)+G_\mathrm{min}(x_i)}{k_BT}\right)$  
  $\displaystyle =$ $\displaystyle \exp \left(-\frac{G_\mathrm{min}(x_i)}{k_BT}\right) \sum_{n=0}^{N} \exp \left(-\frac{G_\mathrm{tot}(x_i,n)-G_\mathrm{min}(x_i)}{k_BT}\right)$  
$\displaystyle Z(x_i)$ $\displaystyle =$ $\displaystyle e^{\displaystyle -\beta G_\mathrm{min}(x_i)} \cdot \tilde{Z}(x_i)$ (H.4)

where

$\displaystyle \tilde{Z}(x_i)=\sum_{n=0}^{N} \exp \left[-\beta\big(G_\mathrm{tot}(x_i,n)-G_\mathrm{min}(x_i)\big)\right]$ (H.5)

is the partition function calculated with an energy offset given by the state of minimum energy at fixed $ x_i$. Equation H.4 can be introduced into Eq. H.2 to obtain
$\displaystyle f(x_i)$ $\displaystyle =$ $\displaystyle -\frac{k_B T}{\Delta x} \ln \frac{Z(x_{i+1})}{Z(x_i)}$  
  $\displaystyle =$ $\displaystyle -\frac{k_B T}{\Delta x} \ln \frac{\displaystyle e^{\displaystyle ...
...1})}{\displaystyle e^{\displaystyle -\beta G_\mathrm{min}(x_i)} \tilde{Z}(x_i)}$  
  $\displaystyle =$ $\displaystyle -\frac{k_B T}{\Delta x}\left[ \ln \frac{e^{\displaystyle -\beta G...
...a G_\mathrm{min}(x_i)}} + \ln \frac{\tilde{Z}(x_{i+1})}{\tilde{Z}(x_i)} \right]$  
$\displaystyle f(x_i)$ $\displaystyle =$ $\displaystyle -\frac{k_B T}{\Delta x}\left[ \beta \Big( G_\mathrm{min}(x_i)-G_\mathrm{min}(x_{i+1}) \Big) + \ln \frac{\tilde{Z}(x_{i+1})}{\tilde{Z}(x_i)} \right]$ (H.6)

where $ G_\mathrm{min}(x_i)$, $ G_\mathrm{min}(x_{i+1})$, $ \tilde{Z}(x_{i+1})$ and $ \tilde{Z}(x_i)$ can be easily calculated numerically.

To sum up, the numerical calculation of the FDC requires the following steps,

  1. Discretize the distance with equidistant values $ x_i$.
  2. For each value of $ x_i$, calculate the energies $ G_\mathrm{tot}(x_i,n)$ for all the values of $ n$ by solving Eq. H.1.
  3. Among the previously calculated energies, find the value of $ n_i^*$ that minimizes the energy for each value of $ x_i$. This energy is called $ G_\mathrm{min}(x_i)$.
  4. For each value of $ x_i$, calculate $ \tilde{Z}(x_i)$ according to Eq. H.5
  5. For each value of $ x_i$, calculate $ f(x_i)$ according to Eq. H.6.

JM Huguet 2014-02-12