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

Questions
,
1

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

2

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.

3

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

4

An example that be problematic:

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

 (*General::ovfl: Overflow occurred in computation.*)
5

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 *)
6

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.

7

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

8

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

9

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*)
10

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

4830rys.png
4830rys.png1060×437 8.46 KB

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 :stuck_out_tongue:
[1]: http://community.wolfram.com//c/portal/getImageAttachment?filename=4830rys.png&userId=476423