Message Boards Message Boards


What algorithm does NSolve use for real solutions of underdetermined system

Posted 11 months ago
11 Replies
2 Total Likes

Hi, I am trying to write up a paper with some results form Mathematica (11.2) among other things. I need to know what algorithm does it use when I solve a system of polynomial equations which is underdetermined (infinitely many solutions). In particular, sometimes Mathematica also find real solutions that are embedded into the complex positive-dimensional components. I couldn't find a mentioned of the specific algorithm that Mathematica is using. I believe this code is written by @DanielLichtblau ? I know this discussion: Also, I have looked into NSolve's documentation, but it does not mention the specific algorithm for underdetermined polynomial systems and in particular for finding real solutions embedded in complex curves (though Mathematica certainly find these!). e.g.

In= NSolve[{x^2 + y^2 - 1} == 0, {x, y}, Reals]
Out= {{x -> -0.997714, y -> 0.0675708}, {x -> -0.335987, y -> -0.941867}}

While this example may be easy enough, Mathematica also finds real solutions for larger systems.

I also wonder if Mathematica finds ALL the real solutions from all the complex curves? Again, to know this, one needs to know what algorithm/method does Mathematica uses.


11 Replies

For sparse linear systems, Solve and NSolve use several efficient numerical methods, mostly based on Gauss factoring with Markowitz products (approximately 250 pages of code). For systems of algebraic equations, NSolve computes a numerical Gröbner basis using an efficient monomial ordering, then uses eigensystem methods to extract numerical roots.

Posted 11 months ago

@FrankKampas, thanks, but that does not address my specific question.

I get this message when I run your code:

NSolve::infsolns: Infinite solution set has dimension at least 1. Returning intersection of solutions with
-((92291 x)/87992)-(121001 y)/175984 == 1.

It is omitted from your quoted result. Did you not get it? Obviously, NSolve adds linear conditions until the system is not underdetermined. (I say, obviously, because that's what the message says.) With that in mind, Frank's response seems to be exactly what is done to find the solutions (at least the Gröbner basis part).

For the following equation, it adds two hyperplanes:

NSolve[x^2 + y^2 + z^2 - 1 == 0, {x, y, z}, Reals]

By the way, you can turn off that method and get parametrized solutions (at least for many algebraic equations):

SetSystemOptions["NSolveOptions" -> "UseSlicingHyperplanes" -> False]
NSolve[x^2 + y^2 + z^2 - 1 == 0, {x, y, z}, Reals]

NSolve::svars: Equations may not give solutions for all "solve" variables.

Out[] = {
 {z ->  ConditionalExpression[-1.` Sqrt[1.` - 1.` x^2 - 1.` y^2],
   -1.` < x < 1.` && -1.` Sqrt[1.` - 1.` x^2] <= y <= Sqrt[1.` - 1.` x^2]]},
 {z -> ConditionalExpression[Sqrt[1.` - 1.` x^2 - 1.` y^2],
   -1.` < x < 1.` && -1.` Sqrt[1.` - 1.` x^2] <= y <= Sqrt[1.` - 1.` x^2]]},
 {x -> -1.`, y -> 0, z -> 0},
 {x -> 1.`, y -> 0, z -> 0}}
Posted 10 months ago

@MichaelRogers, Thanks for your response. No, I did not get that error, though I used to get that error in earlier versions of Mathematica (now I have 11.2). Regardless of that error, I think technically these algorithms shouldn't give real isolated solutions for underdetermined systems. The Groebner basis method is defined over complex variables and gives complex solutions. If the original system has only isolated solutions, the Groebner Basis method gives all isolated complex solutions which would include isolated real solutions as well. So all is good.

However, if the original system has infinitely many solutions (underdetermined), and we add as many random hyperplanes (or linear conditions, as you say) as required to completely determine the system (i.e., only isolated solutions exist), the Groebner basis method will then give us complex isolated solutions for the new and bigger system. That's because you are adding complex hyperplanes (remember, the GB mehod is strictly defined over complex variables).

The fact that the NSolve command gives us real isolated solutions back is suggesting that it may be using some additional method/trick to extract the real isolated solutions from complex solution curves of the original system. So, my question is much more specific than Frank's response. In any case, Mathematica documentation should at least tell us what method it is using, rather than its user figuring out things. Otherwise, I can't use the results generated by it into a scientific publication. Again, I don't need any proprietary information here obviously. I just need to know if the method makes sense and kosher.

What @MichaelRogers wrote is correct at least for the example under discussion. Linear conditions (one in this example) are added to make the system exactly determined for the highest dimensional components. In this instance the intersection happens to contain real valued solutions.

Are you saying there are complex solutions that are omitted? The complex solutions of the system constructed are real numbers. Why is that surprising?

To see complex isolated solutions, try

NSolve[x^2 + y^2 == -1, {x, y, z}]
Posted 10 months ago

@DanielLichtblau, Thanks for your message! For this particular example, it would be straightforward, sure. That's why I had made a comment in my original question "While this example may be easy enough, Mathematica also finds real solutions for larger systems." How does Mathematica find isolated real solutions out of underdetermined system of polynomial equations when the number of equations/variables are larger? The GB method assumes complex variables. For underdetermined systems, it is very tricky to extract real isolated solutions which would now be embedded in complex positive dimensional components. I can send you a more complicated problem if you want. Or I think you would already know such examples anyway.

It does not take much work to find an example for which no real solutions are found. A simple modification of the example under scrutiny will suffice.

In[193]:= NSolve[{x^2 + y^2 - 1/100} == 0, {x, y}, Reals]

During evaluation of In[193]:= NSolve::infsolns: Infinite solution set has dimension at least 1. Returning intersection of solutions with -((92291 x)/87992)-(121001 y)/175984 == 1.

Out[193]= {}

All solutions found using a hyperplane had nontrivial imaginary parts.

I am not sure what is the goal here, but if it is to determine whether or not an underdetermined system has real solutions, this approach using NSolve is not going to do that.

One approach you might try is to find critical points for some distance-like function, maybe a random weighted sum of squares or distance to a random point (the randomization being useful to avoid degeneracy situations). Greg Reid and coauthors have done work along these lines. Below are a few links that might be of use.

Posted 10 months ago

@DanielLichtblau, perfect! That's the answer I was looking for. My goal is indeed to find all real solutions of underdetermined polynomial systems. Since NSolve was giving some results for certain cases, I got confused. But since you confirmed that NSolve is not going to do that for me, that's all I wanted to know. A suggestion though, it may be good to mention this in the Mathematica documentation so that people don't get confused.

Thanks again.

Posted 10 months ago

@MichaelRogers, No, I am not saying that. Nor it is surprising that the NSolve command solves this small and well-known system accurately when Reals is specified or not. The main concern is as follows. The system x^2 + y^2 - 1 ==0 is not the 'same' if x, y \in C as opposed to x, y \in R. The NSolve command has to use different methods based on if the option 'Reals' is specified. If this option is not specified, then I can understand that the NSolve command uses the GB method and then it would spit out some representation of the complex positive dimensional solution space. But if 'Reals' is specified, the GB method can't be used as it is. I want to know what method/algorithm would NSolve use for underdetermined systems when Reals is specified.

Posted 10 months ago


 NSolve[x^2 + y^2 - 1 == 0, {x, y}, Reals],
 TraceInternal -> True

It does use GroebnerBasis. I do not see why NSolve cannot use complex methods and then verify "all variables, parameters, constants, and function values are restricted to be real."

Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
or Discard

Group Abstract Group Abstract