MathGroup Archive 2010

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

Search the Archive

Re: FindRoots?

  • To: mathgroup at smc.vnet.net
  • Subject: [mg112134] Re: FindRoots?
  • From: "David Park" <djmpark at comcast.net>
  • Date: Wed, 1 Sep 2010 06:28:19 -0400 (EDT)

Three comments:

1) This appears to arise from new capabilities added to Root in Version 7,
capabilities that may have slipped under the radar screen of many users,
certainly under mine.

2) Nice work!

3) WRI can't provide a specific command for everything and often users will
have to add their own definitions, but 1-dimensional root searching of real
functions seems one place a conveniently packaged routine would be useful.


David Park
djmpark at comcast.net
http://home.comcast.net/~djmpark/  



From: Daniel Lichtblau [mailto:danl at wolfram.com] 

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




  • Prev by Date: Re: FindRoots?
  • Next by Date: Re: FindRoots?
  • Previous by thread: Re: FindRoots?
  • Next by thread: Re: FindRoots?