MathGroup Archive 2007

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

Search the Archive

Re: Re: Debug of FindRoot

  • To: mathgroup at smc.vnet.net
  • Subject: [mg78949] Re: [mg78882] Re: Debug of FindRoot
  • From: DrMajorBob <drmajorbob at bigfoot.com>
  • Date: Fri, 13 Jul 2007 06:09:38 -0400 (EDT)
  • References: <f72c9h$d75$1@smc.vnet.net> <26979812.1184251892863.JavaMail.root@m35>
  • Reply-to: drmajorbob at bigfoot.com

Rather than 5 starting points for 3 roots, the following code looks at 4 
starting points and finds 4 roots. It finds the root after 10, which  
Miguel seemed to want, and it feeds bracketed ranges to FindRoot, rather 
than single starting points. I also used Cases[plot,_Line,Infinity] to  
eliminate number pairs from various Plot options.

(f[z_] = Numerator[Together[BesselJ[0, z]/BesselJ[1, z] - z/(1/2)]]);
Print[StyleForm[
   "the function numerator is\nf[z] = " <> ToString[f[z]],
   FontColor -> Blue]]
plot = Plot[f[x], {x, 0.1, 11}]
Print[StyleForm[
    "bracketed roots (from the Plot) in the form {mid,left,right},\n\
where f[left]*f[right]<=0 and mid=(left+right)/2",
    FontColor -> Blue]];
brackets =
  Cases[Cases[
     Cases[plot, _Line, Infinity], {(x_)?NumberQ, (y_)?NumberQ},
     Infinity] //. {a___, {x1_, y1_}, {x2_, y2_}, b___} /;
      y1 y2 <= 0 :> {a, {(x1 + x2)/2, x1, x2}, b}, {_, _, _}]
Print[StyleForm["the roots are", FontColor -> Blue]];
(FindRoot[f[x] == 0, {x, Sequence @@ #}, WorkingPrecision -> 30,
     PrecisionGoal -> 10] &) /@ brackets

Bobby

On Thu, 12 Jul 2007 04:10:58 -0500, dimitris <dimmechan at yahoo.com> wrote:

> The reason is clear from the following plot
>
> Plot[BesselJ[0, z]/BesselJ[1, z], {z, 0.1, 12}]
>
> and as well
>
> << "NumericalMath`BesselZeros`"
>
> BesselJZeros[1,4]
> {3.83171,7.01559,10.1735,13.3237}
>
> That is the zeros of the denominator (at the roots
> of BesselJ[1,z] ) causes problems to Mathematica.
>
> Work as follows instead
>
> (*copy paste the following as a whole cell; select the
> cell and press Shift+Enter*)
>
> Print[StyleForm["the function", FontColor -> Blue]]
> f[z_] = Numerator[Together[BesselJ[0, z]/BesselJ[1, z] - z/(1/2)]]
> plot = Plot[f[x], {x, 0.1, 10}];
> Print[StyleForm["the points used by the Plot function", FontColor ->
> Blue]];
>   Short[points = Cases[plot, {(x_)?NumberQ, (y_)?NumberQ}, Infinity]]
> Print[StyleForm["find where the function changes sign", FontColor ->
> Blue]];
>   seeds = Position[Apply[Times, Partition[points[[All,2]], 2, 1],
> {1}], x_ /; x <= 0]
> Print[StyleForm["between this points in x axis there is a change in
> sign of f[x]", FontColor -> Blue]];
>   samples = Extract[Partition[points[[All,1]], 2, 1], seeds]
> Print[StyleForm["the roots, at last!", FontColor -> Blue]];
>   (FindRoot[f[x] == 0, {x, #1[[1]], #1[[2]]}, WorkingPrecision -> 30,
> PrecisionGoal -> 10] & ) /@ samples
>
> (*Above I use the points used by the Plot command in order
> to plot the function. It is found between what points in the
> requested
> interval the function changes sign. This points are used as starting
> points for the secant method*)
>
> Regards
> Dimitris
>
>
>
>             Miguel       :
>> To resolve one of the heat equations it is necesary to calculate the
>> solution of z for BesselJ[0,z]/BesselJ[1,z]==z/Bi, where Bi is the
>> Biot number (equal to 0.5, for example).
>>
>> 1.- Plot[{BesselJ[0,z]/BesselJ[1,z],z/Bi},{z,0.001,12}].
>>
>> >From this plot I deduce the ranges, more or less, {1,4,7,10}.
>>
>> 2.- FindRoot[BesselJ[0,z]/BesselJ[1,z]==z/Bi,{z,#}]&/@{1,4,7,10}
>> {{z->0.940771},{z->3.95937},{z->0.940771},{z->3.95937}}
>>
>> I dont understand the reason. With others differents intervals it
>> works fine.
>
>
>



-- 

DrMajorBob at bigfoot.com


  • Prev by Date: annoying documentation in 6 (rant)
  • Next by Date: Re: Palette to write out Global Names to a Selected NB
  • Previous by thread: Re: Debug of FindRoot
  • Next by thread: The Wolfram Blog Blogged: Analog Clock