Selection of finite element spaces in plasticity computations

General Discussion
1

Hi, everyone

I’m working on a plastic computation script testGauss.edp (8.5 KB) using FreeFEM.

I use two finite‑element spaces"P1-P0"
func Pk = P1;
func PkDC = P0;

The algorithm is convergent. However, when I switch to P1-P1 or P2-P2
func Pk = P1;
func PkDC = P1;

or
func Pk = P2;
func PkDC = P2;

The iterationis divergence. I suspect this is because the incremental update of StrainPlasticity4 (which stores the plastic strain history) may lose inheritance across Gauss points.

StrainPlasticity4 = [
StrainPlasticity4[0]*(1-logicPlastApex) + ...,
...
];
My questions:
How can I correctly implement higher-order elements (e.g., P2-P2)？Should I keep the plastic strain as P0 even if displacement is upgraded to P2? Or is there a standard way in FreeFEM to store history variables directly at Gauss points using arrays instead of FE functions?

Any suggestions would be greatly appreciated!