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