%% LyX 2.1.3 created this file.  For more info, see http://www.lyx.org/.
%% Do not edit unless you really know what you are doing.
\documentclass[english]{article}
\usepackage[T1]{fontenc}
\usepackage[latin9]{inputenc}
\usepackage{babel}
\usepackage{amsfonts}
\usepackage{amsmath}
\usepackage[squaren]{SIunits}
\usepackage{graphicx}
\usepackage{hyperref}
\usepackage{epsfig}
\usepackage{subfig}
\usepackage{epstopdf}
\DeclareGraphicsExtensions{.eps}

\title{Details on the expressions used in the various terms of the modified Rutherford equation (MRE).\\}

\begin{document}
\author{O. Sauter, H. van den Brand, S. Nowak, M. Preynas}
\date{July 2015}
\maketitle

\section{Introduction}

The modified Rutherford equation (MRE), starting from the Rutherford equation [Rutherford1973], can be written as follows:
\begin{align}
  \label{MRE}
 \frac{dw}{dt} &= \frac{\rho_s}{\tau_R} \biggl[ \sum{\rho_s \,\Delta'_{ii} } \biggr],
\end{align}
where $\rho_s$ represents the radial location, in meters, of the rational surface where $q=m/n$ for the m, n mode considered (m, n, being the poloidal and toroidal mode numbers). $\tau_R$ is defined as [SauterPoP1997]:
\begin{align}
  \label{tauR}
 \tau_R &= \frac{\mu_0\,\rho_s^2}{1.22\,\eta_{neo}},
\end{align}
where $\eta_{neo}$ is the neoclassical resitivity at the rational surface using [SauterBS1999;2002]. The definition of $\rho_s$ will be discussed in the next Section.

The various terms contributing to the island growth rate in Eq. (\ref{MRE}), $\Delta'_{ii}$ correspond to the various contribution: classical ($\Delta'_{class}$), perturbed bootstrap ($\Delta'_{bs}$), current drive ($\Delta'_{cd}$), etc. They are detailed in the following sections.

\section{On the radial variable}

The radial variable has historically been ``$r$'', the minor radius (or averaged minor radius $r=(R_{\text{max}}-R_{\text{min}})/2$) since the various terms have been derived from large aspect ratio calculations. However, when various shapes and realistic aspect ratio tokamaks are being compared, and in particular gradients, we prefer a variable based on flux surface. The closest to a radial variable is $\rho_V \sim \sqrt{V}$, but $\rho_\Phi\sim \sqrt{\Phi}$ is good as well, where $V$ is the volume enclosed in a flux surface and $\Phi$ the toroidal flux. In addition the ``simple'' relation between $q$ and $B_p$ is often used to simplify or modify the various terms, since it is a relation valid for large aspect ratio:
\begin{align}
  \label{qB0Bp}
 q &= \frac{\tilde{\rho}\,\tilde{B}_\phi}{\tilde{R}\,\tilde{B}_p}
\end{align}
Note that all these variables need to be defined in realistic aspect ratio and shape, hence the ' $\tilde{.}$ '. For $\tilde{B}_\phi=B_0$ and $\tilde{R}=R_0$, it is a natural choice. Note that it is important to keep $R_0$ and $B_0$ related, thus $B_0$ is usually the vacuum field at $R=R_0$, to keep $R_0 \cdot B_0=cst$ relevant for the case studied. Following [HintonHazeltine1986] we can define:
\begin{align}
 B_{po} &= \frac{1}{R_0}\frac{d\psi}{d\rho_\Phi},
\end{align}
with $\rho_\Phi=\sqrt{\Phi/(\pi B_0)}$. And using $q=(1/2\pi)\;d\Phi/d\psi$ we obtain:
\begin{align}
 q &= \frac{\rho_\Phi\,B_0}{R_0\,B_{po}}
\end{align}
It is often useful to have a normalized radial variable, which value lies between 0 and 1, in particular for interpolations. We therefore also define $\hat{\rho}$:
\begin{align}
\hat{\rho}  &\equiv \rho_{\Phi,norm.} = \frac{\rho_\Phi}{\rho_\Phi(edge)} = \frac{\rho_\Phi}{\rho_{\Phi,a}}
\end{align}
where $_a$ stands for the value at the plasma boundary which has a minor radius noted $a$. The minor radius of the plasma boundary is defined as:
\begin{align}
a  &= \frac{R_{\rm max}-R_{\rm min}}{2},
\end{align}
where $R_{\rm max}$ and $R_{\rm min}$ are the maximum and minimum values of $R$ of the last closed flux surface, or more precisely at $Z=Z_{axis}$. For the modified Rutherford equation (MRE) in particular, but actually in many situations to compare with experimental results, one needs a ``flux surface'' radial variable close to the configuration space minor radius. This is obtained by using either $\hat{\rho}_V=\sqrt{V/V_a}$ or $\hat{\rho}$ multiplied by the last closed flux surface minor radius $a$. Since we want to relate $B_0$ and $B_p$, as we shall see below, it is better to use $\hat{\rho}$ and Eq. (\ref{qB0Bp}). For MHD-related equations, it is useful since the position of the rational surface $q=m/n$ is important and should be between 0 and $a$, thus we define $\bar{\rho}$ as:
\begin{align}
\bar{\rho} &= \hat{\rho} \, a \;[\rm m],
\end{align}
and the value at the rational surface often written $\rho_s$ in published papers which we write here as:
\begin{align}
\bar{\rho}_{mn} &= \hat{\rho}(q=m/n) \, a \;[\rm m].
\end{align}
This allows also to define the flux surface inverse aspect ratio with the same radial variable:
\begin{align}
\epsilon_{mn} &= \frac{\bar{\rho}_{mn}}{R_0}.
\end{align}
In this way we define $B_p$ such that we can replace it with:
\begin{align}
  \label{Bp}
  B_{p} &= \frac{\bar{\rho}\,B_0}{R_0\,q}.
\end{align}
Since all the radial variables ``$\rho$'' are related to the toroidal flux, if not otherwise stated, we omit the subscript ``$_\Phi$'' in the following.

\section{$\Delta_{BS}'$}

The bootstrap term has been combined essentially using [Fitzpatrick1995] and the bootstrap formulas [SauterBS1999;2002]. In addition, following the work by [Poli2002] showing that the effective $w_d$ for ions is much larger than for electrons, we can separate the fraction of bootstrap current due to electrons ($f_{bse}$) from the one due to ions ($1-f_{bse}$). These lead to the equation first published in [SauterPPCF2002, Eq(1)]:
\begin{align}
  \label{rsdeltapbs_1}
  \bar{\rho}_{mn}\Delta_{BS}'=a_2\,\bar{\rho}_{mn}\,\beta_p\,|L_{bs}|\,\frac{L_q}{-L_p}\,\biggl(\frac{f_{bse}\,w}{w^2 + w_{d,e}^2}\;+\;\frac{(1-f_{bse})\,w}{w^2 + w_{d,i}^2} \biggr).
\end{align}
In this equation, $a_2$ is a ``free'' parameter which is known only in large aspect ratio, circular shape. It has been taken around 2.6 since [SauterPoP1997]. This was to match the observed saturated islands, assuming $\bar{\rho}_{mn}\Delta'_0=-m$. Thus they are somewhat related. The term $w_{d,e}$ should be taken as in Eq.(6) of [SauterPPCF2002]:
\begin{align}
 w_{d,e} &= \biggl[5.1\,\biggl(\frac{1}{\epsilon\, s\,n}\biggr)^{(1/2)}\biggr]^{(4/3)} \; \biggl(\frac{\chi_\perp}{\chi_\parallel}\biggr)^{(1/3)},
\end{align}
with $\chi_\parallel=v_{te}R_0L_q/n$ in [$m^2/s$] and $\chi_\perp$ from local value or $a^2/\tau_E$, using a scaling law for energy confinement time. The $w_{d,i}$ is taken as:
\begin{align}
 w_{d,i} &= \sqrt{28}\,\rho_b = \sqrt{28}\,\sqrt{\epsilon}\,\rho_p,
\end{align}
with $\rho_b$ the banana width and $\rho_p$ the ion poloidal Larmor radius $v_{ti}/\omega_{ci}$, where the thermal velocities are taken here with $\sqrt{2}$.

Note that the ion part of the bootstrap current, mainly due to ion density gradient is relatively small, this is why $f_{bse}=1$ is used most of the time. Note also that $L_{bs}\approx L_{31}\approx f_t$ can be used, with $f_t=1.46 \sqrt{\epsilon}$ the large aspect ratio value used in [Fitzpatrick1995].

The value of $\beta_p$ in Eq.(\ref{rsdeltapbs_1}) has been taken either as the global value, especially when comparing with experimental results, or as the local value, in particular when coupled with a transport code. Since there are many definitions for $\beta_p$, it is sometimes better to replace it. Using above equations and $B_{p}$ (Eq.(\ref{Bp})) in $\beta_p$, we have:
\begin{align}
  \label{betapbar}
\beta_p &= \frac{2\mu_0\,p}{B_{p}^2} = \frac{2\mu_0\,p_{mn}\,R_0^2q_{mn}^2}{\bar{\rho}_{mn}^2B_0^2},
\end{align}
which leads to introducing:
\begin{align}
  \label{def_s}
  s=\frac{\rho}{q}\,\frac{dq}{d\rho}.
\end{align}
Thus
\begin{align}
  s = \frac{\rho}{L_q},
\end{align}
which can be introduced in Eq. (\ref{rsdeltapbs_1}) with Eq. (\ref{betapbar}):
\begin{align}
  \label{rsdeltapbs_nobetap}
  \bar{\rho}_{mn}\Delta_{BS}' = a_2\,\frac{2\mu_0R_0^2q_{mn}^2}{s\,B_0^2}\frac{|L_{bs}|\,p_{mn}}{-L_p}\,\biggl(\frac{f_{bse}\,w}{w^2 + w_{d,e}^2}\;+\;\frac{(1-f_{bse})\,w}{w^2 + w_{d,i}^2} \biggr).
\end{align}
It is now interesting to introduce the bootstrap current density itself. In fact the term $p/-L_p$ comes from $j_{bs}=<j_{bs}\cdot B>/B_0 = (p/B_{po})\;(|L_{bs}|/-L_p)$, thus we obtain:
\begin{align}
  \label{jbs}
  \frac{p\,|L_{bs}|}{-L_p} &= \frac{\bar{\rho}\,B_0}{q\,R_0}\;j_{bs},
\end{align}
which we introduce in Eq. (\ref{rsdeltapbs_nobetap})
\begin{align}
  \label{rsdeltapbs_jbs1}
  \bar{\rho}_{mn}\Delta_{BS}'=a_2\,\bar{\rho}_{mn}\,\frac{2\mu_0R_0q_{mn}}{s_{mn}\,B_0}\,j_{bs,mn}\,\biggl(\frac{f_{bse}\,w}{w^2 + w_{d,e}^2}\;+\;\frac{(1-f_{bse})\,w}{w^2 + w_{d,i}^2} \biggr).
\end{align}
It is also sometimes written as:
\begin{align}
  \label{rsdeltapbs_jbs2}
  \bar{\rho}_{mn}\Delta_{BS}' = a_2\,\bar{\rho}_{mn}^2\,\frac{2\mu_0q_{mn}}{s_{mn}\,\epsilon_{mn}\,B_0}\,j_{bs,mn}\,\biggl(\frac{f_{bse}\,w}{w^2 + w_{d,e}^2}\;+\;\frac{(1-f_{bse})\,w}{w^2 + w_{d,i}^2} \biggr),
\end{align}
where again the subscript $_{mn}$ means that the term is evaluated at $q=m/n$. If $f_{bse}$ is evaluated from the actualy profiles, then it should also be the local value $f_{bse,mn}$. The local bootstrap current density $j_{bs,mn}$ represents the equilibrium bootstrap current density, that is the ``maximum'' perturbed current density that can be driven by the presence of an island, if the profiles are fully flattened. The terms in parenthesis reduce this value if the island is not too large and the profiles not fully flattened. If the MRE is coupled to a transport code, which leads to flattened profiles in the presence of an island, then the local $j_{bs}$ from these profiles will be reduced and close to zero. However in the MRE one should use the bootstrap current density one would obtain if there was no island, which we can call $j_{bs0}$. Thus Eq.(\ref{rsdeltapbs_jbs1}) can be written:
\begin{align}
  \label{rsdeltapbs_jbs1b}
  \bar{\rho}_{mn}\Delta_{BS}' = a_2\,\bar{\rho}_{mn}\,\frac{2\mu_0R_0q_{mn}}{s_{mn}\,B_0}\,j_{bs0,mn}\,\biggl(\frac{f_{bse}\,w}{w^2 + w_{d,e}^2}\;+\;\frac{(1-f_{bse})\,w}{w^2 + w_{d,i}^2} \biggr).
\end{align}
In practice, $j_{bs0,mn}$ can be evaluated from a linear interpolation of $j_{bs}$ obtained at one full width away, on each side, from $\bar{\rho}_{mn}$ where the profiles are not altered, thus using:
\begin{align}
  j_{bs0,mn} = \frac{1}{2}\;\bigl[j_{bs}(\bar{\rho}_{mn}-w) + j_{bs}(\bar{\rho}_{mn}+w) \bigl].
\end{align}


%\frac{\mu_{0}}{B_{0}}\frac{q_{mn}}{s_{mn}\epsilon_{mn}}\frac{\bar{\rho}_{tor}^{2}}{w}j_{bs}


\section{$\Delta_{classical}'$}

The classical $\Delta'$ should be calculated from an equilibrium and effective $q$ profile. However it is very sensitive to derivatives of $q$. For simulations using the MRE, another approach is used, assuming the island is relatively large and the ``average'' effect of equilibrium current density leads to a value of:
\begin{align}
  \label{deltap0_m}
  \bar{\rho}_{mn}\Delta_{0}' = - m
\end{align}
It is important to note that this leads to two ``side'' effects. First the saturated island width is essentially related to $\Delta_{BS}'/\Delta_{0}'$, thus the coefficients used in $\bar{\rho}_{mn}\Delta_{BS}'$ would be different if another value for $\Delta_{0}'$ was used. Second, the characteristic increase of the island width with time is also related to this assumption. This is why sometimes, keeping $\bar{\rho}_{mn}\Delta_{0}' = - m$, the value of $\tau_R$ (effective resistive time) needs to be modified to match the experimental observation.
If the island is driven due to an unstable $q$ profile, that is due to a positive $\Delta_{classical}'$ at $w=0$ and if one wants to simulate the full time evolution, from small island size, the term needs to be modified following the stabilizing effect of the modification of the current density by the island [White1982] and the observation as in TCV [Reimerdes2002]:
\begin{align}
  \label{deltap0_alpha}
  \bar{\rho}_{mn}\Delta_{class}' = \bar{\rho}_{mn}\Delta_{0}' - \alpha \frac{w}{\bar{\rho}_{mn}},
\end{align}
where $\bar{\rho}_{mn}\Delta_{0}'$ takes a positive value and $\alpha$ is always positive. From [White1982] we can take as first approximation their lowest order term:
\begin{align}
  \label{alpha}
  \alpha = \frac{m^2}{\hat{\rho}_{mn}^2} - \frac{s}{\hat{\rho}_{mn}} = \frac{m^2}{\hat{\rho}_{mn}^2}\;(1 - \frac{s\,\hat{\rho}_{mn}}{m}) = \frac{m^2\,a^2}{\bar{\rho}_{mn}^2}\;(1 - \frac{s\,\bar{\rho}_{mn}}{m\,a}),
\end{align}
while ensuring that $\alpha$ stays positive (for example taking $\alpha={\rm max}[Eq.(\ref{alpha}),m]$). However this is for $w$ between 0 and the ``classical'' saturation value $\bar{\rho}_{mn}\Delta_{0}' \; \bar{\rho}_{mn}/\alpha$.
At larger $w$, $\Delta_{class}'$ should saturate as well, thus we define:
\begin{align}
  \label{rsdeltap_class}
  \bar{\rho}_{mn}\Delta_{class}' = \bar{\rho}_{mn}\Delta_{0}' - \frac{(\bar{\rho}_{mn}\Delta_{0}' + m)\;w}{\sqrt{w^2 + \bigl((\bar{\rho}_{mn}\Delta_{0}' + m)\,\bar{\rho}_{mn}/\alpha \bigr)^2}},
\end{align}
written such that if $\bar{\rho}_{mn}\Delta_{0}' =  -m$ it is well defined and we recover Eq.(\ref{deltap0_m}), such that we follow Eq.(\ref{deltap0_alpha}) at small island width and we recover Eq.(\ref{deltap0_m}) at large island width.

\section{$\Delta_{GGJ}'$}

From [Lutjens2001] and following the notation of [SauterPPCF2002]
\begin{align}
  \label{rsdeltap_GGJ1}
  \bar{\rho}_{mn}\Delta_{GGJ}' = - \; \bar{\rho}_{mn}\,\frac{6D_R}{\sqrt{w^2 + 0.2 \,w_{d,e}^2}},
\end{align}
where the resistive interchange term $D_R$ can be obtained from CHEASE output ([CHEASE1996]), for example, or approximated by:
\begin{align}
  \label{DR}
  D_R = \frac{\epsilon^2\,\beta_p}{s}\frac{L_q}{-L_p} (1-\frac{1}{q^2})
\end{align}
Thus we get:
\begin{align}
  \label{rsdeltap_GGJ_betap}
  \bar{\rho}_{mn}\Delta_{GGJ}' &= - \; \beta_p\,\frac{6\epsilon^2\,\bar{\rho}_{mn}^2}{q_{mn}^2\,s_{mn}^2(-L_p)}\,\frac{(q_{mn}^2-1)}{\sqrt{w^2 + 0.2 \,w_{d,e}^2}}.
\end{align}
Note that $\beta_p\,L_q/-L_p$ expresses the fact that the term is also proportional to the pressure gradient and indeed it should be kept similar as the bootstrap term, as shown in [SauterPPCF2002] and [Sauter1997], since at large island the GGJ term leads to a ``simple'' reduction of the bootstrap term (proportional to $p'/w$). Therefore one should either use Eq. (\ref{rsdeltapbs_1}) with Eqs. (\ref{rsdeltap_GGJ1} and \ref{DR}). On the other hand, if $j_{bs}$ is introduced, which is useful when coupling with a transport code and when ECCD is involved, using Eq. (\ref{rsdeltapbs_jbs1b}), then one should also introduce it here. But first, as we did for the BS term, we can replace $\beta_p$:
\begin{align}
  \label{rsdeltap_GGJ_nobetap1}
  \bar{\rho}_{mn}\Delta_{GGJ}' &= - \; \frac{12\mu_0R_0^2q_{mn}^2\;p_{mn}}{\bar{\rho}_{mn}\,B_0^2}\,\frac{(D_R/\beta_p)}{\sqrt{w^2 + 0.2 \,w_{d,e}^2}} \\
  \label{rsdeltap_GGJ_nobetap2}
  &= - \; \frac{12\mu_0R_0^2}{B_0^2}\,\frac{\epsilon_{mn}^2}{s_{mn}^2}\,\frac{p_{mn}}{-L_p}\,\frac{(q_{mn}^2 - 1)}{\sqrt{w^2 + 0.2 \,w_{d,e}^2}}
\end{align}
Now we can also introduce $j_{bs}$, following Eqs. (\ref{betapbar}, \ref{jbs}) we first write:
\begin{align}
  \label{rsdeltap_GGJ2}
  \beta_p = \frac{2\mu_0R_0q}{\bar{\rho}\,B_0}\,\frac{-L_p}{|L_{bs}|}\;j_{bs},
\end{align}
and then obtain either:
\begin{align}
  \label{rsdeltap_GGJ3}
  \bar{\rho}_{mn}\Delta_{GGJ}' = - \; \frac{12\mu_0R_0q_{mn}}{B_0}\,(\frac{D_R}{\beta_p})\,\frac{-L_p}{|L_{bs}|}\,\frac{j_{bs0,mn}}{\sqrt{w^2 + 0.2 \,w_{d,e}^2}},
\end{align}
or using Eq. (\ref{DR}):
\begin{align}
  \label{rsdeltap_GGJfinal}
  \bar{\rho}_{mn}\Delta_{GGJ}' = - \; \frac{12\mu_0R_0\,\bar{\rho}_{mn}}{q_{mn}\,B_0}\,\frac{\epsilon^2}{s^2}\,(q_{mn}^2-1)\,\frac{j_{bs0,mn}/|L_{bs}|}{\sqrt{w^2 + 0.2 \,w_{d,e}^2}}.
\end{align}
Note that the term $\epsilon^2$ expresses the fact that the GGJ term is actually proportional to $\beta$ rather than $\beta_p$ and thus is much smaller than the BS term except for tight aspect ratio tokamaks as discussed for example in [Buttery2002]. For MAST-type tokamaks, we might also need to use $D_R$ from CHEASE instead of the approximation.

\section{$\Delta_{CD}'$}

For this term, we shall follow mainly the fitted functions, notations and results given in [Sauter2004] and [DeLazzari2009]:
\begin{align}
  \label{rsdeltap_cd_Icd}
  \bar{\rho}_{mn}\Delta_{cd}' = - \; \frac{16\mu_0\bar{\rho}_{mn}}{\pi s_{mn}\,B_p}\,\frac{I_{cd}}{w_{dep}^2}\,N_{cd}(w/w_{dep},D)\,G(w/w_{dep}, \bar{\rho}_{dep})\,M(w/w_{dep},D),
\end{align}
with $I_{cd}$ the total driven current by the beam considered deposited at $\bar{\rho}_{dep}$ with a deposition width $w_{dep}$. Hence the relation we shall use between $j_{cd}$ and $I_{cd}$:
\begin{align}
  \label{Icd_jcd}
  I_{cd} = \frac{\pi\,\sqrt{\pi}}{2}\,w_{dep}\,\bar{\rho}_{dep}\,j_{cd0}.
\end{align}
Note that in this case using $\rho_\Phi$ would be more consistent, since it is the surface which counts between $j_{cd}$ and $I_{cd}$, however the small differences are included in the free coefficients in front of each terms. The function $N_{cd}(w/w_{dep},D)$ (or $\eta_{aux}$) determines the dependence of $\Delta'_{cd}$ on the island width, the function $G_{cd}(w/w_{dep}, \bar{\rho}_{dep})$ the misalignement with respect to $\bar{\rho}_{mn}$ and $M(w/w_{dep},D)$ the effect of modulation relative to the on-time fraction $D$. We shall assume CW deposition here, thus $D=1$ and $M(w/w_{dep},D)=1$.

We can replace $B_p$ in the same way as before, in order to not depend on its effective definition, yielding:
\begin{align}
  \label{rsdeltap_cd_Icd_noBp}
  \bar{\rho}_{mn}\Delta_{cd}' = - \; \frac{16\mu_0R_0q_{mn}}{\pi s_{mn}\,B_0}\,\frac{I_{cd}}{w_{dep}^2}\,N_{cd}(w/w_{dep})\,G_{cd}(w/w_{dep}, \bar{\rho}_{dep}).
\end{align}
This equations would be used in particular with the relations where $\beta_p$ and $B_p$ have been substituted without introducing $j_{bs}$, thus Eqs. (\ref{rsdeltapbs_nobetap}, \ref{rsdeltap_class}, \ref{rsdeltap_GGJ_nobetap1}/\ref{rsdeltap_GGJ_nobetap2}). If one wants to relate to the driven current density, introducing Eq.(\ref{Icd_jcd}):
\begin{align}
  \label{rsdeltap_cd_jcd_noBp}
  \bar{\rho}_{mn}\Delta_{cd}' = - \; \frac{8\sqrt{\pi}\mu_0R_0q_{mn}}{s_{mn}\,B_0}\,\frac{\bar{\rho}_{dep}\,j_{cd0}}{w_{dep}}\,N_{cd}(w/w_{dep})\,G_{cd}(w/w_{dep}, \bar{\rho}_{dep}),
\end{align}
where $j_{cd0}$ marks the fact that it is not $j_{cd}$ at the $q=m/n$ surface but the value at $\bar{\rho}_{dep}$ since the misalignement effect is already accounted for in $G_{cd}(w/w_{dep}, \bar{\rho}_{dep})$ (defined below).

Since the criteria for NTM stabilization includes either $I_{cd}$, hence Eq. (\ref{rsdeltap_cd_Icd_noBp}), or $j_{cd}/j_{bs}$, we still need to introduce $j_{bs}$ here as well. Using Eq. (\ref{jbs}) we can write:
\begin{align}
  \label{rsdeltap_cd_jcd_jbs}
  \bar{\rho}_{mn}\Delta_{cd}' = - \; \frac{8\sqrt{\pi}\mu_0R_0^2q_{mn}^2}{s_{mn}\,B_0^2}\,\frac{\bar{\rho}_{dep}}{\bar{\rho}_{mn}w_{dep}}\,\frac{p\,|L_{bs}|}{(-L_p)}\,\frac{j_{cd0}}{j_{bs0,mn}}\,N_{cd}(w/w_{dep})\,G_{cd}(w/w_{dep}, \bar{\rho}_{dep}).
\end{align}
The function $N_{cd}(w/w_{dep})$ is taken from [Perkins2003, Sauter2004], assuming CW local deposition:
\begin{align}
  \label{Ncd}
  N_{cd}(w/w_{dep}) = \frac{0.25}{1+\frac{2}{3}\,\bigl(\frac{w}{w_{dep}}\bigr)^2}.
\end{align}
The misalignement term is taken from [DeLazzari2009] (Fig. 2). The misalignement mainly leads to a significant decrease when the deposition is outside the island, thus $|\bar{\rho}_{dep}-\bar{\rho}_{mn}| \le {\rm max}(w,w_{dep})/2$. Therefore, defining $x_{norm}=(\bar{\rho}_{dep}-\bar{\rho}_{mn})/{\rm max}(w,w_{dep})$, one obtains Fig. 2 of [DeLazzari2009] for $G_{cd}$, which can be approximated by: 
\begin{align}
  \label{Gcd}
  G_{cd}(x_{norm},w/w_{dep}) = &(1.+G_{coeff})\;\frac{(1-\tanh[\frac{(0.75.*x_{norm}-0.3)}{0.2}])}{(1-\tanh [\frac{-0.3}{0.2}])+2\,x_{norm}^3}\; \nonumber \\
  -\;&G_{coeff}\,e^{-x_{norm}^2},
\end{align}
with $G_{coeff}=0.6$. The function is constructed such as to be one at $\bar{\rho}_{dep}=\bar{\rho}_{mn}$ and zero at large value of $x_{norm}$. The first function yields the drop due to misalignement, with a small negative value when the deposition is just outside. The second function yields the slower return to zero, as shown in Fig.\ref{fig:Gcd}. A better fit can be obtain and a finer dependence on $w/w_{dep}$ could be recovered as well. However the main effect is the drop from 1 to below 0.2 within $x_{norm} \le 0.5$.
%aaa=0.5;ytest4=(1+aaa)./((1-tanh(-0.25/0.15))+2.*xnorm.^3).*(1-tanh((0.75.*xnorm-0.25)/0.15)) - aaa.*exp(-4.*(xnorm).^2/2.4.^2);
%aaa=0.6;ytest4=(1+aaa)./((1-tanh(-0.3/0.2))+2.*xnorm.^3).*(1-tanh((0.75.*xnorm-0.3)/0.2)) - aaa.*exp(-4.*(xnorm).^2/2.^2);

\begin{figure}
  \center
  \includegraphics[width=8cm] {Gcdtest4.eps}
\caption{$G_{cd}$ from Eq.(\ref{Gcd}).}
\label{fig:Gcd}
\end{figure}


\begin{thebibliography}{10}

\bibitem{Rutherford1973} Rutherford1973: P. H. Rutherford, Phys. Fluids {\bf 16} (1973) 1903.
\bibitem{HintonHazeltine1986} HintonHazeltine1986: F. L. Hinton and R. D. Hazeltine, Rev. Mod. Phys. {\bf 48} (1976) 239.
\bibitem{Fitzpatrick1995} Fitzpatrick1995: R. Fitzpatrick, Phys. Plasmas {\bf 2} (1995) 825.
\bibitem{SauterBS1999;2002} SauterBS1999;2002: O. Sauter, C. Angioni and Y. R. Lin-Liu, {\it Neoclassical conductivity and bootstrap current formulas for general axisymmetric equilibria and arbitrary collisionality regime}, Phys. Plasmas {\bf 6} (1999) 2834; ibidem, {\it Erratum}, {\bf 9} (2002) 5140.
\bibitem{Poli2002} Poli2002: Poli E et al,{\it Reduction of the Ion Drive and $\rho_\theta^*$ Scaling of the Neoclassical Tearing Mode}, Phys. Rev. Lett. {\bf 88} (2002) 075001.
\bibitem{SauterPPCF2002} SauterPPCF2002: O. Sauter {\it et al}, {\it Marginal $\beta$-limit for neoclassical tearing modes in JET H-mode discharges}, Plasma Phys. Control. Fusion {\bf 44} (2002) 1999.
\bibitem{SauterPoP1997} SauterPoP1997: O. Sauter {\it et al}, {\it Beta limits in long-pulse tokamak discharges}, Phys. Plasmas {\bf 4} (1997) 1654.
\bibitem{White1982} White1982: R. B. White and D. A. Monticello, M. N. Rosenbluth, and B. V. Waddell, {\it Saturation of the tearing mode}, Phys. Fluids {\bf 20} (1977) 800.
\bibitem{Reimerdes2002} Reimerdes2002: H. Reimerdes {\it et al}, {\it From Current-Driven to Neoclassically Driven Tearing Modes}, Phys. Rev. Lett. {\bf 88} (2002) 105005.
\bibitem{Lutjens2001} Lutjens2001: H. L\"utjens, J.-F. Luciani and X. Garbet, {\it Curvature effects on the dynamics of tearing modes in tokamaks}, Phys. Plasmas {\bf 8} (2001) 4267.
\bibitem{CHEASE1996} CHEASE1996: H. L\"utjens, A. Bondeson, O. Sauter, {\it The CHEASE code for toroidal MHD equilibria}, Computer Physics Communications {\bf 97} (1996) 219.
\bibitem{Buttery2002} Buttery2002: R. J. Buttery {\it et al},{\it Neoclassical Tearing Physics in the Spherical Tokamak MAST}, Phys. Rev. Lett. {\bf 88} (2002) 125005.
\bibitem{Sauter2004} Sauter2004: O. Sauter, {\it On the contribution of local current density to neoclassical tearing mode stabilization}, Phys. Plasmas {\bf 11} (2004) 4808.
\bibitem{DeLazzari2009} DeLazzari2009: D. De Lazzari and E. Westerhof, {\it On the merits of heating and current drive for tearing mode stabilization}, Nucl. Fusion {\bf 49} (2009) 075002; ibidem, {\it Erratum}, {\bf 50} (2010) 079801.
\bibitem{Perkins2003} F. W. Perkins and R. W. Harvey, Bull. Am. Phys. Soc. {\bf 48} (2003) QP1-025.
  
\end{thebibliography}
\end{document}
