I would like to derive Beltrami–Michell equation. I consider linearised elasticity theory, and assume that the material of interest is an isotropic elastic solid with the constitutive relation $$\tau = \lambda (\text{Tr } \varepsilon)I + 2\mu \varepsilon,$$ and I considered the static problem in which the balance of linear momentum reduces to $$0 = \text{div } \tau + f,$$ where $f =_{def} \rho_R b.$ The goal is to show that the compatibility conditions for the linearised strain $$\text{curl} (\text{curl } \varepsilon)^T = 0,$$ the constitutive relation and the balance of linear momentum imply that stress tensor $\tau$ solves Beltrami–Michell equation.
I used the identity $$0=\text{curl} (\text{curl } \varepsilon)^T = [\Delta \text {Tr }\varepsilon − \text{div} (\text{div } \varepsilon)]I + {\nabla(\text{div } \varepsilon) + [\nabla(\text{div } \varepsilon)]^T } − \nabla(\nabla \text{Tr }\varepsilon) − \Delta \varepsilon,$$ and substituted $$\varepsilon = \frac{1+\nu}{E} \tau + \frac{\nu}{E} \text{Tr }\tau.$$ But I wasn't able to reduce the identity to $$\Delta\tau+\frac{1}{1+\nu}+\nabla (\nabla(\text{Tr } \tau)) = −(\nabla f + \nabla^T f) −\frac{\nu}{1-\nu}(\text{div } f) I.$$
The trick is to use that the trace of a zero tensor is also zero.
Quick definition of operators: $$ \mathrm{div} \, \mathbf{f} = \nabla \cdot \mathbf{f} \, , \\ \mathrm{Div} \, \boldsymbol{\sigma} = \boldsymbol{\sigma} \cdot \nabla \, , \\ \mathrm{D} \, \mathbf{f} = \mathbf{f} \otimes \nabla \, , \\ \mathrm{curl} \, \mathbf{f} = \nabla \times \mathbf{f} \, , \\ \mathrm{Curl} \, \boldsymbol{\sigma} = -\boldsymbol{\sigma} \times \nabla \, , \\ \Delta \boldsymbol{\sigma} = \mathrm{Div} \mathrm{D} \, \boldsymbol{\sigma} \, , \\ \mathrm{hess} \, \lambda = \mathrm{D} \nabla \, \lambda \, . \\ $$
We have static equilibrium, the inverse constitutive equation, and the compatiblity condition $$ -\mathrm{Div} \, \boldsymbol{\sigma} = \mathbf{f} \, , \\ \boldsymbol{\varepsilon} = \dfrac{1}{E} [ (1+\nu) \boldsymbol{\sigma} - (\nu \, \mathrm{tr} \, \boldsymbol{\sigma}) \boldsymbol{I}] \, , \\ \mathrm{Inc} \, \boldsymbol{\varepsilon} = \mathrm{curl} \mathrm{Curl} \,\boldsymbol{\varepsilon} = 0 \, . $$ The incompatability operator can be re-written as $$ \mathrm{Inc} \, \boldsymbol{\varepsilon} = 2 \,\mathrm{sym}(\mathrm{D} \mathrm{Div} \, \boldsymbol{\varepsilon}) - \Delta \boldsymbol{\varepsilon} - \mathrm{hess} (\mathrm{tr} \, \boldsymbol{\varepsilon}) + (\Delta \mathrm{tr} \, \boldsymbol{\varepsilon} - \mathrm{div} \mathrm{Div} \, \boldsymbol{\varepsilon}) \boldsymbol{I} \, . $$ To make things simpler, we multiply the inverse constitutive relation first with $E$ and divide by $(1+\nu)$, before applying the incompatibility operator to get $$ \dfrac{E}{1+\nu} \mathrm{Inc} \,\boldsymbol{\varepsilon} = \mathrm{Inc} \, \boldsymbol{\sigma} - \dfrac{\nu}{1+\nu} \mathrm{Inc} [( \mathrm{tr} \, \boldsymbol{\sigma}) \boldsymbol{I}] = 0 \, . $$ We now observe that clearly $\mathrm{tr}(\mathrm{Inc}\boldsymbol \,{\varepsilon}) = 0$, since $\boldsymbol{\varepsilon}$ is compatible. Using that on the previous equation with the identities $$ \mathrm{tr} (\mathrm{Inc} \, \boldsymbol{\sigma}) = \Delta \mathrm{tr} \, \boldsymbol{\sigma} - \mathrm{div} \mathrm{Div} \, \boldsymbol{\sigma} \, , \\ \mathrm{tr} (\mathrm{Inc} [ (\mathrm{tr} \, \boldsymbol{\sigma}) \boldsymbol{I} ]) = 2 \Delta \mathrm{tr} \, \boldsymbol{\sigma} \, , $$ allows us to derive $$ (1-\nu) \Delta \mathrm{tr}\, \boldsymbol{\sigma} - (1+\nu) \mathrm{div} \mathrm{Div} \, \boldsymbol{\sigma} = 0 \, , \qquad (1-\nu) \Delta \mathrm{tr}\, \boldsymbol{\sigma} = (1+\nu) \mathrm{div} \mathrm{Div} \, \boldsymbol{\sigma} = -(1+\nu) \mathrm{div} \, \mathbf{f} \, . $$ Now, putting it all together we get the Beltrami-Michell equations $$ -\Delta \boldsymbol{\sigma} - \dfrac{1}{1+\nu} \mathrm{hess} (\mathrm{tr} \, \boldsymbol{\sigma}) = 2 \mathrm{sym} \mathrm{D} \, \mathbf{f} + \dfrac{\nu}{1-\nu} (\mathrm{div} \, \mathbf{f}) \boldsymbol{I} \, . \\ $$
For more details check out: https://browse.arxiv.org/pdf/2402.00480.pdf