3
$\begingroup$

I have two polynomials as function of $wa$ and $wb$ , I am going to show those polynomials.

This is the expression for $GS65$:

GS65=((-1 + wb)^4 (1 + 2 wb) (1 + 4 wb^2))/( 1 + 2 (-1 + 
 wb) wb) - (wa (1 - wb) wb^4 (-3 + 2 wb) (5 + 
  4 (-2 + wb) wb) (-1 + wa^3 (1 - wb)^4 - 6 wa^5 wb^5 + 
  6 wa^6 wb^6 + (1 - wa)^6 (1 - wb)^5 (1 + 5 wb) + 
  wa^2 (1 - wb)^3 (-1 + wa wb + 4 wa^2 wb) + 
  wa^4 (1 - wb) wb^2 (-6 + (-4 + 6 wa) wb + (-5 + 4 wa) wb^2 + 
     wa (5 + 6 wa) wb^3) + (1 - wa)^5 (1 - 
     wb)^2 ((5 + wa) (1 - wb)^4 + 6 wb^2 + 4 (1 - wb) wb^2 + 
     5 (1 + 15 wa) (1 - wb)^2 wb^2 + (1 - wb)^3 (-1 + 
        31 wa wb)) + (1 - wa)^4 (1 - wb) (5 wa (1 - wb)^5 + 
     4 (1 + 4 wa) wb^3 + 30 wa (1 - wb) wb^3 + 
     50 wa (1 + 4 wa) (1 - wb)^2 wb^3 + 
     30 wa (1 - wb)^3 wb (-1 + 6 wa wb) + 
     5 (1 - wb)^4 (-1 + 2 wa + 11 wa wb + 6 wa^2 wb)) + 
  wa^3 (1 - wb)^2 wb (-4 + 4 wa wb + 
     wa^2 wb (6 + 4 wb + 5 wb^2)) + (1 - 
     wa)^3 ((1 + 3 wa + 6 wa^2) wb^4 + 40 wa^2 (1 - wb) wb^4 + 
     50 wa^2 (2 + 3 wa) (1 - wb)^2 wb^4 + 
     10 wa^2 (1 - wb)^5 (1 + 5 wb) + 
     50 wa^2 (1 - wb)^3 wb^2 (-3 + 7 wa wb) + 
     2 wa (1 - wb)^4 (-5 + 3 wa - 25 wb + 25 wa wb + 75 wa wb^2 + 
        75 wa^2 wb^2)) + (1 - wa) wa^2 ((1 + 3 wa) (1 - wb)^4 + 
     wa^2 (5 + wa) wb^6 + 3 wa^2 (1 - wb) wb^4 (-25 + 27 wa wb) + 
     wa (1 - wb)^2 wb (-16 + (-30 + 22 wa) wb + (-50 + 
           34 wa) wb^2 + 5 wa (11 + 15 wa) wb^3) + (1 - 
        wb)^3 (-3 + 7 wa wb + 
        2 wa^2 wb (8 + 15 wb + 25 wb^2))) + (1 - 
     wa)^2 wa^2 (10 wa wb^5 + 10 wa (5 + 3 wa) (1 - wb) wb^5 + 
     25 wa (1 - wb)^2 wb^3 (-8 + 11 wa wb) + 
     2 (1 - wb)^3 (-3 + (-20 + 11 wa) wb + 5 (-10 + 7 wa) wb^2 + 
        25 wa (3 + 4 wa) wb^3) + (1 - wb)^4 (3 + 
        2 wa (3 + 20 wb + 50 wb^2)))))/((-1 + (1 - wa) (1 - wb) + 
  wa wb) (1 + 2 (-1 + wb) wb))

The other polynomial is $GS56$:

