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

Re: FindRoots?

• To: mathgroup at smc.vnet.net
• Subject: [mg112103] Re: FindRoots?
• From: Daniel Lichtblau <danl at wolfram.com>
• Date: Tue, 31 Aug 2010 04:14:28 -0400 (EDT)

David Park wrote:
> Install Ted Ersek's package RootSearch from:
> 
> http://library.wolfram.com/infocenter/MathSource/4482/ 
> 
> I don't know why Wolfram doesn't acquire the rights to this package and
> incorporate it as part of regular Mathematica. Finding all the roots of a
> real function on a 1-dimensional domain is the most common case of root
> finding, and it is the one thing regular Mathematica is very poor at. Ted's
> package is quite robust and returns all the roots in order.
> 
> 
> David Park
> djmpark at comcast.net
> http://home.comcast.net/~djmpark/  
> 
> 
> 
> From: Sam Takoy [mailto:sam.takoy at yahoo.com] 
> 
> Hi,
> 
> Is there a command for numerically finding all roots of a function in a 
> given interval. If not, what's the best way of accomplishing this task?
> 
> Many thanks in advance,
> 
> Sam

If I instead define RootSearch to use Reduce, as below, it handles quite 
well all the examples from the Examples section of 
RootSearchExamples.nb. I show them below. I modified one to use 
Piecewise instead of a condition, and for the last I use Reduce directly 
to get an expanded form of the solution set.

I think Mathematica is quite good at this sort of thing, actually.

Daniel Lichtblau
Wolfram Research

-------------------------------

In[33]:= RootSearch[
   eqn_, {x_, lo_, hi_}] := {ToRules[Reduce[{eqn, lo <= x <= hi}]]}

In[34]:= f1[x_] := Sin[4 x] - (x + 1)/8
RootSearch[f1[x] == 0, {x, -13, 12}]

