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