Solve returns different answers to the same question

I am trying to use solve() to solve a system of 4 equations for 4 variables, using the following code.

syms x1 x2 x3 x4
vars = [x1 x2 x3 x4];
eqns = [ x1 == 30/7 - 1/(x3+x4), x2==30/7 - 1/(x3+x4), x3== 10/3 - 1/(x1+x2+x4), x4==10/3 - 1/(x1+x2+x3) ];
C = solve(eqns, vars);

Then vpa(C.x1) returns:

ans =
                                           -1.6666666666666666666666666666667
                                           -1.6666666666666666666666666666667
   4.1316999594422134359804216941161 - 2.5330953622691428389795980612384e-39i
    4.343827578149360853360217735709 + 9.3940138280508998241066580657853e-40i
 0.095901033836997139230789141603439 + 1.6227020204672685655672638821872e-39i

Which is the correct solution. However, when the same problem is reformulated using:

eqns = [ x1 == 10/3 - 1/(x2+x3+x4), x2==10/3 - 1/(x1+x3+x4), x3== 30/7 - 1/(x1+x2), x4==30/7 - 1/(x1+x2) ];

Note that although the names of the variables and the order they are entered is different, this system is identical to the original. Then vpa(C.x4) (which should contain the same solutions as C.x1) instead returns:

ans =
  -1.6666666666666666666666666666667
  -1.6666666666666666666666666666667
    4.343827578149360853360217735709
 0.095901033836997139230789141603439
   4.1316999594422134359804216941161

This is presenting a problem for me since somehow it seems to be losing the imaginary portion of each of the complex solutions, only in the second scenario. Why is this happening and how can I fix it so that it is returning the full complex solutions?

Answer 1

[Walter Roberson]

It is possible to demonstrate that one of the roots you are seeing imaginary portions for has a closed form value of

 10/63 * sqrt(3) * sqrt(499) * cos(atan(63/1596949 * sqrt(3) * sqrt(203226157))/3 + pi/6) - 110/63 - 10/63 * sqrt(499) * sin(atan(63/1596949 * sqrt(3) * sqrt(203226157))/3 + pi/6)

That involves two arctan of the same expression 63*sqrt(3)*sqrt(203226157)*(1/1596949) . You can evaluate that expression and see that it is a number between 0 and 1; therefore the arctan is real-valued. Therefore the expression as a whole is real valued.

The imaginary components you see do not exist and are due to round-off error. If you increase the number of digits in the vpa then the imaginary components will get smaller and smaller.

... actually I need to revise that closed form solution, as I was not paying attention to the way the variables get renamed. I will come back to this with a correction; I need to do something else for now.

[Walter Roberson]

Okay, for x4 for the second equation, the three complicated roots that are showing as complex valued turn out to be real valued. When

 A = 2495*cos(2*arctan(63*sqrt(3)*sqrt(203226157)*(1/1596949))*(1/3))*(1/1323)
 B = 2495*sin(2*arctan(63*sqrt(3)*sqrt(203226157)*(1/1596949))*(1/3))*sqrt(3)*(1/1323)
 C = 55*sqrt(499)*cos((1/3)*arctan(63*sqrt(3)*sqrt(203226157)*(1/1596949)))*(1/1323)
 D = 55*sqrt(499)*sin((1/3)*arctan(63*sqrt(3)*sqrt(203226157)*(1/1596949)))*sqrt(3)*(1/1323)
 E = 20/7

then the three solutions are

 [-A - B + C - D + E
  -A + B + C + D + E
   2 * A - 2 * C + E]

As above, the 63*sqrt(3)*sqrt(203226157)*(1/1596949) of the arctan argument is between 0 and 1, so the arctan is real and so the expressions are real.

[Phil]

Thank you so much for this detailed answer.

The script I am writing involves solving a bunch of different "versions" of these equations. Some are identical (like in the example) but many differ in the constant term and the level of interdependence. Just out of curiosity, how likely do you think it is that an equation of this form to has an actual complex solution?

I hesitate to simply take the real portions of the answers since I wouldn't want an actual complex answer to slip through the cracks.