5
$\begingroup$

In v11.3 NSolve returns a bunch of redundant solutions for the following set of six equations:

eqns={
  0 == -0.07 s1 + 0.3 (0.75 i1 + s1) (1 - 0.05 (i1 + i2 + i3 + s1 + s2 + s3)),
  0 == -0.28 i1,
  0 == -0.07 s2 - 0.5 i1 s2 - 0.5 i2 s2 - 0.5 i3 s2 + 0.8 (0.75 i2 + s2) (1 - 0.05 (i1 + i2 + i3 + s1 + s2 + s3)), 
  0 == -0.28 i2 + 0.5 i1 s2 + 0.5 i2 s2 + 0.5 i3 s2,
  0 == -0.07 s3 - 0.4 i1 s3 - 0.4 i2 s3 - 0.4 i3 s3 + 0.7 (0.75 i3 + s3) (1 - 0.05 (i1 + i2 + i3 + s1 + s2 + s3)), 
  0 == -0.28 i3 + 0.4 i1 s3 + 0.4 i2 s3 + 0.4 i3 s3};

unks={s1, i1, s2, i2, s3, i3};

eq = NSolve[eqns, unks];
Length[eq]
(* 30 *)

For example, {s1 -> 0, i1 -> 0, s2 -> 18.25, i2 -> 0, s3 -> 0, i3 -> 0} shows up in eq twice, {s1 -> 0, i1 -> 0, s2 -> 0, i2 -> 0, s3 -> 18., i3 -> 0} shows up four times, etc.

In v11.2, NSolve returns twelve unique solutions, as does Solve in v11.3:

eq = Solve[eqns, unks];
Length[eq]
(* 12 *)

Giving NSolve a Method such as "EndomorphismMatrix", "CompanionMatrix", "Legacy", "Aberth", "JenkinsTraub" as discussed here results in the proper set of twelve solutions. Edit: Even non-existent methods give twelve. However Method->"Homotopy" gives 30.

Does v11.3 use a new Method for NSolve that produces this seeming bug? What would the best workaround be?

Edit:

$Version gives "11.3.0 for Mac OS X x86 (64-bit) (March 7, 2018)" so something must have changed since "11.3.0 for Mac OS X x86 (64-bit) (January 22, 2018)" in @MichaelE2's comment.

I've reported it to WRI.

BTW, this example surfaced in a general function I wrote, so any work-around should be generally applicable not just deal with this particular example.

$\endgroup$
7
  • $\begingroup$ I get only 12 in V11.3 MacOS: i.stack.imgur.com/FPM3x.png $\endgroup$
    – Michael E2
    Dec 31, 2018 at 6:34
  • $\begingroup$ If the system is always polynomial, consider Method -> "Legacy". Since NSolve started as a polynomial solver, maybe you're not losing anything. $\endgroup$
    – Michael E2
    Dec 31, 2018 at 15:33
  • $\begingroup$ @MichaelE2 Good idea but unfortunately the system isn't necessarily polynomial. $\endgroup$
    – Chris K
    Dec 31, 2018 at 15:36
  • 2
    $\begingroup$ (1) The sparse homotopy method can overcount the actual number and fail to recognize that several solutions are actually the same. Offhand I do not recall if there is a good way to detect multiplicity vs overcount, other than maybe by difficulty of obtaining convergence for solutions with multiplicity (which of course makes matters difficult, since they might not look the same). $\endgroup$ Jan 1, 2019 at 18:13
  • 2
    $\begingroup$ (2) I believe this method is only implemented in machine precision, hence using anything else, even low arbitrary precision, will cause a different method to be used. (3) Likewise, invoking a nonexistent method will cause the endomorphism matrix method to be used instead. $\endgroup$ Jan 1, 2019 at 18:14

3 Answers 3

3
$\begingroup$

This is an extended comment. Actually, the answers only appear to be the same. For instance, eq[[6]] and eq[[7]] actually are

{{s1 -> 0, i1 -> 0, s2 -> 18.249999999999982`, i2 -> 0, s3 -> 0, i3 -> 0}}
{{s1 -> 0, i1 -> 0, s2 -> 18.249999999988777`, i2 -> 0, s3 -> 0, i3 -> 0}}

More generally, Union@eq does not reduce the number of solutions, indicating that there are small differences among them. Perhaps, NSolve performs its analysis on the equations in various orders, yielding solutions that differ in some cases only by roundoff.

Incidentally,

Solve[Rationalize[eqns, 0], unks]

gives the twelve desired solutions exactly.

$\endgroup$
1
  • $\begingroup$ Yes this issue seems particular to NSolve. $\endgroup$
    – Chris K
    Dec 31, 2018 at 14:58
2
$\begingroup$
$Version

(* "11.3.0 for Mac OS X x86 (64-bit) (March 7, 2018)" *)

Rationalize the equations

