variational derivative Let $\Omega \subset \mathbb{R}^n,\ n=1,2 \mbox{ or } 3$.
Define the following energy
$$E=\int_{\Omega} \frac{1}{\varepsilon}\left[f(u)+\frac{\varepsilon^2}{2}|\gamma(n)\nabla u|^2\right]\,dx$$
where, 
$$n=\frac{\nabla u}{|\nabla u|}$$
this problem talking about the derivation of Cahn-Hilliard System, $\gamma$ defines the anisotropy function, we take it as general as we want e.g. 
$$\gamma(n)=1+\Gamma(n)$$ 
and $f \in C^3(\Omega)$.
the point is that the article talk about the variational derivative of $E$ w.r.t. $u$ and he get the following
$$\frac{\delta E}{\delta u}=f'(u)-\varepsilon^2 \nabla \cdot \textbf{m}$$
where,
$$\textbf{m}= \gamma^2 \textbf{p}+\gamma\ |\textbf{p}|\ \textbf{P} \nabla_n \gamma$$
with $$\textbf{p}=\nabla u\ ,\ \textbf{P}=|\nabla u|\ \frac{\partial n}{\partial \textbf{p}}$$ 
i start working in the variational derivative i get the following 
\begin{eqnarray*}
\delta E &=&\int_{\Omega} \frac{1}{\varepsilon}\left[f'(u) \delta u+\varepsilon^2 (\gamma(n) \nabla u)\delta(\gamma(n) \nabla u)\right]\,dx\\\\
&=& \int_{\Omega}\frac{1}{\varepsilon} \left[ f'(u) \delta u+\varepsilon^2(\gamma(n) \nabla u)\left(\delta(\gamma(n)) \nabla u+\gamma(n)\delta(\nabla u)\right)\right]\,dx
\end{eqnarray*}
my problem is that i don't know how to continue to get $\textbf{m}$. any help?
$\textbf{Remark :}$ This problem mentioned in "S. Torabi, S. Wise, J. Lowengrub , A New Method for Simulating Strongly Anisotropic Cahn-Hilliard Equations"
 A: $$E=\int_{\Omega} \left[f(u)+\frac{\varepsilon^2}{2}|\gamma(n)\nabla u|^2\right]dx$$
The variational/functional/Gateaux derivative is given by
$$\delta E(u)(v)=\frac{d}{d\delta} E(u+\delta v)|_{\delta=0}$$
for all $v \in C_c^1(0,T)$. Hence
$$E(u+\delta v) = \int_\Omega  f(u+\delta v) + \frac{\epsilon^2}{2} 
\left|\gamma(n(u+\delta v))(\nabla u + \delta \nabla v) \right|^2 dx$$
and taking the derivative
$$\frac{d}{d\delta}E(u+\delta v) =  \int_\Omega  f'(u+\delta v)v + \epsilon^2 
\gamma(n(u+\delta v)) (\nabla u + \delta\nabla v ) \Big(\gamma(n(u+\delta v))\nabla v + \gamma'(n(u+\delta v)) n'(u+\delta v) v (\nabla u + \delta \nabla v)\Big)dx$$
and setting $\delta=0$
$$\frac{d}{d\delta}E(u+\delta v)|_{\delta=0} =  \int_\Omega  f'(u)v + \epsilon^2 
\gamma(n(u)) \nabla u \Big(\gamma(n(u))\nabla v + \gamma'(n(u)) n'(u) v \nabla u\Big)dx$$
and integration by parts
$$\frac{d}{d\delta}E(u+\delta v)|_{\delta=0} =  \int_\Omega  f'(u)v - \epsilon^2 
\nabla \cdot \Big(\gamma(n(u)) \nabla u \gamma(n(u))\Big) v  + \epsilon^2 
\gamma(n(u)) \nabla u \gamma'(n(u)) n'(u)  \nabla u  v dx$$
and by the fundamental theorem of calculus of variations
$$\frac{d}{d\delta}E(u+\delta v)|_{\delta=0} =  f'(u) - \epsilon^2 
\nabla \cdot \Big(|\gamma(n(u))|^2 \nabla u \Big)   + \epsilon^2 
\gamma(n(u))  \gamma'(n(u)) n'(u)  |\nabla u|^2     $$
