Group Abstract Group Abstract

Message Boards Message Boards

0
|
6.2K Views
|
9 Replies
|
0 Total Likes
View groups...
Share
Share this post:

Compute accurately Exp[x^2]*Erfc[x] in a wide range of x?

Leslaw Bieniasz
Posted 9 years ago

Hi,

I need to compute accurately (with let's say at least 25 significant digits) the values of several expressions involving exponential function, in particular the expression Exp[x^2]*Erfc[x], in a wide range of x (ideally the range should correspond to the range of extended precision variables on modern PCs). Unfortunately, MATHEMATICA seems to fail performing such a task, owing to underflow or overflow errors in Exp[ ], even if I increase $MaxExtraPrecision considerably. I need to note that the above expression converges to zero very slowly (with increasing x), much more slowly than Exp[-x]. It's a pity that there seems to be no special function implemented in MATHEMATICA, for representing the above expression, despite the fact that the expression occurs frequently in many applications. A similar problem arises with expressions such as Exp[-x]BesselI[0,x] or Exp[-x]BesselI[1,x] . Is there any way to overcome the above difficulty?

Leslaw
POSTED BY: Leslaw Bieniasz
9 Replies
Sort By:

Now we testing simple case and verify are they correct.

DIGITPRECISION = 100;
ORDER = 20;
NUMBER = 10^7;
sol2 = Table[ser2 = Normal@Series[Exp[x^2]*Erfc[x], {x, Infinity, k}]; 
s2 = Block[{$MaxExtraPrecision = Infinity}, N[ser2 /. x -> NUMBER, DIGITPRECISION]]; 
ANSWER = s2, {k, 1, ORDER}]; m = Last[Counts[sol2]]; 
If[TrueQ[m >= 3], "The required precision is obtained", "Not enough ,increase ORDER"]
(*The required precision is obtained*)
sol2

enter image description here

In red column last digits (...4046517361) repeated 6 times ,Oder 15 it's enough to precision 100.

You can become more secure and write in last "code":

 If[TrueQ[m >= 5], "The required precision is obtained", "Not enough ,increase ORDER"]

Ps. My English is bad sorry :P

POSTED BY: Mariusz Iwaniuk

Maybe like so:

    DIGITPRECISION = 100000;
    ORDER = 30;
    NUMBER = 10^4900;
    sol2 = Table[ ser2 = Normal@Series[Exp[x^2]*Erfc[x], {x, Infinity, k}]; 
    s2 = Block[{$MaxExtraPrecision = Infinity}, N[ser2 /. x -> NUMBER, DIGITPRECISION]]; 
    ANSWER = s2, {k, 1, ORDER}]; m = Last[Counts[sol2]];
    If[TrueQ[m >= 3], "The required precision is obtained", "Not enough ,increase ORDER"]

   (*The required precision is obtained*)

   ANSWER
   (* 5.6418958354...34303*10^-4901*)
POSTED BY: Mariusz Iwaniuk

I have tried that Normal@Series approach, and it seems that I can get results up to x=10^5000 by this method (although I am not sure how to verify that they are correct). I wonder why taking the logarithm is supposed to be advantageous over just calculating the asymptotic series for Exp[x^2]*Erfc[x] ?

Les?aw
POSTED BY: Leslaw Bieniasz

Hi,

The example posted by Mariusz exactly touches the problem, showing the overflow. As I wrote, I need to calculate Exp[x^2]*Erfc[x] in the entire range of x corresponding to extended precision variables. The range is up to about 10^4900. Therefore, the values such as x=10^7 or x=10^8 are far too small.

Les?aw
POSTED BY: Leslaw Bieniasz

Exp[-x] does not strike me as problematic.

In[1168]:= N[Exp[-1000], 20]

Out[1168]= 5.0759588975494567653*10^-435

If this is not the sort of example you had in mind then (once again) I note that it would be helpful to have specific examples.
POSTED BY: Daniel Lichtblau

As a power series expansions: Testing:

Block[{$MaxExtraPrecision = Infinity},  N[Exp[x^2]*Erfc[x] /. x -> 10^7, 100]]
(* 5.64189583547753466000161712821652064212873489131706443451983934426010
2089725651175753054274046517361*10^-8 *)

n = 14; (*The result agrees to 100 decimal places*)
ser = Normal@Series[Log[Exp[x^2]*Erfc[x]], {x, Infinity, n}]; s = Block[{$MaxExtraPrecision = Infinity}, 
N[ser /. x -> 10^7, 100]]; Exp[s]
(* 5.64189583547753466000161712821652064212873489131706443451983934426010
208972565117575305427404651736*10^-8 *)
POSTED BY: Mariusz Iwaniuk

An example that be problematic:

 Block[{$MaxExtraPrecision = Infinity}, N[Exp[x^2]*Erfc[x] /. x -> 10^8, 20]]

 (*General::ovfl: Overflow occurred in computation.*)
POSTED BY: Mariusz Iwaniuk

Thanks. Yes, I assume real-valued x. I doubt if expansion of Log[1+y] as a series of y can be useful for any computing, as the series converges extremely slowly. Power series of x, well, yes, the expansions are well known, but my question was more about how to force MATHEMATICA into a truly multiprecision mode of operation in this case. It's also frustrating why there is no special function implemented for this expression. As far as I can see, another similar expression known as the Dawson integral has already been included in MATHEMATICA. So, why not including Exp[x^2]*Erfc[x] ? I might reformulate my question in another way: why one cannot compute Exp[-x] with arbitrary precision for arbitrary x? This appears impossible due to the existing limitation on the range of real values (if I understand the situation). Les?aw
POSTED BY: Leslaw Bieniasz

Examples of what is actually wanted would be quite helpful. Is there a restriction to real-valued x?

It might help to take the logarithm, and perhaps expand Log[Erfc[x]] as a series in Erf[x] (depending on specifics of x). But it looks like the precision needed far exceeds the size of inputs so maybe a power series of the original expression should be used for large x.
POSTED BY: Daniel Lichtblau
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Tag limit exceeded
Note: Only the first five people you tag will receive an email notification; the other tagged names will appear as links to their profiles.
Reply Preview
Attachments
Remove
or Discard