Let $\Omega \subseteq \mathbb{R}^n$ be a nice convex domain (bounded,open,connected, with smooth boundary...)
Let $u:\Omega \to \mathbb{R} \in W^{1,\infty}(\Omega)$.
Define $c_u=\text{ess} \sup_{\Omega}{\|Du\|_{op}}. \, \,$ More precisely, for every $x \in \Omega$ we consider $Du(x)=\big(\partial_1 u(x),\partial_2 u(x),...,\partial_n u(x)\big)$ (all the derivatives are weak of course).
Now we take the operator norm (this is merely the usual Euclidean norm) of $Du(x)$, and so we get a function $x \to \|Du(x)\|_{op} $ from $\Omega$ to $\mathbb{R}$. $\, \,c_u$ is defined to be its essential supremum.
Question: Is it true that every $u \in W^{1,\infty}(\Omega)$ is $c_u$-Lipschitz?
Remarks:
$(1)\,$ The statement clearly holds in the case where $u$ is continuously differentiable everywhere,by the mean value theorem. (It even holds if $u$ is only differentiable everywhere, when we take $c_u$ to be the actual supremum, instead of the essential supremum).
$(2) \,$ It is a well-known theorem that every $u \in W^{1,\infty}(\Omega)$ is $\|Du\|_{\infty}$-Lipschitz where $\|Du\|_{\infty}$ is defined as $\sum_{i=1}^n \text{ess} \sup_{\Omega}|\partial_i u|$.
The classic proofs (see Evans PDE book for instance, theorem 4, pg 279) goes as follows:
Let $\epsilon >0$, and let $\eta_{\epsilon}$ be the usual mollifier, and define $u^{\epsilon}(x)=u * \eta_{\epsilon}(x)= \int_{\Omega} \eta_{\epsilon}(x-y)u(y)$.
It can be proven that $ \partial_i u^{\epsilon}=\partial_i u * \eta_{\epsilon}$, hence (since $\int_{\mathbb{R}^n} \eta_{\epsilon}=1$, and $\partial_i u(y) \le \text{ess} \sup_{\Omega}|\partial_i u|$ a.e)
$|\partial_i u^{\epsilon}(x)| \le \text{ess} \sup_{\Omega}|\partial_i u|$ for all $x \in \Omega$. Thus,
$$ (1) \, \, \|Du^{\epsilon}(x)\|_{op}=\sqrt{\sum_{i=1}^n |\partial_i u^{\epsilon}(x)|^2} \le \sum_{i=1}^n |\partial_i u^{\epsilon}(x)| \le \sum_{i=1}^n \text{ess} \sup_{\Omega}|\partial_i u|=\|Du\|_{\infty}$$
holds for every $x \in \Omega$.
From here, the proof proceeds as follows:
$$ u^{\epsilon}(x) - u^{\epsilon}(y)=\int_0^1 \frac{d}{dt}u^{\epsilon}(tx+(1-t)y)dt=\int_0^1 Du^{\epsilon}(tx+(1-t)y)\big((x-y)\big)dt \Rightarrow$$
$$ |u^{\epsilon}(x) - u^{\epsilon}(y)| \le \int_0^1 \|Du^{\epsilon}(tx+(1-t)y)\|_{op}\cdot\|x-y\|\big)dt \stackrel{(1)}{\le} \|Du\|_{\infty} \cdot\|x-y\|$$
Now letting $\epsilon \to 0$ (recalling $u^{\epsilon} \to u$ a.e), we conclude that
$$ |u(x) - u(y)| \le \|Du\|_{\infty} \cdot\|x-y\|$$
as required.
The key point is that in estimate $(1)$ we had to pass to the essential supremum of each component separately.
If we could prove $\text{ess} \sup_{\Omega}{\|Du^{\epsilon}\|_{op}} \le \text{ess} \sup_{\Omega}{\|Du\|_{op}}$, then the above proof would produce the stronger result.
Since $\partial_i u^\epsilon = \partial_i u * \eta^\epsilon$, you have that $D_v u^\epsilon = D_v u * \eta^\epsilon$ for every vector $v$. Therefore, you have $$|D_vu^\epsilon(x)| \le \text{ess} \sup |D_v u|\le \|Du\|_{op}|v|.$$ By dividing by $|v|$ and taking supremum over $x$, you obtain the wanted bound.