GS56=(wa (1 - wb) (-1 + wb)^4 (1 + 2 wb) (1 + 4 wb^2) (-1 + 
 wa^3 (1 - wb)^4 - 6 wa^5 wb^5 + 
 6 wa^6 wb^6 + (1 - wa)^6 (1 - wb)^5 (1 + 5 wb) + 
 wa^2 (1 - wb)^3 (-1 + wa wb + 4 wa^2 wb) + 
 wa^4 (1 - wb) wb^2 (-6 + (-4 + 6 wa) wb + (-5 + 4 wa) wb^2 + 
    wa (5 + 6 wa) wb^3) + (1 - wa)^5 (1 - 
    wb)^2 ((5 + wa) (1 - wb)^4 + 6 wb^2 + 4 (1 - wb) wb^2 + 
    5 (1 + 15 wa) (1 - wb)^2 wb^2 + (1 - wb)^3 (-1 + 
       31 wa wb)) + (1 - wa)^4 (1 - wb) (5 wa (1 - wb)^5 + 
    4 (1 + 4 wa) wb^3 + 30 wa (1 - wb) wb^3 + 
    50 wa (1 + 4 wa) (1 - wb)^2 wb^3 + 
    30 wa (1 - wb)^3 wb (-1 + 6 wa wb) + 
    5 (1 - wb)^4 (-1 + 2 wa + 11 wa wb + 6 wa^2 wb)) + 
 wa^3 (1 - wb)^2 wb (-4 + 4 wa wb + 
    wa^2 wb (6 + 4 wb + 5 wb^2)) + (1 - 
    wa)^3 ((1 + 3 wa + 6 wa^2) wb^4 + 40 wa^2 (1 - wb) wb^4 + 
    50 wa^2 (2 + 3 wa) (1 - wb)^2 wb^4 + 
    10 wa^2 (1 - wb)^5 (1 + 5 wb) + 
    50 wa^2 (1 - wb)^3 wb^2 (-3 + 7 wa wb) + 
    2 wa (1 - wb)^4 (-5 + 3 wa - 25 wb + 25 wa wb + 75 wa wb^2 + 
       75 wa^2 wb^2)) + (1 - wa) wa^2 ((1 + 3 wa) (1 - wb)^4 + 
    wa^2 (5 + wa) wb^6 + 3 wa^2 (1 - wb) wb^4 (-25 + 27 wa wb) + 
    wa (1 - wb)^2 wb (-16 + (-30 + 22 wa) wb + (-50 + 
          34 wa) wb^2 + 5 wa (11 + 15 wa) wb^3) + (1 - 
       wb)^3 (-3 + 7 wa wb + 
       2 wa^2 wb (8 + 15 wb + 25 wb^2))) + (1 - 
    wa)^2 wa^2 (10 wa wb^5 + 10 wa (5 + 3 wa) (1 - wb) wb^5 + 
    25 wa (1 - wb)^2 wb^3 (-8 + 11 wa wb) + 
    2 (1 - wb)^3 (-3 + (-20 + 11 wa) wb + 5 (-10 + 7 wa) wb^2 + 
       25 wa (3 + 4 wa) wb^3) + (1 - wb)^4 (3 + 
       2 wa (3 + 20 wb + 50 wb^2)))))/((-1 + (1 - wa) (1 - wb) + 
 wa wb) (1 + 2 (-1 + wb) wb))

Given the values for $GS65$ and $GS56$, I want to get $wa$ and $wb$. If I use NSolve for $GS65=0.75$ and $GS56=0.25$, Mathematica shows the correct result which is $wa=0.5$ and $wb=0.5$. At this point everything has worked as expected, but when I use NSolve with the values $GS65=0.827531$ and $GS56=0.327531$ Mathematica doesn't get a solution, and I know it should be $wa=0.6$ and $wb=0.5$.

I waited for 24 hours and Mathematica didn't finish. I saw that the Mathematica kernel was using 50% of my computer processing power, but it still didn't finish.

WHat can I do in this case? I guess waiting more time won't work so I don't know what to do.

All these variables ($GS65$, $GS56$, $wa$, $wb$) are probabilities and therefore fall in the range 0 to 1.

This is the code I used to solve the problem.

This works fine:

NSolve[GS65==0.75 && GS56==0.25 && 0<wa<1 && 0<wb<1, {wa,wb}]

But this runs for more than 24 hours:

NSolve[GS65==0.827531 && GS56==0.327531 && 0<wa<1 && 0<wb<1, {wa,wb}]

I know the answer for this case should be $wa=0.6$ and $wb=0.5$.

$\endgroup$
2
  • $\begingroup$ Thanks for your answer. I am going to search information about GroebnerBase. I have used the find tool to check if you got $wa=0.6$ but I couldn't find that solution, and I guess you didn't use the fact that $wa$ is >0 and <1. Do you think using that information you could get the right solution using that Groebner base? $\endgroup$ – George Dec 29 '15 at 12:35
  • $\begingroup$ I edited the answer and explicitly solved 5 of 6 pairs of Groebner bases equations, The last one seems is causing the problem (not finished within few minutes), I think because coefficients in these equations are very large. Fortunatelly, your solution is present in the answer. $\endgroup$ – Acus Dec 29 '15 at 13:33
5
$\begingroup$

The problem seems to be that NSolve needs to work at better than machine precision to solve the equations. I was able to get the desired solution with the following:

