Solution (source code)

= Solution

In the <Voigt elastic approximation>, every local <strain> is $E$. The corresponding local <stress> is $\sigma_V(x)=D_\varepsilon W(E,x)$. Differentiating the averaged <strain energy density> under the integral gives
$$
\boxed{DW_V(E)=\langle D_\varepsilon W(E,x)\rangle=\langle\sigma_V\rangle.}
$$
The uniform <strain> need not generate an equilibrated local <stress>; this calculation concerns the stated approximation, not an assertion that it solves the heterogeneous boundary problem.

For the <Reuss elastic approximation>, use the <convex conjugate>
$$
W^*(S,x)=\sup_e\{S:e-W(e,x)\}.
$$
Under differentiable strictly convex duality, the local <strain> produced by uniform <stress> $S$ is $e_R(x)=D_SW^*(S,x)$. Set $F(S)=\langle W^*(S,x)\rangle$. Averaging the local <strains> gives $E=DF(S)$. At the maximizing <stress> in $W_R(E)=F^*(E)$, the stationarity condition is precisely $DF(S)=E$. Consequently
$$
\begin{aligned}
DW_R(E):H
&=D_E\bigl(S(E):E-F(S(E))\bigr):H\\
&=S:H+DS(E)[H]:\bigl(E-DF(S)\bigr)=S:H.
\end{aligned}
$$
Thus \b[the derivative of each approximate energy is its own mean stress], including
$$
\boxed{DW_R(E)=S=\langle\sigma_R\rangle.}
$$
For positive <linear elasticity>, $W=e:C(x)e/2$ and $W^*=S:C(x)^{-1}S/2$, giving the useful check $C_V=\langle C\rangle$ and $C_R=\langle C^{-1}\rangle^{-1}$. More general convex duality gives the same conjugacy relation with <subgradients> when differentiability or strict convexity is absent.