eqns = {0 == -0.07 s1 + 
       0.3 (0.75 i1 + s1) (1 - 0.05 (i1 + i2 + i3 + s1 + s2 + s3)), 
     0 == -0.28 i1, 
     0 == -0.07 s2 - 0.5 i1 s2 - 0.5 i2 s2 - 0.5 i3 s2 + 
       0.8 (0.75 i2 + s2) (1 - 0.05 (i1 + i2 + i3 + s1 + s2 + s3)), 
     0 == -0.28 i2 + 0.5 i1 s2 + 0.5 i2 s2 + 0.5 i3 s2, 
     0 == -0.07 s3 - 0.4 i1 s3 - 0.4 i2 s3 - 0.4 i3 s3 + 
       0.7 (0.75 i3 + s3) (1 - 0.05 (i1 + i2 + i3 + s1 + s2 + s3)), 
     0 == -0.28 i3 + 0.4 i1 s3 + 0.4 i2 s3 + 0.4 i3 s3} // Rationalize // 
   Simplify;

unks = {s1, i1, s2, i2, s3, i3};

Use arbitrary precision by setting the WorkingPrecision

eq = NSolve[eqns, unks, WorkingPrecision -> 7]

(* {{i1 -> 0, s1 -> 0, s2 -> -62.9200 - 11.5510 I, i2 -> -15.1800 + 38.4020 I, 
  s3 -> 79.3500 + 14.4388 I, i3 -> 15.2500 - 38.7566 I}, {i1 -> 0, s1 -> 0, 
  s2 -> -62.9200 + 11.5510 I, i2 -> -15.1800 - 38.4020 I, 
  s3 -> 79.3500 - 14.4388 I, i3 -> 15.2500 + 38.7566 I}, {i1 -> 0, s1 -> 0, 
  s2 -> 18.25000, i2 -> 0, s3 -> 0, i3 -> 0}, {i1 -> 0, s1 -> 0, s2 -> 0, 
  i2 -> 0, s3 -> 18.00000, i3 -> 0}, {i1 -> 0, s1 -> 14.21852, s2 -> 0, 
  i2 -> 0, s3 -> 0.700000, i3 -> 0.4148148}, {i1 -> 0, s1 -> 15.33333, 
  s2 -> 0, i2 -> 0, s3 -> 0, i3 -> 0}, {i1 -> 0, s1 -> 0, s2 -> 0, i2 -> 0, 
  s3 -> 0.700000, i3 -> 9.41499}, {i1 -> 0, s1 -> 14.30667, s2 -> 0.560000, 
  i2 -> 0.4666667, s3 -> 0, i3 -> 0}, {i1 -> 0, s1 -> 0, s2 -> 0.560000, 
  i2 -> 10.60545, s3 -> 0, i3 -> 0}, {i1 -> 0, s1 -> 0, s2 -> 0, i2 -> 0, 
  s3 -> 0.700000, i3 -> -1.714995}, {i1 -> 0, s1 -> 0, s2 -> 0.560000, 
  i2 -> -1.245448, s3 -> 0, i3 -> 0}, {i1 -> 0, s1 -> 0, s2 -> 0, i2 -> 0, 
  s3 -> 0, i3 -> 0}} *)

Verifying solutions

And @@ (eqns /. eq // Flatten)

(* True *)

NSolve returned the expected 12 solutions

Length[eq]

(* 12 *)
$\endgroup$
4
  • $\begingroup$ Thanks! It seems that the WorkingPrecision->7 alone does the trick; Rationalize alone results in 31 solutions. $\endgroup$
    – Chris K
    Dec 31, 2018 at 14:57
  • $\begingroup$ @ChrisK You probably want WorkingPrecision -> $MachinePrecision (~ 15.95) or higher, unless that doesn't work. Not sure why Bob chose 7, which is less than single-precision. $\endgroup$
    – Michael E2
    Dec 31, 2018 at 15:37
  • $\begingroup$ @MichaelE2 - WorkingPrecision -> 7] was chosen to show that even low-precision arbitrary-precision is better than machine precision. If you remove the WorkingPrecision option, And @@ (eqns /. eq // Flatten) evaluates to False, i.e., the results are not sufficiently precise to demonstrably satisfy the equations. Any specified WorkingPrecision above 7 (including $MachinePrecision) is clearly more precise; however, it will clutter the output. $\endgroup$
    – Bob Hanlon
    Dec 31, 2018 at 15:55
  • 1
    $\begingroup$ Thanks. But I'd assert that the machine precision calculations are more precise than the WorkingPrecision -> 7 ones. The reason you get True is that the rounding error estimate is taken into consideration, given that the initial numbers are accurate to only 7 digits. The maximum of the residual errors (i.e. the right-hand sides of the OP's eqns) at machine precision is 2.37442*10^-9, whereas for 7 digits, Mathematica estimates a bound of 3.16306 and the actual error is 0.655472, which results in an arbtrary-precision zero of about 0``-0.5. $\endgroup$
    – Michael E2
    Dec 31, 2018 at 16:43
1
$\begingroup$

Adding a non-existent Method->"Foo" as here fixes the problem, as does Method->"EndomorphismMatrix".

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.