Let $u v\in W^{1,1}(\Omega)$ and consider the gradient functional $$ J(u) = \int_\Omega |\nabla u(x)| \, dx $$ and the perturbation $u_\varepsilon = u + \varepsilon v$, with $\varepsilon > 0$.
On the set $A = \{ x \in \Omega \mid \nabla u(x) \neq 0 \}$, we have differentiability, and $$ \frac{1}{\varepsilon} \left( |\nabla u + \varepsilon \nabla v| - |\nabla u| \right) \to \frac{\nabla u}{|\nabla u|} \cdot \nabla v \quad \text{as } \varepsilon \to 0^+ $$
On the set $B = \{ x \in \Omega \mid \nabla u(x) = 0 \}$, we have $$ |\nabla u + \varepsilon \nabla v| = \varepsilon |\nabla v| \quad \Rightarrow \quad \frac{|\nabla u + \varepsilon \nabla v| - |\nabla u|}{\varepsilon} = |\nabla v| $$
Therefore, we obtain define the directional derivative $$ \delta J(u)[v] = \lim_{\varepsilon \to 0^+} \frac{J(u + \varepsilon v)-J(u)}{\varepsilon} = \int_{\{\nabla u \ne 0\}} \frac{\nabla u}{|\nabla u|} \cdot \nabla v \, dx + \int_{\{\nabla u = 0\}} |\nabla v| \, dx $$
Question: what about the case $u,v\in BV(\Omega)$
Let now $u,v \in BV(\Omega)$. Recall that the total variation is defined as $$ J(u) = |Du|(\Omega)= \sup \left\{ \int_\Omega u \, \operatorname{div} \phi \, dx : \phi \in C_c^1(\Omega; \mathbb{R}^n),\ |\phi(x)| \leq 1 \right\} $$ where $Du$ is a vector-valued Radon measure.
How we compute $$ \delta J(u)[v] = \lim_{\varepsilon \to 0^+} \frac{|D(u + \varepsilon v)|(\Omega) - |Du|(\Omega)}{\varepsilon} $$
Any reference or answer is welcome.
Recall that for a function $u\in BV(\Omega)$ there is a sequence $u_n\in W^{1,1}(\Omega)$ such that $\|u_n-u\|_{L^1(\Omega)}\to 0$ and $$|D u_n|(\Omega)\to |D u|(\Omega)$$ Moreover, for $u\in W^{1,1}(\Omega)$ we have $u\in BV(\Omega)$ and $$|D u|(\Omega)= \int_\Omega |\nabla u(x)|dx$$
