I am trying to find the root of a certain expression in Mathematica
version 7:
expr = 110.52499999999998 + (300. - 135.52499999999998/(1 - x)) (1 -
x) - 300. x - 135.52499999999998 Log[1 - x] + 135.52499999999998
Log[x]
It appears to plot fine, for example using Plot[expr, {x, 0, 1}]. The
plot shows that there should be a root at about x=0.85. However, when
I try to find this root, using for example the following:
FindRoot[expr, {x, 0.5}]
I get an error message:
"FindRoot::lstol: The line search decreased the step size to within
tolerance specified by AccuracyGoal and PrecisionGoal but was unable
to find a sufficient decrease in the merit function. You may need
more than MachinePrecision digits of working precision to meet these
tolerances."
and it prints a seemingly incorrect (according to the qualitative form
of the plot) result: {x -> 0.344678}. Only if I use for example
FindRoot[expr, {x, 0.7}]
do I get the seemingly "correct" root: {x -> 0.849823}.
Can you help me see why the FindRoot is getting stuck at {x ->
0.344678} when I use starting values far away from 0.7 or 0.8? I will
ultimately want to find the roots of many similar functions, which may
have more than one "actual" root, so it would be helpful if I could
see why FindRoot[expr, {x, 0.5}] does not give {x -> 0.849823}. (also
when I tried NSolve[expr==0,x], Mathematica will not solve it.)
Thank you,
Andrew DeYoung
Carnegie Mellon University
please try:
In[18]:= expr =
110.52499999999998 + (300. - 135.52499999999998/(1 - x))*(1 - x) -
300.*x -
135.52499999999998*Log[1 - x] + 135.52499999999998*Log[x]
FindRoot[expr == 0, {x, 0.1, 0.9}]
Out[18]= 110.52499999999998 + (300. -
135.52499999999998/(1 - x))*(1 - x) - 300.*x -
135.52499999999998*Log[1 - x] + 135.52499999999998*Log[x]
Out[19]= {x -> 0.849823089232009}
Saludos,
MATTHIAS BODE
S 17.35775=B0, W 066.14577=B0
2'740 m
AMSL.
> Date: Tue, 18 Jan 2011 05:51:50 -0500
> From: adey...@andrew.cmu.edu
> Subject: Using FindRoot on an equation involving Log terms
> To: math...@smc.vnet.net
The problem is likely the direction it uses to begin the search.
x = 0.655 is a minimum. If you give a initial value above that it finds your
root.
If you give a initial value below that it will go the other way which
eventually diverges when x approaches zero.
You can take a better look using the following two lines of code:
<< Optimization`UnconstrainedProblems`
FindRootPlot[expr, {x, 0.6}, PlotRange -> {{0, 1}, {-100, 100}}]
(The output graph comes out bugged but you can still see what you need to
see: the steps it took).
The colors are explained in
http://reference.wolfram.com/mathematica/tutorial/UnconstrainedOptimizationPlottingSearchData.html
By the way, FindInstance works:
FindInstance[expr == 0, x, Reals] // N
Cheers,
Gabriel
The plot also shows that 0.5 is not in the basin of attraction for
Newton's method to find the root, and that it will send one in the
direction of the local max which is around x=0.35
Other methods might work better. For example, if you provide two start
points that bracket a root, then Brent's method will be used and it will
successfully find this root.
In[84]:= FindRoot[expr, {x, 0.5, .9}]
Out[84]= {x -> 0.849823}
> FindRoot[expr, {x, 0.7}]
>
> do I get the seemingly "correct" root: {x -> 0.849823}.
The zone of attraction for the root is to the right of the local min,
which appears to be somewhere around x=0.65.
> Can you help me see why the FindRoot is getting stuck at {x ->
> 0.344678} when I use starting values far away from 0.7 or 0.8? I will
> ultimately want to find the roots of many similar functions, which may
> have more than one "actual" root, so it would be helpful if I could
> see why FindRoot[expr, {x, 0.5}] does not give {x -> 0.849823}. (also
> when I tried NSolve[expr==0,x], Mathematica will not solve it.)
>
> Thank you,
>
> Andrew DeYoung
> Carnegie Mellon University
Solve or NSolve with restrictions on the domain will handle this (in
version 8).
In[82]:= Solve[{expr == 0, 0 <= x <= 1}, x]
Out[82]= {{x -> 0.849823}}
Daniel Lichtblau
Wolfram Research
It needn't be this verbose, but here's an approach that works well when
the plot is well-behaved:
Clear[expr]
expr[x_] =
110.52499999999998 + (300. - 135.52499999999998/(1 - x)) (1 - x) -
300. x - 135.52499999999998 Log[1 - x] +
135.52499999999998 Log[x];
p0 = Plot[expr@x, {x, 0, 1}];
pt1 = First@Cases[p0, Line[pts_List] :> pts, Infinity];
p1 = ListPlot[pt1, PlotRange -> All];
pt2 = pt1 /. {x_, y_} -> {x, Abs@y};
p2 = ListPlot@pt2;
sorted = SortBy[pt2, Last];
x0 = sorted[[1, 1]]
root = x /. FindRoot[expr@x, {x, x0}]
0.849313
0.849823
expr /@ {x0, root}
{-0.235046, -2.84217*10^-14}
expr'[root]
461.91
Notice the large derivative and how far off (vertically) is the plot's
closest point to the root.
highLow = {root, #} & /@ pt1[[{1, -1}, -1]];
Show[p0, p1, p2, Graphics@{Thick, Red, Line@highLow}]
Bobby
On Tue, 18 Jan 2011 04:51:50 -0600, Andrew DeYoung
<adey...@andrew.cmu.edu> wrote:
> Hi,
>
> I am trying to find the root of a certain expression in Mathematica
> version 7:
>
> expr = 110.52499999999998 + (300. - 135.52499999999998/(1 - x)) (1 -
> x) - 300. x - 135.52499999999998 Log[1 - x] + 135.52499999999998
> Log[x]
>
> It appears to plot fine, for example using Plot[expr, {x, 0, 1}]. The
> plot shows that there should be a root at about x=0.85. However, when
> I try to find this root, using for example the following:
>
> FindRoot[expr, {x, 0.5}]
>
> I get an error message:
>
> "FindRoot::lstol: The line search decreased the step size to within
> tolerance specified by AccuracyGoal and PrecisionGoal but was unable
> to find a sufficient decrease in the merit function. You may need
> more than MachinePrecision digits of working precision to meet these
> tolerances."
>
> and it prints a seemingly incorrect (according to the qualitative form
> of the plot) result: {x -> 0.344678}. Only if I use for example
>
> FindRoot[expr, {x, 0.7}]
>
> do I get the seemingly "correct" root: {x -> 0.849823}.
>
> Can you help me see why the FindRoot is getting stuck at {x ->
> 0.344678} when I use starting values far away from 0.7 or 0.8? I will
> ultimately want to find the roots of many similar functions, which may
> have more than one "actual" root, so it would be helpful if I could
> see why FindRoot[expr, {x, 0.5}] does not give {x -> 0.849823}. (also
> when I tried NSolve[expr==0,x], Mathematica will not solve it.)
>
> Thank you,
>
> Andrew DeYoung
> Carnegie Mellon University
>
>I am trying to find the root of a certain expression in Mathematica
>version 7:
>expr = 110.52499999999998 + (300. - 135.52499999999998/(1 - x)) (1 -
>x) - 300. x - 135.52499999999998 Log[1 - x] + 135.52499999999998
>Log[x]
>It appears to plot fine, for example using Plot[expr, {x, 0,
>1}]. The plot shows that there should be a root at about
>x=0.85. However, when I try to find this root, using for
>example the
>following:
>FindRoot[expr, {x, 0.5}]
>I get an error message:
<snip>
>FindRoot[expr, {x, 0.7}]
>do I get the seemingly "correct" root: {x -> 0.849823}.
>
>Can you help me see why the FindRoot is getting stuck at {x ->
>0.344678} when I use starting values far away from 0.7 or 0.8?
With only one starting point, FindRoot uses Newton's method to
find the root. Newton's method simply is not guaranteed to
converge for arbitrary starting points.
If you use the EvaluationMonitor as follows:
In[23]:= {res, {evx}} =
Reap[FindRoot[expr == 0, {x, .5}, EvaluationMonitor :> Sow[x]]];
evx[[;; 5]]
Out[24]= {0.5,0.0682211,0.456822,0.0377306,0.414913}
you can see with 0.5 as a starting point, the algorithm is not
converging as you want. Contrast that with starting at 0.7
In[25]:= {res, {evx}} =
Reap[FindRoot[expr == 0, {x, .7}, EvaluationMonitor :> Sow[x]]];
evx[[;; 5]]
Out[26]= {0.7,1.36516,0.766516,0.917897,0.793223}
With Newton's method success won't occur if the starting point
is to far away from the root.
There are a couple of ways around this difficulty. First you can
bound the root and do:
In[27]:= FindRoot[expr == 0, {x, .5, .9}]
Out[27]= {x->0.849823}
With two initial points which bound the root, FindRoot uses the
secant method which will always converge on the desired root
when there is only one root in the interval.
Another approach would be to use NMimimze as follows:
In[28]:= NMinimize[{Abs[expr], -1 < x < 1}, x]
Out[28]= {5.86936*10^-7,{x->0.849823}}
NMinimize uses different algorithms than FindRoot which are
often more robust for finding a root in this manner. However,
this will not always be the case without some tweaking of the
settings used by NMimimize.
Use the form of FindRoot that avoids use of derivatives
FindRoot[expr, {x, 0.5, 0.9}]
{x -> 0.849823}
Alternatively,
expr2 = Rationalize[expr] // FullSimplify
275 - 600*x - (5421/20)*ArcTanh[1 - 2*x]
Reduce[expr2 == 0, x, Reals] // N // ToRules
{x -> 0.849823}
Bob Hanlon
---- Andrew DeYoung <adey...@andrew.cmu.edu> wrote:
=============
Hi,
I am trying to find the root of a certain expression in Mathematica
version 7:
expr = 110.52499999999998 + (300. - 135.52499999999998/(1 - x)) (1 -
x) - 300. x - 135.52499999999998 Log[1 - x] + 135.52499999999998
Log[x]
It appears to plot fine, for example using Plot[expr, {x, 0, 1}]. The
plot shows that there should be a root at about x=0.85. However, when
I try to find this root, using for example the following:
FindRoot[expr, {x, 0.5}]
I get an error message:
"FindRoot::lstol: The line search decreased the step size to within
tolerance specified by AccuracyGoal and PrecisionGoal but was unable
to find a sufficient decrease in the merit function. You may need
more than MachinePrecision digits of working precision to meet these
tolerances."
and it prints a seemingly incorrect (according to the qualitative form
of the plot) result: {x -> 0.344678}. Only if I use for example
FindRoot[expr, {x, 0.7}]
do I get the seemingly "correct" root: {x -> 0.849823}.
Can you help me see why the FindRoot is getting stuck at {x ->
For instance:
expr = 110.52499999999998 + (300. - 135.52499999999998/(1 - x)) (1 -
x) - 300. x - 135.52499999999998 Log[1 - x] +
135.52499999999998 Log[x];
<< Optimization`UnconstrainedProblems`
FindRootPlot[expr, {x, 0.28}, PlotRange -> {{0, 1}, {-100, 100}}]
Bobby
On Wed, 19 Jan 2011 04:26:20 -0600, Gabriel Landi <gtl...@gmail.com>
wrote:
> Hey Andrew,
>
> The problem is likely the direction it uses to begin the search.
> x = 0.655 is a minimum. If you give a initial value above that it finds
> your
> root.
> If you give a initial value below that it will go the other way which
> eventually diverges when x approaches zero.
>
> You can take a better look using the following two lines of code:
>
> << Optimization`UnconstrainedProblems`
>
> FindRootPlot[expr, {x, 0.6}, PlotRange -> {{0, 1}, {-100, 100}}]
>
>
> (The output graph comes out bugged but you can still see what you need to
> see: the steps it took).
>
> The colors are explained in
> http://reference.wolfram.com/mathematica/tutorial/UnconstrainedOptimizationPlottingSearchData.html
>
>
> By the way, FindInstance works:
>
>
>
> FindInstance[expr == 0, x, Reals] // N
>
>
> Cheers,
>
>
> Gabriel
>
>
>
> On Tue, Jan 18, 2011 at 8:51 AM, Andrew DeYoung
> <adey...@andrew.cmu.edu>wrote:
>
Clear[expr]
expr[x_] =
110.52499999999998 + (300. - 135.52499999999998/(1 - x)) (1 - x) -
300. x - 135.52499999999998 Log[1 - x] +
135.52499999999998 Log[x] // Rationalize // FullSimplify
g[x_] = x - expr[x]/expr'[x] // FullSimplify;
Plot[{expr@x, g[x]}, {x, 0, 1}]
g has poles near .34 and .66, where expr'[x] is 0 at local minima/maxima.
For Newton's method to work reliably all by itself, g should be a
contraction, meaning Abs[g'[x]] < 1, in an interval containing both the
root and the starting point.
That's not true at most x in this case, as this plot will show:
Plot[{1, Abs[g'[x]]}, {x, 0, 1}, PlotRange -> {0, 2}]
Hence, you can't expect Newton to do a lot of good (except by luck) unless
the starting point is between .8 and .95 (from the graph), or more
accurately:
FindRoot[Abs[g'[x]] == 1, {x, .8}]
{x -> 0.795996}
and
FindRoot[Abs[g'[x]] == 1, {x, .95}]
{x -> 0.942448}
FindRoot is smarter than just "Newton's method", but we'd at least want to
start above the local minimum near .66.
FindRoot[expr@x == 0, {x, .7}]
{x -> 0.849823}
Bobby
On Wed, 19 Jan 2011 04:26:31 -0600, Daniel Lichtblau <da...@wolfram.com>
wrote:
> Andrew DeYoung wrote:
>> Hi,
>>
>> I am trying to find the root of a certain expression in Mathematica
>> version 7:
>>
>> expr = 110.52499999999998 + (300. - 135.52499999999998/(1 - x)) (1 -
>> x) - 300. x - 135.52499999999998 Log[1 - x] + 135.52499999999998
>> Log[x]
>>
>> It appears to plot fine, for example using Plot[expr, {x, 0, 1}]. The
>> plot shows that there should be a root at about x=0.85. However, when
>> I try to find this root, using for example the following:
>>
>> FindRoot[expr, {x, 0.5}]
>>
>> I get an error message:
>>
>> "FindRoot::lstol: The line search decreased the step size to within
>> tolerance specified by AccuracyGoal and PrecisionGoal but was unable
>> to find a sufficient decrease in the merit function. You may need
>> more than MachinePrecision digits of working precision to meet these
>> tolerances."
>>
>> and it prints a seemingly incorrect (according to the qualitative form
>> of the plot) result: {x -> 0.344678}. Only if I use for example
>
> The plot also shows that 0.5 is not in the basin of attraction for
> Newton's method to find the root, and that it will send one in the
> direction of the local max which is around x=0.35
>
> Other methods might work better. For example, if you provide two start
> points that bracket a root, then Brent's method will be used and it will
> successfully find this root.
>
> In[84]:= FindRoot[expr, {x, 0.5, .9}]
> Out[84]= {x -> 0.849823}
>
>
>> FindRoot[expr, {x, 0.7}]
>>
>> do I get the seemingly "correct" root: {x -> 0.849823}.
>
> The zone of attraction for the root is to the right of the local min,
> which appears to be somewhere around x=0.65.
>
>
>> Can you help me see why the FindRoot is getting stuck at {x ->
>> 0.344678} when I use starting values far away from 0.7 or 0.8? I will
>> ultimately want to find the roots of many similar functions, which may
>> have more than one "actual" root, so it would be helpful if I could
>> see why FindRoot[expr, {x, 0.5}] does not give {x -> 0.849823}. (also
>> when I tried NSolve[expr==0,x], Mathematica will not solve it.)
>>
>> Thank you,
>>
>> Andrew DeYoung
>> Carnegie Mellon University
>
> Solve or NSolve with restrictions on the domain will handle this (in
> version 8).
>
> In[82]:= Solve[{expr == 0, 0 <= x <= 1}, x]
> Out[82]= {{x -> 0.849823}}
>
> Daniel Lichtblau
> Wolfram Research
>
Thank you very much for all your time and wonderful help! This
discussion was immensely helpful to me.
Thank you,
Andrew
Carnegie Mellon University