Stress in residual operator contracted on wrong index

[btalamini]

The functional-based solid mechanics operator uses a total Lagrangian setting, so it uses the first Piola-Kirchhoff stress in the residual calculation. However, in the kernel here the stress is contracted with the shape function derivatives on its first index, which is the wrong one.

Needed: fix this in a dimensionally independent way, so both scalar and vector-valued primal field problems work.