Principal Stresses Program - help sought

[trs32596]

Hi there,

Having trouble creating a program to solve the principal stress equation. This is done by finding the roots of the third order polynomial.

Here is an example of input data:
σx = 60MPa σy =-30MPa σz =-20MPa τxy =40MPa Tyz = 0, Tzx = 0.

Using the normal equation:
σ^3 −(σx +σy +σz )⋅σ^2 +(σxσy +σxσz +σyσz −τxy^2 −τyz^2 −τzx^2 )⋅σ −(σxσyσz +2*τxy*τyz*τzx −σx*τyz^2 −σy*τzx^2 −σz*τxy^2)=0
produces roots of -45.21, -20.0, and 75.21MPa.

My goal is to create a simple program to do this with the input data.

Code:
EXPORT PrincipalStressCoeffs(X,Y,Z,L,M,N)
BEGIN
LOCAL C1:= (X+Y+Z);
LOCAL C2:= (X*Y+X*Z+Y*Z-L^2-M^2-N^2);
LOCAL C3:= (X*Y*Z+2*L*M*N-X*M^2-Y*N^2-Z*L^2);
LOCAL r = CAS( proot([1, -C1, C2, -C3]) );
RETURN(r);
END;

This code ends up returning imaginary numbers and not the same roots. Any help would be much appreciated.