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"
$$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 $$