Services & Resources / Wolfram Forums
-----
 /
MathGroup Archive
2002
*January
*February
*March
*April
*May
*June
*July
*August
*September
*October
*November
*December
*Archive Index
*Ask about this page
*Print this page
*Give us feedback
*Sign up for the Wolfram Insider

MathGroup Archive 2002

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Find many Roots

  • To: mathgroup at smc.vnet.net
  • Subject: [mg32181] Re: [mg32098] Find many Roots
  • From: Daniel Lichtblau <danl at wolfram.com>
  • Date: Fri, 4 Jan 2002 05:03:56 -0500 (EST)
  • References: <200112240832.DAA15698@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

Yoram Pollack wrote:
> 
> Hello,
> I am trying to find all roots of the Equation: Tan[x]==1/x, {x,0,12}].
> What is the best way?
> Thanks
> Yoram

Several nice methods have been posted. I thought I would draw attention
to a similar question from last June.

http://library.wolfram.com/mathgroup/archive/2001/Jun/msg00226.html

The biggest difference is that the present function has poles. One can
instead find zeroes of x*Sin[x]-Cos[x] to avoid this issue.

My own reply to the June post (URL below) gave one method for attacking
such problems, using the Cauchy integral formula to count roots in an
interval, bisection to isolate them, and CIF again to approximate
individually isolated roots.

http://library.wolfram.com/mathgroup/archive/2001/Jun/msg00226.html

Since the code was not too long I just include it below. The actual note
goes into more detail regarding some of the advantages, disadvantages,
and potential extensions of this approach.


epsilon = 10^(-1);
Off[NIntegrate::slwcon]

segmentIntegrate[func_,a_,b_,offset_,z_, f2_:1] := 
  Module[{t1,t2,t3,integrand=f2*D[func,z]/func},
    t1 = NIntegrate[
        Evaluate[(integrand/.z->z-I*offset)-
          (integrand/.z->z+I*offset)], {z,a,b},
		 PrecisionGoal->3, MaxRecursion->10];
    t2 = NIntegrate[Evaluate[integrand], {z,a+I*offset,a-I*offset},
        PrecisionGoal->3, MaxRecursion->10];
    t3 = NIntegrate[Evaluate[integrand], {z,b-I*offset,b+I*offset},
        PrecisionGoal->3, MaxRecursion->10];
    (t1+t2+t3)/(2*Pi*I)]

isolateRoots[func_,lower_,upper_,z_,eps_] := 
  Module[{numroots,avg=(lower+upper)/2,leftroots},
    numroots = Round[segmentIntegrate[func,lower,upper,eps,z]];
    If[numroots==0, Return[{}]];
    If[numroots==1, Return[{Interval[{lower,upper}]}]];
    leftroots = isolateRoots[func,lower,avg,z,eps];
    If[Length[leftroots]==numroots,Return[leftroots]];
    Join[leftroots,isolateRoots[func,avg,upper,z,eps]]
    ]

refineRoots[func_,intervals_,eps_,z_]:=
  Map[Re[segmentIntegrate[func,#[[1,1]],#[[1,2]],eps,z,z]]&, intervals]


Now we try it out.

f[x_] := x*Sin[x]-Cos[x]

InputForm[Timing[rootintervals =
  isolateRoots[f[z],0,12,z,epsilon]]]

This gives some useless messages regarding loss of precision, followed
by

Out[46]//InputForm= 
{0.03999999999999915*Second, {Interval[{0, 3}], Interval[{3, 6}], 
  Interval[{6, 9}], Interval[{9, 12}]}}

In[47]:= Timing[rts=refineRoots[f[z],rootintervals,epsilon,z]]
Out[47]= {0.02 Second, {0.860334, 3.42562, 6.4373, 9.52933}}

We see that the residuals are not fairly small.

In[53]:= InputForm[f[x] /. x->rts]
Out[53]//InputForm= 
{-3.984138952084493*^-8, 3.80665711174899*^-8, 3.887810040570372*^-8, 
 7.042598237916309*^-8}

We could of course do better either by higher precision quadrature or,
more sensibly, by taking these as starting values for a high
precision/accuracy FindRoot.

Other replies to the June query gave code to implement entirely
different methods e.g. polynomial approximation, interval search, and
maybe others as well; these all could be applied to this example as
well.


Daniel Lichtblau
Wolfram Research


  • Prev by Date: Re: ORDINARY DIFFERENTIAL EQUATION
  • Next by Date: Re: Front end problems!
  • Previous by thread: Re: Solutions that are not solutions
  • Next by thread: Re: Find many Roots