Hi, I have had some issues getting reasonable results from a VUMAT of a Burgers model (Kelvin + Maxwell units). I have cut out the section in the VUMAT that causes the problem (attached). So now the model is much simpler (just a Kelvin unit without damage) to find what is wrong. The problem happens during unloading/relaxation otherwise loading step completes perfectly. As shown in the attached figure, mises stress increases during the relaxation period although there are no loads applied on the model during the relaxation (only confining pressure that has been kept constant during the whole analysis). This increase in Mises stress due to the development of some unwanted/parasitic lateral stresses prevents the delayed strain to recover correctly. I have run this model with different types of elements, loads (with and without confining pressure), with and without scaling factor, and always this problem happens but with different intensities. I tried some iterative techniques to find the best current stress (by changing trail stress) for the calculation of the delayed strain, but none have worked. In the attached vumat, I have used the back stress as the trial stress (theta=0) to keep it simple for review as the choice of the trial stress had no effect on the results. The issue becomes much more prominent when there is damage in the model so I cannot neglect it. I will be grateful if anyone could help. The input file and the vumat are attached.

P.S. for the unloading if I drop the load amplitude to a larger value (0.05 -0.2), no additional lateral stresses occur and thus the strain recovers correctly. 

Edit: there were some minor typos in the figure that are corrected

Abaqus UMAT  VUMAT Viscoelastic ​​​​​​​

​​​​​​​​​​​​​​