Re: LogErfc
- To: mathgroup at smc.vnet.net
- Subject: [mg122635] Re: LogErfc
- From: Daniel Lichtblau <danl at wolfram.com>
- Date: Fri, 4 Nov 2011 05:59:53 -0500 (EST)
- Delivered-to: l-mathgroup@mail-archive0.wolfram.com
- References: <201111030844.DAA15151@smc.vnet.net>
On 11/03/2011 03:44 AM, Shawn Garbett wrote:
> I need to have a LogErfc function to be able to compute another
> function in a numerically stable manner. I found the function in the
> GSL, I found it in R, and even in another system, there's a copy of code to do
> it floating around. However, I need it in Mathematica and Log[Erfc[x]]
> doesn't cut in due to underflow. Cody solved this in 1969, and I could
> code his solution, but I was hoping that someone had already solved
> this in Mathematica. Anyone?
>
If you are interested mostly in x either positive or not too negative,
could do Pade approximations at zero and infinity, and use whichever is
appropriate in given ranges.
padeZero[x_] = PadeApproximant[Log[Erfc[x]], {x,0,{5,5}}];
padeInfinity[x_] = Log[1/x] +
PadeApproximant[Log[Erfc[x]]-Log[1/x], {x,Infinity,{5,5}}];
logErfcSmallApprox[x_] := Block[{z},padeZero[z]/.z->x]
logErfcLargeApprox[x_] := Block[{z},padeInfinity[z]/.z->x]
Quick check:
InputForm[Table[{Log[Erfc[x]],logErfcSmallApprox[x],logErfcLargeApprox[x]},
{x, -3.5, 10.5, 1}]]
Out[238]//InputForm=
{{0.69314680901069, 1.0493788756124784, -14.112164948457105 +
3.141592653589793*I}, {0.6929436838471712, 0.716213867787085,
-7.805322341465522 + 3.141592653589793*I},
{0.6760545027339606, 0.6761720084686005, -3.370687193890008 +
3.141592653589793*I}, {0.41903914777555973, 0.41903914813762605,
-0.46255109569808817 + 3.141592653589793*I},
{-0.7350111298370846, -0.7350111299422797, -0.46255109569808817},
{-3.384492089551553, -3.3844975258920105, -3.370687193890008},
{-7.806815272727264, -7.8072778682726955, -7.805322341465522},
{-14.112437402148174, -14.118997281238098, -14.112164948457105},
{-22.34976830343594, -22.390749679580356, -22.34969815365446},
{-32.543008907376965, -32.70519933060154, -32.54298605103614},
{-44.705670189533116, -45.18860445062361, -44.70566137269984},
{-58.845967473872314, -60.03403469956622, -58.845963615640876},
{-74.96923568857265, -77.5233971405431, -74.96923382750941},
{-93.07912219224849, -98.05396915740401, -93.07912122240553},
{-113.1782250431496, -122.17524249202071, -113.17822450502089}}
There is clearly an overlap range where both approximations are good to
three or so places. If this is not enough, could use higher degree Pade
approximations.
For sufficiently negative x you can approximate by Log[2].
Daniel Lichtblau
Wolfram Research
- References:
- LogErfc
- From: Shawn Garbett <shawn.garbett@gmail.com>
- LogErfc