How to extract only real solutions from the result of Solve

Question

For Example, I have the following polynomial equations:

$$ \left\{ \begin{array}{ll} \dfrac{2025 (208 x+5 (y-3446))}{52 (y+90)^2}+\dfrac{300 (8 x-21 y+1920) (y-80)}{7 (x+30)^3}+\frac{300 (80 x-21 z+4500) (z-100)}{7 (x+30)^3}-\dfrac{1521 (-539 x+50 z+15425)}{539 (z+39)^2} =0\\[30pt] \frac{2025 (445 y-208 z+49410) (z-45)}{52 (y+90)^3}-\dfrac{2025 (x-85) (208 x+5 (y-3446))}{52 (y+90)^3}-\frac{300 (8 x-21 y+1920)}{7 (x+30)^2}-\frac{1521 (-539 y+90 z+19680)}{539 (z+39)^2} =0\\[30pt] \frac{2025 (445 y-208 z+49410)}{52 (y+90)^2}+\frac{300 (80 x-21 z+4500)}{7 (x+30)^2}+\frac{1521 (x-25) (539 x-25 (2 z+617))}{539 (z+39)^3}+\frac{1521 (y-30) (539 y-30 (3 z+656))}{539 (z+39)^3}=0 \\ \end{array} \right.$$

The mathematica expression:

{{(300 (1920 + 8 x - 21 y) (-80 + y))/(7 (30 + x)^3) + (
   2025 (208 x + 5 (-3446 + y)))/(52 (90 + y)^2) + (
   300 (4500 + 80 x - 21 z) (-100 + z))/(7 (30 + x)^3) - (
   1521 (15425 - 539 x + 50 z))/(
   539 (39 + z)^2)}, {-((300 (1920 + 8 x - 21 y))/(7 (30 + x)^2)) - (
   2025 (-85 + x) (208 x + 5 (-3446 + y)))/(52 (90 + y)^3) + (
   2025 (49410 + 445 y - 208 z) (-45 + z))/(52 (90 + y)^3) - (
   1521 (19680 - 539 y + 90 z))/(
   539 (39 + z)^2)}, {-((2025 (49410 + 445 y - 208 z))/(
    52 (90 + y)^2)) - (300 (4500 + 80 x - 21 z))/(7 (30 + x)^2) - (
   1521 (-25 + x) (539 x - 25 (617 + 2 z)))/(539 (39 + z)^3) - (
   1521 (-30 + y) (539 y - 30 (656 + 3 z)))/(539 (39 + z)^3)}}

How to extract only those real solutions as new Rule? I mean a general approach to such solutions with mixed real and complex solutions in a rule.

Answer 1

eqns = {(300 (1920 + 8 x - 21 y) (-80 + y))/(7 (30 + x)^3) + (2025 (208 x + 
           5 (-3446 + y)))/(52 (90 + y)^2) + (300 (4500 + 80 x - 21 z) (-100 +
            z))/(7 (30 + x)^3) - (1521 (15425 - 539 x + 
           50 z))/(539 (39 + z)^2) == 
     0, -((300 (1920 + 8 x - 21 y))/(7 (30 + x)^2)) - (2025 (-85 + 
           x) (208 x + 5 (-3446 + y)))/(52 (90 + y)^3) + (2025 (49410 + 
           445 y - 208 z) (-45 + z))/(52 (90 + y)^3) - (1521 (19680 - 539 y + 
           90 z))/(539 (39 + z)^2) == 
     0, -((2025 (49410 + 445 y - 208 z))/(52 (90 + y)^2)) - (300 (4500 + 
           80 x - 21 z))/(7 (30 + x)^2) - (1521 (-25 + x) (539 x - 
           25 (617 + 2 z)))/(539 (39 + z)^3) - (1521 (-30 + y) (539 y - 
           30 (656 + 3 z)))/(539 (39 + z)^3) == 0} // Simplify;

vars = Variables[Level[eqns, {-1}]];

allSolns = Solve[eqns, vars] // Simplify;

Length[allSolns]

47

Since the exact results contain Root objects they are neither explicitly Real nor explicitly Complex.

ro = {z /. allSolns[[47]]};

{Cases[ro, _Real], Cases[ro, _Complex]}

{{}, {}}

The root objects need to be converted to numbers to be explicitly Complex or Real. Machine precision calculations are too inaccurate when converting the root objects. Specifying a precision will require Mathematica to maintain sufficient precision in the calculations.

realSolns = Cases[allSolns,
   _?(FreeQ[N[#, 16], Complex] &)];

Length[realSolns]

7

realSolns // N[#, 16] & // Grid[#, Alignment -> Left] &

Answer 2

Version 9. Specifying a domain may speed things up or may greatly slow things down.

Cases[N[Solve[{...}, {x,y,z}], 30], {x->_Real, y->_Real, z->_Real}] Result in 28 seconds.

N[Solve[{...}, {x,y,z}, Reals], 30] Stopped it after 6 hours with no result.

Answer 3

Setting expr equal to the OP's expression, the equations are given by Flatten@expr == 0.

expr = {{(300 (1920+8 x-21 y) (-80+y))/(7 (30+x)^3)+(2025 (208 x+5 (-3446+y)))/(52 (90+y)^2)+(300 (4500+80 x-21 z) (-100+z))/(7 (30+x)^3)-(1521 (15425-539 x+50 z))/(539 (39+z)^2)},{-((300 (1920+8 x-21 y))/(7 (30+x)^2))-(2025 (-85+x) (208 x+5 (-3446+y)))/(52 (90+y)^3)+(2025 (49410+445 y-208 z) (-45+z))/(52 (90+y)^3)-(1521 (19680-539 y+90 z))/(539 (39+z)^2)},{-((2025 (49410+445 y-208 z))/(52 (90+y)^2))-(300 (4500+80 x-21 z))/(7 (30+x)^2)-(1521 (-25+x) (539 x-25 (617+2 z)))/(539 (39+z)^3)-(1521 (-30+y) (539 y-30 (656+3 z)))/(539 (39+z)^3)}};
eqns = Flatten@expr == 0;

Somewhat unexpectedly the quickest way to get solutions was with NSolve and setting WorkingPrecision. Using the default MachinePrecision was slower.

nrealsol = NSolve[
   eqns, {x, y, z}, Reals, 
   WorkingPrecision -> 16]; // AbsoluteTiming
NSolve[
   eqns, {x, y, z}, Reals]; // AbsoluteTiming
(*
  {2.780016, Null}
  {5.659495, Null}
*)

Both were faster than Solve without a domain specification. Solve[eqns, {x, y, z}, Reals] took longer than I was willing to wait, which was much less than the six hours of Bill. We can get the real solution by picking the ones for which {x, y, z} is an element of the Reals. Beware, just using N on the exact solutions without specifying a desired precision results in too much error.

AbsoluteTiming[
 sol = Solve[Thread[Flatten@expr == 0], {x, y, z}]; 
 realsol = Pick[sol, {x, y, z} \[Element] Reals /. # & /@ sol];
 N@N[realsol, $MachinePrecision] // Sort
 ]
(*
  {14.680846, {{x -> -264.933, y -> -516.956, z -> 176.326},
               {x -> -104.422, y -> 219.102,  z -> -146.78},
               {x -> -78.0687, y -> 74.0697,  z -> 131.459},
               {x -> 28.8908,  y -> 31.8585,  z -> -43.8798},
               {x -> 77.1941,  y -> 125.209,  z -> 506.349},
               {x -> 184.305,  y -> -229.324, z -> -238.648},
               {x -> 618.239,  y -> -1340.34, z -> 616.787}}}
*)

Check:

Sort[{x, y, z} /. N[realsol, 16]] == Sort[{x, y, z} /. nrealsol]
(* True *)

Answer 4

Given the equation :

eqns = {(300 (1920 + 8 x - 21 y) (-80 + y))/(7 (30 + x)^3) + (2025 (208 x + 
         5 (-3446 + y)))/(52 (90 + y)^2) + (300 (4500 + 80 x - 21 z) (-100 +
         z))/(7 (30 + x)^3) - (1521 (15425 - 539 x + 50 z))/(539 (39 + z)^2) == 0,
         -((300 (1920 + 8 x - 21 y))/(7 (30 + x)^2)) - (2025 (-85 + x) (208 x + 5 
         (-3446 + y)))/(52 (90 + y)^3) + (2025 (49410 + 445 y - 208 z) (-45 + z))/
         (52 (90 + y)^3) - (1521 (19680 - 539 y + 90 z))/(539 (39 + z)^2) == 0,
         -((2025 (49410 + 445 y - 208 z))/(52 (90 + y)^2)) - (300 (4500 + 80 x - 21 z))/
         (7 (30 + x)^2) - (1521 (-25 + x) (539 x - 25 (617 + 2 z)))/(539 (39 + z)^3) - 
         (1521 (-30 + y) (539 y - 30 (656 + 3 z)))/(539 (39 + z)^3) == 0} // Simplify;

Using NSolve to obain Real values and Complex having Real parts, thus making filtering simpler:

allSolns = NSolve[eqns, WorkingPrecision -> 8] // Simplify;

Extracting out the Real roots:

Realroots = Cases[allSolns, {_ -> _Real, _ -> _Real, _ -> _Real}] 

{{x -> 618.23891, y -> -1340.3394, z -> 616.78706}, {x -> 184.30453, 
    y -> -229.32371, z -> -238.64759}, {x -> 77.194088, y -> 125.20903, 
    z -> 506.34885}, {x -> -78.068738, y -> 74.069706, 
    z -> 131.45885}, {x -> -104.42229, y -> 219.10230, 
    z -> -146.77983}, {x -> -264.93286, y -> -516.95612, 
    z -> 176.32593}, {x -> 28.890843, y -> 31.858477, z -> -43.879780}}

Answer 5

You can extract the roots and then select real ones:

x^3 - x + 1 == 0 // Solve[#, {x}]& // (x /. #)& // Select[# ∈ Reals &]