GS65[wa_, wb_] = FullSimplify[...];
GS56[wa_, wb_] = FullSimplify[...];
NSolve[
  GS65[wa, wb] == 827531/1000000 && GS56[wa, wb] == 327531/1000000 && 
    0 < wa < 1 && 0 < wb < 1, 
  {wa, wb}, 
  WorkingPrecision -> 100]
{{wa -> 0.60000000000000000000000000000000000000000000000000000000000000000000000000000, 
  wb -> 0.50000000000000000000000000000000000000000000000000000000000000000000000000000}, 
 {wa -> 0.43646749272763910788150807887685400254767003351976708419711793418129855151965, 
  wb -> 0.43646749272763910788150807887685400254767003351976708419711793418129855151965}}

Less than 100 digits of working precision might work. I didn't have the time to experiment with lower values.

$\endgroup$
9
  • $\begingroup$ Well I use NSolve because I remember Solve didn't finish too, Do you think it would be possible to get the solution using Solve? Well now it was just great that you find a method to solve my problem but as I don't understand how you arrived at your conclusions the next time I face this problem (I have to solve a lot more polynomials that are even larger) maybe I won't know how to solve them. I mean, How did you get WorkingPrecision 100? $\endgroup$ – George Dec 29 '15 at 14:07
  • $\begingroup$ Another question is, Why did you enter the numbers this way? 827531/1000000 $\endgroup$ – George Dec 29 '15 at 14:10
  • $\begingroup$ @George. I didn't do anything with Solve because I regarded this as a numerics problem. I think Solve would have a lot of trouble with the 2nd set of your equations. I choose 100 for value of working precision by intuition. As I noted in my answer, if I had time, I would tried some smaller values. Numerics is an art, not a science. $\endgroup$ – m_goldberg Dec 29 '15 at 14:14
  • 1
    $\begingroup$ @George. I used rationals for the constants because high working precision needs the problem expressed with high-precision values. It is futile to attempt to compute at high precision when the inputs are low precision. $\endgroup$ – m_goldberg Dec 29 '15 at 14:20
  • $\begingroup$ Goldberg now it worked in my computer, it took some minutes but it worked. THank you very much for all your help, was great. :) $\endgroup$ – George Dec 29 '15 at 14:28
4
$\begingroup$

A similar approach using Groebner basis using DegreeReversedLexographic monomial ordering. It tends to yeild polynomials of lower degree.

gs65 = Numerator[Together[GS65 - Rationalize[0.827531, 0.0000001]]] //FullSimplify // Expand

gs56 = Numerator[Together[GS56 - Rationalize[0.327531, 0.000001]]] //FullSimplify // Expand

gb = GroebnerBasis[{gs65, gs56}, {wa, wb}, MonomialOrder -> DegreeReverseLexicographic]

sols = NSolve[gb == ConstantArray[0, Length[gb]] && 0 < wa < 1 && 0 < wb < 1, {wa,wb}, WorkingPrecision -> 20]

Out[]={{wa->0.59999827733375774661,wb->0.49999918333272064412},{wa->0.43646936814548225095,wb->0.43646808801540740031}}
$\endgroup$
3
  • $\begingroup$ I don't have the knowledge to understand your answer but I will try to understand it. It looks very helpful. Thank you very much. $\endgroup$ – George Dec 29 '15 at 14:44
  • $\begingroup$ Your functions were rational but to solve f(x,y)=a it can be turned into a root problem by solving f(x,y)-a=0, then only the numerator if this expression matters. The Groebner basis is a set of polynomials that have the same roots, (zero set) but in a simplified form, much like Gaussian elimination simplifies a set of linear equations. $\endgroup$ – John McGee Dec 29 '15 at 15:55
  • $\begingroup$ Interesting. Just a side note, one of the key point of this approach is the Method of GroebnerBasis can't be default i.e. Lexicographic. $\endgroup$ – xzczd Dec 30 '15 at 3:07
3
$\begingroup$

(1) One method for obtaining real solution is to solve numerically and discard explicitly complex ones. This might in some cases require using nondefault setting for WorkingPrecision.

(2) I believe this particular example also suffers from the use of a method in version 10 that might not be well suited, or sufficiently well tuned, to figure out zero imaginary parts. Hence I revert to the "old" (versions 4-9) default method below.