Out[35]= {{x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -8.3482893782997216496}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -8.1289315341182707392}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -6.8629511940895252882}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -6.4714706256130258036}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -5.3539165396979340111}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -4.8374622685220237568}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -3.8363831890811367726}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -3.21161777565464349626}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -2.3149160596763252331}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -1.58922633877627429127}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, -0.79190196071777069054}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 0.032351189053103634282}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 0.73087710353047892233}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 1.65538160232955244046}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 2.2515547207396359640}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 3.2828240626757282491}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 3.7673870629720211202}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 4.9206885635879098258}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 5.2724915094753945872}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 6.5961054327905341135}]}, {x ->
    Root[{1 - 8 Sin[4 #1] + #1 &, 6.7398184613474433969}]}}

In[36]:= f2[x_] := If[x < 3/2, Sqrt[3/2 - x] - Exp[-4 x], 2* x Exp[-x]]
RootSearch[f2[x] == 0, {x, -25, 25}]

Out[37]= {{x ->
    Root[{-Sqrt[2] +
        E^(4 #1) Sqrt[3 - 2 #1] &, -0.055200564731044670903}]}, {x ->
    Root[{-Sqrt[2] + E^(4 #1) Sqrt[3 - 2 #1] &,
      1.49999385548561363888}]}}

In[38]:= RootSearch[Exp[x - \[Pi]] == 1 - \[Pi] + x, {x, -3, 4}]

Out[38]= {{x -> \[Pi]}}

In[39]:= f4[x_] := (1.5 - Erf[x] - Erf[2 - x]) Exp[-Abs[x - 1]]
RootSearch[f4[x] == 0, {x, -400, 200}]

During evaluation of In[39]:= Reduce::ratnz: Reduce was unable to solve 
the system with inexact coefficients. The answer was obtained by solving 
a corresponding exact system and numericizing the result. >>

Out[40]= {{x -> 0.517891}, {x -> 1.48211}}

In[41]:= f5[x_] :=
  3.2 - 3*Cos[x] - Exp[-(x - 8 \[Pi])^2] - Exp[-(x - 14 \[Pi])^2]
s5 = RootSearch[f5[x] == 0, {x, 0, 30 \[Pi]}]

During evaluation of In[41]:= Reduce::ratnz: Reduce was unable to solve 
the system with inexact coefficients. The answer was obtained by solving 
a corresponding exact system and numericizing the result. >>

Out[42]= {{x -> 24.5432}, {x -> 25.7223}, {x -> 43.3928}, {x ->
    44.5718}}

In[43]:= f6[x_] :=
  Piecewise[{{ Sqrt[Abs[x] - 1.2], x <= 1.4}}, Indeterminate]
RootSearch[f6[x] == 0, {x, -8, 8}]

During evaluation of In[43]:= Reduce::ratnz: Reduce was unable to solve 
the system with inexact coefficients. The answer was obtained by solving 
a corresponding exact system and numericizing the result. >>

Out[44]= {{x -> -1.2}, {x -> 1.2}}

In[45]:= RootSearch[Zeta[1/2 + I y] == 0, {y, 0, 237}]

Out[45]= {ToRules[
   C[1] \[Element] Integers && 1 <= C[1] <= 100 &&
    y == -(1/2) I (-1 + 2 ZetaZero[C[1]])]}

In[48]:= explicitsol =
  Resolve[Exists[C[1],
    Reduce[{Zeta[1/2 + I y] == 0, 0 <= y <= 237}, y]]]

Out[48]= y == Im[ZetaZero[1]] || y == Im[ZetaZero[2]] ||
  y == Im[ZetaZero[3]] || y == Im[ZetaZero[4]] ||
  y == Im[ZetaZero[5]] || y == Im[ZetaZero[6]] ||
  y == Im[ZetaZero[7]] || y == Im[ZetaZero[8]] ||
  y == Im[ZetaZero[9]] || y == Im[ZetaZero[10]] ||
  y == Im[ZetaZero[11]] || y == Im[ZetaZero[12]] ||
  y == Im[ZetaZero[13]] || y == Im[ZetaZero[14]] ||
  y == Im[ZetaZero[15]] || y == Im[ZetaZero[16]] ||
  y == Im[ZetaZero[17]] || y == Im[ZetaZero[18]] ||
  y == Im[ZetaZero[19]] || y == Im[ZetaZero[20]] ||
  y == Im[ZetaZero[21]] || y == Im[ZetaZero[22]] ||
  y == Im[ZetaZero[23]] || y == Im[ZetaZero[24]] ||
  y == Im[ZetaZero[25]] || y == Im[ZetaZero[26]] ||
  y == Im[ZetaZero[27]] || y == Im[ZetaZero[28]] ||
  y == Im[ZetaZero[29]] || y == Im[ZetaZero[30]] ||
  y == Im[ZetaZero[31]] || y == Im[ZetaZero[32]] ||
  y == Im[ZetaZero[33]] || y == Im[ZetaZero[34]] ||
  y == Im[ZetaZero[35]] || y == Im[ZetaZero[36]] ||
  y == Im[ZetaZero[37]] || y == Im[ZetaZero[38]] ||
  y == Im[ZetaZero[39]] || y == Im[ZetaZero[40]] ||
  y == Im[ZetaZero[41]] || y == Im[ZetaZero[42]] ||
  y == Im[ZetaZero[43]] || y == Im[ZetaZero[44]] ||
  y == Im[ZetaZero[45]] || y == Im[ZetaZero[46]] ||
  y == Im[ZetaZero[47]] || y == Im[ZetaZero[48]] ||
  y == Im[ZetaZero[49]] || y == Im[ZetaZero[50]] ||
  y == Im[ZetaZero[51]] || y == Im[ZetaZero[52]] ||
  y == Im[ZetaZero[53]] || y == Im[ZetaZero[54]] ||
  y == Im[ZetaZero[55]] || y == Im[ZetaZero[56]] ||
  y == Im[ZetaZero[57]] || y == Im[ZetaZero[58]] ||
  y == Im[ZetaZero[59]] || y == Im[ZetaZero[60]] ||
  y == Im[ZetaZero[61]] || y == Im[ZetaZero[62]] ||
  y == Im[ZetaZero[63]] || y == Im[ZetaZero[64]] ||
  y == Im[ZetaZero[65]] || y == Im[ZetaZero[66]] ||
  y == Im[ZetaZero[67]] || y == Im[ZetaZero[68]] ||
  y == Im[ZetaZero[69]] || y == Im[ZetaZero[70]] ||
  y == Im[ZetaZero[71]] || y == Im[ZetaZero[72]] ||
  y == Im[ZetaZero[73]] || y == Im[ZetaZero[74]] ||
  y == Im[ZetaZero[75]] || y == Im[ZetaZero[76]] ||
  y == Im[ZetaZero[77]] || y == Im[ZetaZero[78]] ||
  y == Im[ZetaZero[79]] || y == Im[ZetaZero[80]] ||
  y == Im[ZetaZero[81]] || y == Im[ZetaZero[82]] ||
  y == Im[ZetaZero[83]] || y == Im[ZetaZero[84]] ||
  y == Im[ZetaZero[85]] || y == Im[ZetaZero[86]] ||
  y == Im[ZetaZero[87]] || y == Im[ZetaZero[88]] ||
  y == Im[ZetaZero[89]] || y == Im[ZetaZero[90]] ||
  y == Im[ZetaZero[91]] || y == Im[ZetaZero[92]] ||
  y == Im[ZetaZero[93]] || y == Im[ZetaZero[94]] ||
  y == Im[ZetaZero[95]] || y == Im[ZetaZero[96]] ||
  y == Im[ZetaZero[97]] || y == Im[ZetaZero[98]] ||
  y == Im[ZetaZero[99]] || y == Im[ZetaZero[100]]

In[49]:= N[explicitsol]

Out[49]= y == 14.1347 || y == 21.022 || y == 25.0109 || y == 30.4249 ||
   y == 32.9351 || y == 37.5862 || y == 40.9187 || y == 43.3271 ||
  y == 48.0052 || y == 49.7738 || y == 52.9703 || y == 56.4462 ||
  y == 59.347 || y == 60.8318 || y == 65.1125 || y == 67.0798 ||
  y == 69.5464 || y == 72.0672 || y == 75.7047 || y == 77.1448 ||
  y == 79.3374 || y == 82.9104 || y == 84.7355 || y == 87.4253 ||
  y == 88.8091 || y == 92.4919 || y == 94.6513 || y == 95.8706 ||
  y == 98.8312 || y == 101.318 || y == 103.726 || y == 105.447 ||
  y == 107.169 || y == 111.03 || y == 111.875 || y == 114.32 ||
  y == 116.227 || y == 118.791 || y == 121.37 || y == 122.947 ||
  y == 124.257 || y == 127.517 || y == 129.579 || y == 131.088 ||
  y == 133.498 || y == 134.757 || y == 138.116 || y == 139.736 ||
  y == 141.124 || y == 143.112 || y == 146.001 || y == 147.423 ||
  y == 150.054 || y == 150.925 || y == 153.025 || y == 156.113 ||
  y == 157.598 || y == 158.85 || y == 161.189 || y == 163.031 ||
  y == 165.537 || y == 167.184 || y == 169.095 || y == 169.912 ||
  y == 173.412 || y == 174.754 || y == 176.441 || y == 178.377 ||
  y == 179.916 || y == 182.207 || y == 184.874 || y == 185.599 ||
  y == 187.229 || y == 189.416 || y == 192.027 || y == 193.08 ||
  y == 195.265 || y == 196.876 || y == 198.015 || y == 201.265 ||
  y == 202.494 || y == 204.19 || y == 205.395 || y == 207.906 ||
  y == 209.577 || y == 211.691 || y == 213.348 || y == 214.547 ||
  y == 216.17 || y == 219.068 || y == 220.715 || y == 221.431 ||
  y == 224.007 || y == 224.983 || y == 227.421 || y == 229.337 ||
  y == 231.25 || y == 231.987 || y == 233.693 || y == 236.524