Re: Re: Re: Re: How to solve nonlinear equations?

```DrBob wrote:
> Neither of us has answered the original question, which was "How to solve...".
>
> Let's say we need 200 digits. How do you solve it?
>
> Bobby
> [...]
>>"DrBob" <drbob at bigfoot.com> schrieb im Newsbeitrag
>>news:comgo1\$89e\$1 at smc.vnet.net...
>>
>>>The only trouble with THAT is that it doesn't work.
>>>
>>>eqns = {20.428836553499018 - Log[X1] == 468/67*X1 + 5790/1273*X2 -
>>>66257/1273*
>>>        X3 - 24660/1273*X4 - 79150/1273*X5, 17.011108423692498 -
>>>           Log[X2] == 5790/1273*X1 + 6294/1273*X2 - 66257/1273*X3 - 24660/
>>>          1273*X4 - 403.8069285*
>>>            X5, -29.72639373695347 - Log[X3] == -66257/1273*
>>>          X1 - 66257/1273*X2 - 2*
>>>        X3, -26.273726271581616 - Log[X4] == -24660/1273*X1 - 24660/
>>>          1273*X2 + 10.15330715*X4, -38.76695085346396 - Log[X5]
>>>== -79150/
>>>          1273*X1 - 403.8069285*X2 - 10.67374705*X5};
>>>
>>>FindRoot[eqns, {{X1, 1}, {X2, 1}, {X3, 1}, {X4, 1}, {X5, 1}}]
>>>eqns /. Equal -> Subtract /. %
>>>
>>>(FindRoot::lstol error)
>>>{X1 -> 0.5287215008048355 + 0.0078093799604902845*I,
>>>  X2 -> 0.002514442570050682 - 0.007506800892680369*I,
>>>  X3 -> 0.2750175495952261 + 0.0051666476248725625*I,
>>>  X4 -> -1.6157105189511614 -4.778719413301337*^-15*I,
>>>  X5 -> -0.006058532522555838 - 0.003759309847863656*I}
>>>{-0.0000172993518866571 - 5.25106242902848*^-7*I,
>>>  0.0014397436481736747 + 0.000060457451791506855*I,
>>>  -0.23588972784670026 + 0.007297504899863284*I,
>>>  -0.057825246363192695 + 3.1474540831359445*I,
>>>  0.0009303265262765592 + 0.00036856201666024546*I}
>>>
>>>Bobby
>>> [...]
>>>>
>>>>"Wei Wang" <weiwang at baosteel.com> schrieb im Newsbeitrag
>>>>news:cohj9d\$1nr\$1 at smc.vnet.net...
>>>>
>>>>>How to solve the following equations, where X1, X2, X3, X4 and X5 are
>>>>>variables?
>>>>>
>>>>>eqns = {20.428836553499018-ln(X1) ==
>>>>>468/67*X1+5790/1273*X2-66257/1273*X3-24660/1273*X4-79150/1273*X5,
>>>>>17.011108423692498-ln(X2) ==
>>>>>5790/1273*X1+6294/1273*X2-66257/1273*X3-24660/1273*X4-403.8069285*X5, -29.72639373695347-ln(X3)
>>>>>== -66257/1273*X1-66257/1273*X2-2*X3, -26.273726271581616-ln(X4)
>>>>>== -24660/1273*X1-24660/1273*X2+10.15330715*X4, -38.76695085346396-ln(X5)
>>>>>== -79150/1273*X1-403.8069285*X2-10.67374705*X5};
>>>>>
>>>
>>>--
>>>DrBob at bigfoot.com
>>>www.eclecticdreams.net

I'm not sure about the case of complex roots, but there is no real
solution to this system. Showing that took me several steps though there
may be a more enlightened method than what I reproduce below.

eqns = {20.428836553499018 - Log[X[1]] -
(468/67*X[1] + 5790/1273*X[2] - 66257/1273*X[3] -
24660/1273*X[4] - 79150/1273*X[5]),
17.011108423692498 - Log[X[2]] -
(5790/1273*X[1] + 6294/1273*X[2] - 66257/1273*X[3] -
24660/1273*X[4] - 403.8069285X[5]),
-29.72639373695347 - Log[X[3]] -
(-66257/1273*X[1] - 66257/1273*X[2] - 2*X[3]),
-26.273726271581616 - Log[X[4]] -
(-24660/1273*X[1] - 24660/1273*X[2] + 10.15330715*X[4]),
-38.76695085346396 - Log[X[5]] -
(-79150/1273*X[1] - 403.8069285*X[2] - 10.67374705*X[5])};

exprs = eqns /. Equal->Subtract;

First notice that we can replace each X[j] by its exponential. Moreover
as all right hand sides are free of the logs, we know the X[j]'s must be
positive. Hence we can use PowerExpand to in effect make those logs our
new variables. I will also Rationalize though that's not strictly necessary.

newexprs = Rationalize[PowerExpand[exprs /. X[j_]:>Exp[lX[j]]], 0]

I now use some elimination ideals to obtain various expressions that
must vanish. From these certain inequality deductions may be made.

expr34 = GroebnerBasis[newexprs, Exp[Array[lX,5]],
Complement[Array[lX,5],{lX[3],lX[4]}]];

We first investigate what happens if lX[4] is nonnegative.

In[280]:= InputForm[e3 = First[expr34] /. lX[4]->Interval[{0,Infinity}]]
Out[280]//InputForm=
647238762444175200000000*E^lX[3] +
Interval[{22053533331180859235803693, Infinity}] -
323619381222087600000000*lX[3]

The next steps perhaps obvious to many people. We first see that this is
strictly positive if lX[3] is negative.

In[289]:= e3 /. lX[3]->Interval[{-Infinity,0}]
Out[289]= Interval[{22053533331180859235803693, Infinity}]

We now check that the derivative with respect to lX[3] is strictly
positive at for positive lX[3]. This guarantees that e3 cannot get
smaller, hence is positive for all real lX[3].

In[290]:= D[e3,lX[3]] /. lX[3]->Interval[{0,Infinity}]
Out[290]= Interval[{323619381222087600000000, Infinity}]

We conclude at this point that we must require lX[4] to be negative.
Taking this into account we will show that lX[1]<1 by showing that if we
assume otherwise then a certain expression in our list will be strictly
positive.

In[291]:= N[newexprs[[-2]] /. {lX[1]->Interval[{1,Infinity}],
lX[4]->Interval[{-Infinity,0}],
lX[2]->Interval[{-Infinity,Infinity}]}]
Out[291]= Interval[{16.2303, Infinity}]

Similarly we see we must have lX[2]<1.

In[293]:= N[newexprs[[-2]] /. {lX[2]->Interval[{1,Infinity}],
lX[4]->Interval[{-Infinity,0}],
lX[1]->Interval[{-Infinity,Infinity}]}]
Out[293]= Interval[{16.2303, Infinity}]

We now show these each must in fact be negative.

expr5 = First[GroebnerBasis[newexprs, Exp[Array[lX,5]],
Complement[Array[lX,5],{lX[5]}]]];

e15 = expr5 /. {lX[1]->Interval[{0,1}], lX[2]->Interval[{-Infinity,1}]};

First show this is strictly positive for all negative lX[5].

In[297]:= e15 /. lX[5]->Interval[{-Infinity,0}]
Out[297]= Interval[{1022379647864000000, Infinity}]

Next show the derivative with respect to lX[5] is strictly positive for
lX[5]>0.

In[299]:= D[e15,lX[5]] /. lX[5]->Interval[{0,Infinity}]
Out[299]= Interval[{422497213422049953, Infinity}]

This shows e15 cannot vanish for any real lX[5]. We can similarly show that

e25 = {lX[2]->Interval[{0,1}], lX[1]->Interval[{-Infinity,1}]};

cannot vanish for any lX[5]. We now know that {lX[1],lX[2],lX[4]} are
all negative. We move to another elimination ideal.

expr3 = First[GroebnerBasis[newexprs, Exp[Array[lX,5]],
Complement[Array[lX,5],{lX[3]}]]]

We see that this is strictly positive if lX[3] is 3.

In[311]:= (expr3 /. {lX[3]->3,
lX[1]->Interval[{-Infinity,0}], lX[2]->Interval[{-Infinity,0}]}) > 0
Out[311]= True

The derivative with respect to lX[3] is always positive for lX[3]>3.

In[313]:= (D[expr3,lX[3]] /. lX[3]->Interval[{3,Infinity}]) > 0
Out[313]= True

Thus expr3 cannot vanish for lX[3]>=3, so lX[3]<3. We next look at one
of our original expressions when lX[3]<-30. We will see it must in this
case be strictly positive.

In[316]:= (newexprs[[3]] /. {lX[3]->Interval[{-Infinity,-30}],
lX[2]->Interval[{-Infinity,0}], lX[1]->Interval[{-Infinity,0}]}) > 0
Out[316]= True

We now know -30<lX[3]<3. Finally we use all this in the first of newexprs.

In[318]:= (newexprs[[1]] /. {lX[3]->Interval[{-30,3}],
lX[2]->Interval[{-Infinity,0}],
lX[4]->Interval[{-Infinity,0}], lX[5]->Interval[{-Infinity,Infinity}],
lX[1]->Interval[{-Infinity,0}]}) > 0
Out[318]= True

As it is strictly positive there are no real values that make newexprs
simultaneously vanish.

Daniel Lichtblau
Wolfram Research

```

• Prev by Date: Re: Re: A problem of numerical precision
• Next by Date: How to evaluate Exp[I Pi(1+x)]?
• Previous by thread: Re: Re: Re: How to solve nonlinear equations?
• Next by thread: GUIKit and VisibleRect