Gradient of exponential scalar function of a matrix variable 
Write the explicit formula of the gradient of
$$ E(u)= \sum_{i=2}^{n-1} \sum_{j=2}^{m-1} \exp \left( - \frac{(u_{i+1,j} - u_{i,j-1})^2 + (u_{i,j+1} - u_{i-1,j})^2}{2 {\cal E} ^ 2}  \right) $$
with respect to $n \times m$ matrix $u$ , and where ${\cal E} > 0$  is a given
  constant. Show all the steps of your calculations.

Define fourth-order tensors $\Omega_1$ and $\Omega_2$ with components
$$\eqalign{
\Omega_{1ijkl} &= \begin{cases}
+1  &{\rm if}\quad(k = i+1)\;{\rm and}\;(l = j) \\
-1  &{\rm if}\quad(k = i)\;{\rm and}\;(l = j-1) \\
\;\;\,0  &{\rm otherwise}\ \
\end{cases} \\
}$$
and 
$$\eqalign{
\Omega_{2ijkl} &= \begin{cases}
+1  &{\rm if}\quad(k = i)\;{\rm and}\;(l = j+1) \\
-1  &{\rm if}\quad(k = i-1)\;{\rm and}\;(l = j) \\
\;\;\,0  &{\rm otherwise}\ \
\end{cases} \\
}$$
and use it to define the following matrices
$$\eqalign{
A_1 &= \Omega_1:U\quad&\implies A_{1ij}
 = \sum_{k=1}^n\sum_{l=1}^m\Omega_{1ijkl} U_{kl} \\
A_2 &= \Omega_2:U\quad&\implies A_{2ij}
 = \sum_{k=1}^n\sum_{l=1}^m\Omega_{2ijkl} U_{kl} \\
X_1 &= A_1\odot A_1\quad&\implies X_{1ij} = A_{1ij}A_{1ij} \\
X_2 &= A_2\odot A_2\quad&\implies X_{2ij} = A_{2ij}A_{2ij} \\
X_3 &= \frac{-1}{2 {\cal E} ^ 2} (X_1 + X_2) \quad&\implies X_{3ij} = \frac{-1}{2 {\cal E} ^ 2} ( X_{1ij} + X_{2ij}) \\
C &= \exp(X_3)\quad&\implies C_{ij} = \exp(X_{3ij}) \\
\quad&\implies dC = C\odot dX_3 \\
B &= (\frac{-1}{ {\cal E} ^ 2} ((A_1 : \Omega_1 ) + (A_2 : \Omega_2))) \\
}$$
where $(:)$ denotes the double-dot product, $(\odot)$ denotes the element-wise product, and all functions are applied element-wise on their arguments.
$$\eqalign{
m &= {\tt1}:C \\
dm
 &= {\tt1}:dC \\
 &= {\tt1}:(C\odot dX_3) \\
 &= C:dX_3 \\
 &= C:(\frac{-1}{2 {\cal E} ^ 2} (2A_1 \odot dA_1 + 2A_2 \odot dA_2 ) ) \\
 &= C:(\frac{-1}{ {\cal E} ^ 2} (A_1 \odot dA_1 + A_2 \odot dA_2 ) ) \\
 &= C:(\frac{-1}{ {\cal E} ^ 2} (A_1 \odot \Omega_1:dU + A_2 \odot \Omega_2:dU ) ) \\
 &= C:(\frac{-1}{ {\cal E} ^ 2} ((A_1 : \Omega_1 ) + (A_2 : \Omega_2))) \odot dU \\
 &= C:B \odot dU \\
 &= C \odot B : dU \\
\frac{\partial{\cal m}}{\partial U} &= C \odot B \\
}$$
Is it correct? If it is not correct why? Can anybody help me? Tnx
 A: You're almost there. For ease of typing, let's  bury the $\Big(\frac{-1}{2{\cal E}^2}\Big)$ factor in the definition of the $\Omega$ tensors.
Then picking up after the third line of your differential
$$\eqalign{
dm &= C:d(X_1+X_2) \\
   &= C:dX_1 + C:dX_2 \\
   &= C:(2A_1\odot dA_1)         \;+\; C:(2A_2\odot dA_2) \\
   &= (2C\odot A_1):dA_1         \;+\; (2C\odot A_2):dA_2 \\
   &= (2C\odot A_1):\Omega_1:dU \;+\; (2C\odot A_2):\Omega_2:dU \\
\frac{\partial m}{\partial U}
   &= (2C\odot A_1):\Omega_1 \;+\; (2C\odot A_2):\Omega_2 \\
}$$
At this point, pulling the buried factor out recovers a result in terms of your original $\Omega$ tensors
$$\eqalign{
\frac{\partial m}{\partial U}
 &= -\frac{1}{{\cal E}^2}
     \Big(\big(C\odot A_1\big):\Omega_1
        + \big(C\odot A_2\big):\Omega_2\Big) 
\\
}$$
NB: You cannot form products like  $(A_1\odot\Omega_1)$ because the operands in such products must be tensors of the same order and of the same dimension within each order.
So I can't take the Hadamard product of a matrix with a vector, nor can I take the product of a $(3\times 1)$ vector with a $(2\times 1)$ vector. 
The $\Omega$ tensors can be defined in terms of Kronecker delta 
symbols as follows
$$
(\Omega_1)_{ijkl} = \delta_{(i+1)k}\delta_{jl}-\delta_{ik}\delta_{(j-1)l} \\
(\Omega_2)_{ijkl} = \delta_{ik}\delta_{(j+1)l}-\delta_{(i-1)k}\delta_{jl} \\
$$
