Re: Cannot NSolve a system of equations
- To: mathgroup at smc.vnet.net
- Subject: [mg88944] Re: Cannot NSolve a system of equations
- From: murat.koyuncu at gmail.com
- Date: Thu, 22 May 2008 02:34:40 -0400 (EDT)
- References: <g0uahi$4tk$1@smc.vnet.net> <g11rbc$agd$1@smc.vnet.net>
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
suddenly disappear or jump to an improbable value.
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!
- Follow-Ups:
- Re: Re: Cannot NSolve a system of equations
- From: Andrzej Kozlowski <akoz@mimuw.edu.pl>
- Re: Re: Cannot NSolve a system of equations