AbsoluteTiming[
 solns = NSolve[{GS65 - 0.827531, GS56 - 0.327531}, {wa, wb}, 
   Method -> "EndomorphismMatrix"];
 realsolns = Select[{wa, wb} /. solns, FreeQ[#, Complex] &];
 Sort[realsolns]]

(* Out[631]= {6.699682, {{-3.78071749624, 
   0.436467492728}, {-1.94039219538, 
   0.5}, {-0.423358009017, -0.421257113116}, {-0.384257820741, \
-0.384257820742}, {0.02922159727, -0.384257820739}, {0.0302914846213, \
-0.421257113116}, {0.358127822324, -0.384257820741}, {0.382260818155, \
-0.421257113116}, {0.436467492728, 0.436467492728}, {0.600000000001, 
   0.500000000001}, {0.631150406517, 1.4503714644}, {0.64686212397, 
   1.42125711309}, {0.980516028163, 1.45037146452}, {0.984049548683, 
   1.42125711351}, {1.41911218891, 1.42125711383}, {1.45037146479, 
   1.45037146503}, {2.25209053247, 0.436467492728}, {2.93740145723, 
   0.5}}} *)

The expected solution is certainly on the list.

Note that if one has a good idea of where the roots are located, then FindRoot becomes a useful tool for this. In this case it is tricky since there is another real solution not too far away.

FindRoot[{GS65 == 0.827531, GS56 == 0.327531}, {wa, .5}, {wb, .5}]

(* Out[634]= {wa -> 0.436467492728, wb -> 0.436467492728} *)

It seems the dynamics of whatever method FindRoot is using do not favor the desired solution in this example.

$\endgroup$
2
  • $\begingroup$ Indeed, in v9 one only needs to take the inequalities {0 < wa < 1, 0 < wb < 1} out of NSolve. (So the most troublesome part is actually the inequalities?) $\endgroup$ – xzczd Dec 30 '15 at 3:24
  • $\begingroup$ @xzczd Yes, inequalities are hard to deal with. They take one into the realm of real algebraic geometry, which, computationally, tends to be more taxing than the complex case (even though there are now fewer solutions). $\endgroup$ – Daniel Lichtblau Dec 30 '15 at 17:28
1
$\begingroup$

First, simplify you input:

GS65new = Simplify[GS65]

-(((-1 + wb) ((-1 + wb)^3 wb (1 + 2 wb + 4 wb^2 + 8 wb^3) + 15 wa^2 (-1 + wb)^4 wb^5 (-15 - 41 wb + 142 wb^2 - 132 wb^3 + 40 wb^4) - 5 wa^3 (-1 + wb)^3 wb^5 (165 - 389 wb - 708 wb^2 + 2264 wb^3 - 1952 wb^4 + 560 wb^5) + 5 wa^4 (-1 + wb)^2 wb^5 (-255 + 1613 wb - 2402 wb^2 - 984 wb^3 + 4992 wb^4 - 4144 wb^5 + 1120 wb^6) - 5 wa^7 wb^5 (-15 + 244 wb - 1344 wb^2 + 3564 wb^3 - 5166 wb^4 + 4228 wb^5 - 1848 wb^6 + 336 wb^7) - 3 wa^5 wb^5 (-345 + 3812 wb - 14607 wb^2 + 25702 wb^3 - 19648 wb^4 - 878 wb^5 + 12180 wb^6 - 7896 wb^7 + 1680 wb^8) + wa^6 wb^5 (-435 + 5951 wb - 28026 wb^2 + 63316 wb^3 - 75574 wb^4 + 45430 wb^5 - 8148 wb^6 - 4200 wb^7 + 1680 wb^8) - wa (-1 + wb)^3 (-1 + wb^4 - 11 wb^5 + 209 wb^6 - 457 wb^7 + 446 wb^8 - 212 wb^9 + 40 wb^10)))/((1 - 2 wb + 2 wb^2) (-wb + wa (-1 + 2 wb))))

GS56new = Simplify[GS56]

-((wa (-1 + wb)^5 (1 + 2 wb) (1 + 4 wb^2) ((-1 + wb)^5 (1 + 5 wb) - 15 wa (-1 + wb)^4 wb (1 + 5 wb) + 5 wa^2 (-1 + wb)^3 wb (-11 + wb + 70 wb^2) - 5 wa^3 (-1 + wb)^2 wb (17 - 69 wb - 28 wb^2 + 140 wb^3) + 5 wa^6 wb (1 - 14 wb + 56 wb^2 - 84 wb^3 + 42 wb^4) + 3 wa^4 wb (23 - 202 wb + 473 wb^2 - 252 wb^3 - 252 wb^4 + 210 wb^5) - wa^5 wb (29 - 331 wb + 1064 wb^2 - 1176 wb^3 + 210 wb^4 + 210 wb^5)))/((1 - 2 wb + 2 wb^2) (-wb + wa (-1 + 2 wb))))

Then calculate GroebnerBase

(gb = GroebnerBasis[{GS56new - 25/100, GS65new - 75/100}, {wa, 
     wb}]) // Length

4

Taking gb elements which includes both variables we get

NSolve[{gb[[1]] == 0, gb[[3]] == 0}, {wa, wb}]

{{wb -> -0.421257, wa -> -0.421257}, {wb -> -0.421257, wa -> 0.0231239}, {wb -> -0.421257, wa -> 0.365881}, {wb -> -0.421257, wa -> 0.703403 - 0.315787 I}, {wb -> -0.421257, wa -> 0.703403 + 0.315787 I}, {wb -> -0.421257, wa -> 1.15004 - 0.190458 I}, {wb -> -0.421257, wa -> 1.15004 + 0.190458 I}, {wb -> -0.421257, wa -> -0.421257}, {wb -> -0.421257, wa -> 0.0231239}, {wb -> -0.421257, wa -> 0.365881}, {wb -> -0.421257, wa -> 0.703403 - 0.315787 I}, {wb -> -0.421257, wa -> 0.703403 + 0.315787 I}, {wb -> -0.421257, wa -> 1.15004 - 0.190458 I}, {wb -> -0.421257, wa -> 1.15004 + 0.190458 I}, {wb -> -0.0023513 + 0.420374 I, wa -> -0.051519 - 0.00405429 I}, {wb -> -0.0023513 + 0.420374 I, wa -> -0.0023513 + 0.420374 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 0.0793953 + 0.468781 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 0.699276 - 0.192388 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 0.774672 + 0.677983 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 1.17814 - 0.187992 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 1.27259 + 0.290135 I}, {wb -> -0.0023513 + 0.420374 I, wa -> -0.051519 - 0.00405429 I}, {wb -> -0.0023513 + 0.420374 I, wa -> -0.0023513 + 0.420374 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 0.0793953 + 0.468781 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 0.699276 - 0.192388 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 0.774672 + 0.677983 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 1.17814 - 0.187992 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 1.27259 + 0.290135 I}, {wb -> -0.0023513 - 0.420374 I, wa -> -0.051519 + 0.00405429 I}, {wb -> -0.0023513 - 0.420374 I, wa -> -0.0023513 - 0.420374 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 0.0793953 - 0.468781 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 0.699276 + 0.192388 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 0.774672 - 0.677983 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 1.17814 + 0.187992 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 1.27259 - 0.290135 I}, {wb -> -0.0023513 - 0.420374 I, wa -> -0.051519 + 0.00405429 I}, {wb -> -0.0023513 - 0.420374 I, wa -> -0.0023513 - 0.420374 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 0.0793953 - 0.468781 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 0.699276 + 0.192388 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 0.774672 - 0.677983 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 1.17814 + 0.187992 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 1.27259 - 0.290135 I}, {wb -> 0.5, wa -> -1.9389}, {wb -> 0.5, wa -> -0.400713 - 0.332275 I}, {wb -> 0.5, wa -> -0.400713 + 0.332275 I}, {wb -> 0.5, wa -> 0.5}, {wb -> 0.5, wa -> 1.40071 - 0.332275 I}, {wb -> 0.5, wa -> 1.40071 + 0.332275 I}, {wb -> 0.5, wa -> 2.9389}, {wb -> 0.5, wa -> -1.9389}, {wb -> 0.5, wa -> -0.400713 - 0.332275 I}, {wb -> 0.5, wa -> -0.400713 + 0.332275 I}, {wb -> 0.5, wa -> 0.5}, {wb -> 0.5, wa -> 1.40071 - 0.332275 I}, {wb -> 0.5, wa -> 1.40071 + 0.332275 I}, {wb -> 0.5, wa -> 2.9389}, {wb -> 1.00235 - 0.420374 I, wa -> -0.272591 - 0.290135 I}, {wb -> 1.00235 - 0.420374 I, wa -> -0.178142 + 0.187992 I}, {wb -> 1.00235 - 0.420374 I, wa -> 0.225328 - 0.677983 I}, {wb -> 1.00235 - 0.420374 I, wa -> 0.300724 + 0.192388 I}, {wb -> 1.00235 - 0.420374 I, wa -> 0.920605 - 0.468781 I}, {wb -> 1.00235 - 0.420374 I, wa -> 1.00235 - 0.420374 I}, {wb -> 1.00235 - 0.420374 I, wa -> 1.05152 + 0.00405429 I}, {wb -> 1.00235 - 0.420374 I, wa -> -0.272591 - 0.290135 I}, {wb -> 1.00235 - 0.420374 I, wa -> -0.178142 + 0.187992 I}, {wb -> 1.00235 - 0.420374 I, wa -> 0.225328 - 0.677983 I}, {wb -> 1.00235 - 0.420374 I, wa -> 0.300724 + 0.192388 I}, {wb -> 1.00235 - 0.420374 I, wa -> 0.920605 - 0.468781 I}, {wb -> 1.00235 - 0.420374 I, wa -> 1.00235 - 0.420374 I}, {wb -> 1.00235 - 0.420374 I, wa -> 1.05152 + 0.00405429 I}, {wb -> 1.00235 + 0.420374 I, wa -> -0.272591 + 0.290135 I}, {wb -> 1.00235 + 0.420374 I, wa -> -0.178142 - 0.187992 I}, {wb -> 1.00235 + 0.420374 I, wa -> 0.225328 + 0.677983 I}, {wb -> 1.00235 + 0.420374 I, wa -> 0.300724 - 0.192388 I}, {wb -> 1.00235 + 0.420374 I, wa -> 0.920605 + 0.468781 I}, {wb -> 1.00235 + 0.420374 I, wa -> 1.00235 + 0.420374 I}, {wb -> 1.00235 + 0.420374 I, wa -> 1.05152 - 0.00405429 I}, {wb -> 1.00235 + 0.420374 I, wa -> -0.272591 + 0.290135 I}, {wb -> 1.00235 + 0.420374 I, wa -> -0.178142 - 0.187992 I}, {wb -> 1.00235 + 0.420374 I, wa -> 0.225328 + 0.677983 I}, {wb -> 1.00235 + 0.420374 I, wa -> 0.300724 - 0.192388 I}, {wb -> 1.00235 + 0.420374 I, wa -> 0.920605 + 0.468781 I}, {wb -> 1.00235 + 0.420374 I, wa -> 1.00235 + 0.420374 I}, {wb -> 1.00235 + 0.420374 I, wa -> 1.05152 - 0.00405429 I}, {wb -> 1.42126, wa -> -0.150039 - 0.190458 I}, {wb -> 1.42126, wa -> -0.150039 + 0.190458 I}, {wb -> 1.42126, wa -> 0.296597 - 0.315787 I}, {wb -> 1.42126, wa -> 0.296597 + 0.315787 I}, {wb -> 1.42126, wa -> 0.634119}, {wb -> 1.42126, wa -> 0.976876}, {wb -> 1.42126, wa -> 1.42126}, {wb -> 1.42126, wa -> -0.150039 - 0.190458 I}, {wb -> 1.42126, wa -> -0.150039 + 0.190458 I}, {wb -> 1.42126, wa -> 0.296597 - 0.315787 I}, {wb -> 1.42126, wa -> 0.296597 + 0.315787 I}, {wb -> 1.42126, wa -> 0.634119}, {wb -> 1.42126, wa -> 0.976876}, {wb -> 1.42126, wa -> 1.42126}}

which does not match your claim. Similarly, for other case we have:

(gb1 = GroebnerBasis[{GS56new - Rationalize[0.327531, 0], 
     GS65new - Rationalize[0.827531, 0]}, {wa, wb}]) // Length

4

NSolve[{gb1[[1]] == 0, gb[[3]] == 0}, {wa, wb}]

{{wb -> -0.421257, wa -> -0.421257}, {wb -> -0.421257, wa -> 0.0231239}, {wb -> -0.421257, wa -> 0.365881}, {wb -> -0.421257, wa -> 0.703403 - 0.315787 I}, {wb -> -0.421257, wa -> 0.703403 + 0.315787 I}, {wb -> -0.421257, wa -> 1.15004 - 0.190458 I}, {wb -> -0.421257, wa -> 1.15004 + 0.190458 I}, {wb -> -0.384258, wa -> -2.00696}, {wb -> -0.384258, wa -> -0.814158}, {wb -> -0.384258, wa -> -0.0321924 - 0.434251 I}, {wb -> -0.384258, wa -> -0.0321924 + 0.434251 I}, {wb -> -0.384258, wa -> 1.07842 - 0.461507 I}, {wb -> -0.384258, wa -> 1.07842 + 0.461507 I}, {wb -> -0.384258, wa -> 1.62828}, {wb -> -0.00294383 - 0.38345 I, wa -> -0.356037 - 0.0105394 I}, {wb -> -0.00294383 - 0.38345 I, wa -> -0.0682849 - 2.03837 I}, {wb -> -0.00294383 - 0.38345 I, wa -> 0.0271765 - 0.415577 I}, {wb -> -0.00294383 - 0.38345 I, wa -> 0.144308 + 0.378458 I}, {wb -> -0.00294383 - 0.38345 I, wa -> 1.03701 - 0.484021 I}, {wb -> -0.00294383 - 0.38345 I, wa -> 1.0963 + 0.378683 I}, {wb -> -0.00294383 - 0.38345 I, wa -> 1.54464 - 0.0621063 I}, {wb -> -0.00294383 + 0.38345 I, wa -> -0.356037 + 0.0105394 I}, {wb -> -0.00294383 + 0.38345 I, wa -> -0.0682849 + 2.03837 I}, {wb -> -0.00294383 + 0.38345 I, wa -> 0.0271765 + 0.415577 I}, {wb -> -0.00294383 + 0.38345 I, wa -> 0.144308 - 0.378458 I}, {wb -> -0.00294383 + 0.38345 I, wa -> 1.03701 + 0.484021 I}, {wb -> -0.00294383 + 0.38345 I, wa -> 1.0963 - 0.378683 I}, {wb -> -0.00294383 + 0.38345 I, wa -> 1.54464 + 0.0621063 I}, {wb -> -0.0023513 - 0.420374 I, wa -> -0.051519 + 0.00405429 I}, {wb -> -0.0023513 - 0.420374 I, wa -> -0.0023513 - 0.420374 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 0.0793953 - 0.468781 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 0.699276 + 0.192388 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 0.774672 - 0.677983 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 1.17814 + 0.187992 I}, {wb -> -0.0023513 - 0.420374 I, wa -> 1.27259 - 0.290135 I}, {wb -> -0.0023513 + 0.420374 I, wa -> -0.051519 - 0.00405429 I}, {wb -> -0.0023513 + 0.420374 I, wa -> -0.0023513 + 0.420374 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 0.0793953 + 0.468781 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 0.699276 - 0.192388 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 0.774672 + 0.677983 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 1.17814 - 0.187992 I}, {wb -> -0.0023513 + 0.420374 I, wa -> 1.27259 + 0.290135 I}, {wb -> 0.436467, wa -> -3.19878}, {wb -> 0.436467, wa -> -0.376344 - 0.188877 I}, {wb -> 0.436467, wa -> -0.376344 + 0.188877 I}, {wb -> 0.436467, wa -> 0.334027}, {wb -> 0.436467, wa -> 1.36623 - 0.459684 I}, {wb -> 0.436467, wa -> 1.36623 + 0.459684 I}, {wb -> 0.436467, wa -> 2.21998}, {wb -> 0.5, wa -> -1.9389}, {wb -> 0.5, wa -> -0.400713 - 0.332275 I}, {wb -> 0.5, wa -> -0.400713 + 0.332275 I}, {wb -> 0.5, wa -> 0.5}, {wb -> 0.5, wa -> 1.40071 - 0.332275 I}, {wb -> 0.5, wa -> 1.40071 + 0.332275 I}, {wb -> 0.5, wa -> 2.9389}, {wb -> 1.00165 - 0.449612 I, wa -> -0.528678 - 0.0307153 I}, {wb -> 1.00165 - 0.449612 I, wa -> -0.064596 + 0.384441 I}, {wb -> 1.00165 - 0.449612 I, wa -> -0.0492989 - 0.491668 I}, {wb -> 1.00165 - 0.449612 I, wa -> 0.883162 - 2.14683 I}, {wb -> 1.00165 - 0.449612 I, wa -> 0.891719 + 0.347417 I}, {wb -> 1.00165 - 0.449612 I, wa -> 0.986124 - 0.425576 I}, {wb -> 1.00165 - 0.449612 I, wa -> 1.37354 - 0.0410766 I}, {wb -> 1.00165 + 0.449612 I, wa -> -0.528678 + 0.0307153 I}, {wb -> 1.00165 + 0.449612 I, wa -> -0.064596 - 0.384441 I}, {wb -> 1.00165 + 0.449612 I, wa -> -0.0492989 + 0.491668 I}, {wb -> 1.00165 + 0.449612 I, wa -> 0.883162 + 2.14683 I}, {wb -> 1.00165 + 0.449612 I, wa -> 0.891719 - 0.347417 I}, {wb -> 1.00165 + 0.449612 I, wa -> 0.986124 + 0.425576 I}, {wb -> 1.00165 + 0.449612 I, wa -> 1.37354 + 0.0410766 I}, {wb -> 1.00235 - 0.420374 I, wa -> -0.272591 - 0.290135 I}, {wb -> 1.00235 - 0.420374 I, wa -> -0.178142 + 0.187992 I}, {wb -> 1.00235 - 0.420374 I, wa -> 0.225328 - 0.677983 I}, {wb -> 1.00235 - 0.420374 I, wa -> 0.300724 + 0.192388 I}, {wb -> 1.00235 - 0.420374 I, wa -> 0.920605 - 0.468781 I}, {wb -> 1.00235 - 0.420374 I, wa -> 1.00235 - 0.420374 I}, {wb -> 1.00235 - 0.420374 I, wa -> 1.05152 + 0.00405429 I}, {wb -> 1.00235 + 0.420374 I, wa -> -0.272591 + 0.290135 I}, {wb -> 1.00235 + 0.420374 I, wa -> -0.178142 - 0.187992 I}, {wb -> 1.00235 + 0.420374 I, wa -> 0.225328 + 0.677983 I}, {wb -> 1.00235 + 0.420374 I, wa -> 0.300724 - 0.192388 I}, {wb -> 1.00235 + 0.420374 I, wa -> 0.920605 + 0.468781 I}, {wb -> 1.00235 + 0.420374 I, wa -> 1.00235 + 0.420374 I}, {wb -> 1.00235 + 0.420374 I, wa -> 1.05152 - 0.00405429 I}, {wb -> 1.42126, wa -> -0.150039 - 0.190458 I}, {wb -> 1.42126, wa -> -0.150039 + 0.190458 I}, {wb -> 1.42126, wa -> 0.296597 - 0.315787 I}, {wb -> 1.42126, wa -> 0.296597 + 0.315787 I}, {wb -> 1.42126, wa -> 0.634119}, {wb -> 1.42126, wa -> 0.976876}, {wb -> 1.42126, wa -> 1.42126}, {wb -> 1.45037, wa -> -0.630092}, {wb -> 1.45037, wa -> -0.0748663 - 0.462553 I}, {wb -> 1.45037, wa -> -0.0748663 + 0.462553 I}, {wb -> 1.45037, wa -> 1.02769 - 0.432315 I}, {wb -> 1.45037, wa -> 1.02769 + 0.432315 I}, {wb -> 1.45037, wa -> 1.71992}, {wb -> 1.45037, wa -> 3.70918}}

Taking other Groebner base elements your will get more/less roots.

Edit1 As it was correctly noted in the comment below, I try to select all solutions which can be given probability interpretation.

Note that in the code below I removed (command Most[]) one pair of GroebnerBase, which actually cases the problem for NSolve. (The coefficients of Groebner base seems are very large, may be somebody will comment on that)

positiveSols = 
 Select[Flatten[
   NSolve[#, {wa, wb}] & /@ 
    Most[(Part[gb1, #] & /@ Subsets[Range[4], {2}])], 1], 
  MatchQ[({wa, wb} /. #), {_?Positive, _?Positive}] &]

{{wb -> 0.436467, wa -> 0.436467}, {wb -> 0.436467, wa -> 2.25209}, {wb -> 0.5, wa -> 0.6}, {wb -> 0.5, wa -> 2.9374}, {wb -> 1.42126, wa -> 0.646862}, {wb -> 1.42126, wa -> 0.98405}, {wb -> 1.42126, wa -> 1.41911}, {wb -> 1.45037, wa -> 0.63115}, {wb -> 1.45037, wa -> 0.980516}, {wb -> 1.45037, wa -> 1.45037}, {wb -> 0.5, wa -> 0.194035}, {wb -> 0.5, wa -> 2.26841}, {wb -> 1.45037, wa -> 0.63115}, {wb -> 1.45037, wa -> 0.980516}, {wb -> 1.45037, wa -> 1.45037}, {wb -> 1.42126, wa -> 0.646862}, {wb -> 1.42126, wa -> 0.98405}, {wb -> 1.42126, wa -> 1.41911}, {wb -> 0.5, wa -> 0.6}, {wb -> 0.5, wa -> 2.9374}, {wb -> 0.436467, wa -> 0.436467}, {wb -> 0.436467, wa -> 2.25209}}

And at last

probabilitySols = 
 Select[positiveSols, 
   MatchQ[({wa, wb} /. #), {_?(0 <= # <= 1 &), _?(0 <= # <= 1 &)}] &]

{{wb -> 0.436467, wa -> 0.436467}, {wb -> 0.436467, wa -> 0.436467}, {wb -> 0.5, wa -> 0.6}, {wb -> 0.5, wa -> 0.194035}, {wb -> 0.5, wa -> 0.6}}

we see that mentioned solution {wb -> 0.5, wa -> 0.6} indeed is in the list. Once again, this list of solutions is still not complete, because one pair of Groebner base equations was removed.

$\endgroup$
2
  • 2
    $\begingroup$ The question states that wa and wb are probabilities, therefore real and in the range 0-1. It might be helpful to eliminate the non-compliant solutions in your answer. $\endgroup$ – Simon Woods Dec 29 '15 at 12:32
  • $\begingroup$ Agreed. One needs to get all possible solutions of total 4 possible Groebner base equations and select all satisfying these conditions. $\endgroup$ – Acus Dec 29 '15 at 12:56

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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