Re: Re: Cannot NSolve a system of equations

• To: mathgroup at smc.vnet.net
• Subject: [mg88995] Re: [mg88944] Re: Cannot NSolve a system of equations
• From: Andrzej Kozlowski <akoz at mimuw.edu.pl>
• Date: Fri, 23 May 2008 03:07:32 -0400 (EDT)
• References: <g0uahi\$4tk\$1@smc.vnet.net> <g11rbc\$agd\$1@smc.vnet.net> <200805220634.CAA22304@smc.vnet.net>

Note that:

phi = Rationalize[0.74615385]

so rationalizing phi in this way does nothing at all. Moreover, the
fact that your zet is a machine precision number means that the
computation is done with machine precision numbers.
There is, of course, no unique way to rationalize approximate numbers.
When I use Rationalize[ ,0] to force everything in your equations to
be rational, I do not get any solutions which satisfy your condition
that all the values should lie between 0 and 1. But what is more
strange: when simply copy and past your code into a new notebook and
evaluate it with fresh kernel I do not get your solution. This may be
platform dependent (I am using 6.02 on a MacBook).

I think the fact that the computation with machine precision is
numerically unstable is not surprising. What is surprising is that the
answer that you obtained appears quite accurate - even though it was
done essentially with machine precision numbers. That may not be
entirely true, since some of the equations involve rational numbers so
some of the intermediate computations may well have bee done with high

Still, the fact that, even running your code without any changes I get
a different answer is strange, to say the least. I hope this is
verified by other people using Macs.

Andrzej Kozlowski

PS. When I run your code the only answer that I get that satisfies
your conditions that everything lies between 0 and 1 is

{x -> 1., x1 -> 1., x2 -> 1., x3 -> 1., x4 -> 1., x5 -> 1., abar ->
0.160682,
roverw -> 0.}

which is quite far form yours.

As I wrote earlier, if the equations are genuinely rationalized, none
of the solutions satisfies your conditions.

Andrzej Kozlowski

On 22 May 2008, at 15:34, murat.koyuncu at gmail.com wrote:

> I sincerely thank you all for your help. It cleared up a lot of things
> for me. I had tried Rationalize but didn't know about precision
> issues.
>
> Now I have a follow-up question. I still cannot trust my results 100%
> because when I slightly change one of my parameters, solution set
> changes drastically. My theory says that resulting variables should
> between 0 and 1, so I am only interested in that range. And sometimes
> I don't even get any such solutions.
>
> For example, here is my model (very much polished, thanks to your
> suggestions):
>
> Unprotect[In, Out]; Clear[In, Out]; ClearAll["Global`*"];
> phi = Rationalize[0.74615385];
> eta = 7/4;
> alpha = 16/25;
> taubar = 13/100;
> {y1, y2, y3, y4, y5} =
>  Rationalize[{0.235457064, 0.512465374, 0.781779009, 1.109572176,
>    2.360726377}, 0];
> zet = N[5 taubar/(y1^(1 + phi) + y2^(1 + phi) + y3^(1 + phi) + y4^(
>       1 + phi) + y5^(1 + phi))];
> tau1 = zet y1^phi;
> tau2 = zet y2^phi;
> tau3 = zet y3^phi;
> tau4 = zet y4^phi;
> tau5 = zet y5^phi;
> a1 = (1 + phi) tau1;
> a2 = (1 + phi) tau2;
> a3 = (1 + phi) tau3;
> a4 = (1 + phi) tau4;
> a5 = (1 + phi) tau5;
> eqns1 = {x1 == (roverw (1 - tau1) + ((roverw + (1 - x)) y1 -
>          1) ((1 - taubar + (1 - abar)/eta) x + (taubar - tau1)
>           roverw - (1 - taubar)))/(roverw (1 -
>          tau1 + (1 - a1)/eta) - ((1 - taubar + (1 - abar)/eta)
>          x + (taubar - tau1) roverw - (1 - taubar))),
>   x2 == (roverw (1 - tau2) + ((roverw + (1 - x)) y2 -
>          1) ((1 - taubar + (1 - abar)/eta) x + (taubar - tau2)
>           roverw - (1 - taubar)))/(roverw (1 -
>          tau2 + (1 - a2)/eta) - ((1 - taubar + (1 - abar)/eta)
>          x + (taubar - tau2) roverw - (1 - taubar))),
>   x3 == (roverw (1 - tau3) + ((roverw + (1 - x)) y3 -
>          1) ((1 - taubar + (1 - abar)/eta) x + (taubar - tau3)
>           roverw - (1 - taubar)))/(roverw (1 -
>          tau3 + (1 - a3)/eta) - ((1 - taubar + (1 - abar)/eta)
>          x + (taubar - tau3) roverw - (1 - taubar))),
>   x4 == (roverw (1 - tau4) + ((roverw + (1 - x)) y4 -
>          1) ((1 - taubar + (1 - abar)/eta) x + (taubar - tau4)
>           roverw - (1 - taubar)))/(roverw (1 -
>          tau4 + (1 - a4)/eta) - ((1 - taubar + (1 - abar)/eta)
>          x + (taubar - tau4) roverw - (1 - taubar))),
>   x5 == (roverw (1 - tau5) + ((roverw + (1 - x)) y5 -
>          1) ((1 - taubar + (1 - abar)/eta) x + (taubar - tau5)
>           roverw - (1 - taubar)))/(roverw (1 -
>          tau5 + (1 - a5)/eta) - ((1 - taubar + (1 - abar)/eta)
>          x + (taubar - tau5) roverw - (1 - taubar))),
>   x == (x1 + x2 + x3 + x4 + x5)/5,
>   abar == (a1 x1 + a2 x2 + a3 x3 + a4 x4 + a5 x5)/(5 x),
>   roverw == (1 - x) (1 - alpha)/alpha};
> vars = {x, x1, x2, x3, x4, x5, abar, roverw};
> sol = NSolve[N[eqns1, 500], vars];
> InputForm[N[sol]]
>
> One of the solutions I get is between 0 and 1:
> {x -> 0.5004885149810059, x1 -> 0.6893393376625371,
>  x2 -> 0.6337652704581314, x3 -> 0.5731109245369462,
>  x4 -> 0.4917632705087048, x5 -> 0.1144637717387098,
>  abar -> 0.1244737777525507, roverw -> 0.2809752103231842}
>
> So I am happy. But when I drop the last digit of the first parameter,
> i.e. phi = Rationalize[0.7461538];, the new set of solutions do not
> have any results compatible with my assumptions:
>
> {{x -> 0.2949617829869077, x1 -> -1.9491786143565832, x2 ->
> 10.96762641745164, x3 -> 96.0369893781559, x4 -> -73.18097932930044,
>  x5 -> -30.399648937016046, abar -> -5.824073507107748, roverw ->
> 0.39658399706986447},
> {x -> 2.0689803436369942, x1 -> -8.65276583855146, x2 ->
> -18.39904428993495, x3 -> -33.55295774059405, x4 ->
> -69.9736576669751,
>  x5 -> 140.92332725424055, abar -> 2.4458078645586667, roverw ->
> -0.6013014432958093},
> {x -> 0.6213605439915423, x1 -> -12.823085779076337, x2 ->
> 41.57610980587073, x3 -> -6.337885849417286, x4 ->
> -8.665252197121296,
>  x5 -> -10.64308326029811, abar -> -0.7583206767716181, roverw ->
> 0.21298469400475745},
> {x -> 0.6854511826071373 - 0.04362902433398709*I, x1 ->
> 0.59964456896311 + 0.08549874839569205*I,
>  x2 -> 0.6411455792212575 + 0.04079939260558141*I, x3 ->
> 0.6760747168284579 - 0.003707834204090894*I,
>  x4 -> 0.7125762414563436 - 0.05913956075670094*I, x5 ->
> 0.7978148065665167 - 0.28159586771041706*I,
>  abar -> 0.17051392124951886 - 0.016129276403199968*I, roverw ->
> 0.17693370978348527 + 0.024541326187867737*I},
> {x -> 0.6854511826071373 + 0.04362902433398709*I, x1 ->
> 0.59964456896311 - 0.08549874839569205*I,
>  x2 -> 0.6411455792212575 - 0.04079939260558141*I, x3 ->
> 0.6760747168284579 + 0.003707834204090894*I,
>  x4 -> 0.7125762414563436 + 0.05913956075670094*I, x5 ->
> 0.7978148065665167 + 0.28159586771041706*I,
>  abar -> 0.17051392124951886 + 0.016129276403199968*I, roverw ->
> 0.17693370978348527 - 0.024541326187867737*I},
> {x -> 0.9999999999989753, x1 -> 0.9999999999982316, x2 ->
> 0.9999999999985612, x3 -> 0.999999999998809, x4 ->
> 0.9999999999991211,
>  x5 -> 1.0000000000001525, abar -> 0.16068192162412143, roverw ->
> 0.}}
>
> And this is just one of the instances that this happens. I solve the
> same set of equations for many parameter sets. Sometimes I get a
> solution in the [0,1] range, but it is so different than the previous
> one, it doesn't make sense. Or sometimes, I have to tweak parameter
> values (drop a digit here, add another one there) to get a meaningful
> solution. I need a set of results over different parameter values that
> I can compare and contrast, but it is hard to do that when my results
>
> Now I know that this is probably the most Mathematica can do for me,
> but is there a way to make my system more "stable"? Maybe play with my
> equations a little bit? Or put another way, what makes my system so
> unstable?
>
> Again, any suggestions will be very much appreciated. Thanks!
>

• Prev by Date: Re: Re: On which OS is Mathematica best implemented?
• Next by Date: Re: Interpolation and plot doing strange things with mathematica
• Previous by thread: Re: Cannot NSolve a system of equations
• Next by thread: Re: Re: Cannot NSolve a system of equations