Re: Root finding needs higher accuracy

*To*: mathgroup at smc.vnet.net*Subject*: [mg123176] Re: Root finding needs higher accuracy*From*: DrMajorBob <btreat1 at austin.rr.com>*Date*: Fri, 25 Nov 2011 04:57:02 -0500 (EST)*Delivered-to*: l-mathgroup@mail-archive0.wolfram.com*References*: <201111241157.GAA29024@smc.vnet.net>*Reply-to*: drmajorbob at yahoo.com

Here's the function and an attempt using FindRoot: Clear[f] f[x_] = Log[x] + x^2/(2 Exp[1]) - 2 x/Sqrt[Exp[1]] + 1; x /. Quiet@ FindRoot[f[x], {x, 1.0, 3.4}, Method -> "Brent", WorkingPrecision -> 100, PrecisionGoal -> 25] // Timing N[Sqrt@Exp@1 - Last@%, 25] {0.003579, 1.64872127073570073780869588498038437600} -3.557259096004509716622080*10^-11 Looking at the plot: Plot[{f@x, f'[x], f''[x], f'''[x]}, {x, 1, 2}] suggests that it's a double, even TRIPLE root, and if that's true, the solution is easy: root = x /. First@Solve[f'[x] == 0, x] f@root Sqrt[E] 0 We can also solve this: Solve[f''[x] == 0, x] {{x -> -Sqrt[E]}, {x -> Sqrt[E]}} Sqrt[E] is the root. BUT, if we don't know it's a double root or Solve can't find a root of f', THEN what? Here's a method that works for this function (and maybe a few others): Clear[trisect] trisect[{low_, high_}] := Module[{tbl, split}, tbl = Table[{f@x, x}, {x, {low, low + (high - low)/3, low + 2 (high - low)/3, high}}]; split = SplitBy[tbl, Negative]; {split[[1, -1, -1]], split[[-1, 1, -1]]} ] precision = 100; iterations = 80; Timing@Mean@ Block[{$MinPrecision = precision, $MaxPrecision = precision}, NestWhile[trisect, SetPrecision[{1, 2}, precision], Abs[Subtract @@ #] > 10^-25 &, 1]] {0.012396, \ 1.64872127070012814684865080152405879098770003155385410176215449001051\ 1663472116782925272061059040521} Choosing adequate "precision" is a problem in itself, of course. Bobby On Thu, 24 Nov 2011 05:57:06 -0600, HwB <hwborchers at googlemail.com> wrote: > I would like to numerically find the root of the following function > with up to 20 digits. > > f12[x_] := Log[x] + x^2 / (2 Exp[1]) - 2 x / Sqrt[Exp[1]] + 1 > > This problem is known to be difficult for solvers in double precision > arithmetics. I thought it should be easy with Mathematica, but the > following attempts were not successful. > > SetPrecision[ > x /. FindRoot[f12[x], {x, 1.0, 3.4}, Method -> "Brent", > AccuracyGoal -> Infinity, PrecisionGoal -> 20], 16] > # 1.648732212532746 > SetPrecision[ > x /. FindRoot[f12[x], {x, 1.0, 3.4}, Method -> "Secant", > AccuracyGoal -> Infinity, PrecisionGoal -> 20], 16] > # 1.648710202030051 > > The true root obviously is Sqrt[Exp[1]]//N = 1.648721270700128... > > The symbolic solver explicitely says it cannot solve this expression. > What do I need to do to get a much more exact result out of > Mathematica? > > Many thanks, Hans Werner > -- DrMajorBob at yahoo.com

**References**:**Root finding needs higher accuracy***From:*HwB <hwborchers@googlemail.com>