2
$\begingroup$

In the Mathematica 11.1 documentation, they have given the following series expansion of the inverse error function. My question is: how does Mathematica compute it? Does it do an approximation of this series? If yes, how many terms do they include? Thanks.

Series[InverseErf[z], {z, 0, 10}]

(* (Sqrt[π] z)/2 + 1/24 π^(3/2) z^3 + 7/960 π^(5/2) z^5 + (127 π^(7/2) z^7)/80640 +
   (4369 π^(9/2) z^9)/11612160 + O[z]^11 *)
Improve this question
$\endgroup$
0
5
$\begingroup$

For inverse functions like InverseErf[], the usual approach for determining their series coefficients is Lagrangian inversion. (That is, what InverseSeries[] does.)

In particular, Mathematica supports the (partial) Bell polynomials that appear in the series coefficients for the inverse series as BellY[]. (See also the discussion in Charalambides's book.) Using this, along with the general formula for the derivative of Erf[] (formula 7.1.19 in Abramowitz and Stegun), we have

SetAttributes[inverseErfCoefficient, Listable];
inverseErfCoefficient[0] = 0; inverseErfCoefficient[1] = Sqrt[π]/2;
inverseErfCoefficient[m_Integer?Positive] :=
       (Sqrt[π]/2)^m BellY[Table[{(m + k - 2)!, (-2)^k/(2 k) Sqrt[π]/Gamma[1 - k/2]},
                                 {k, 2, m}]]/(m! (m - 1)!)

Test:

inverseErfCoefficient[Range[0, 7]]
   {0, Sqrt[π]/2, 0, π^(3/2)/24, 0, (7 π^(5/2))/960, 0, (127 π^(7/2))/80640}

Sum[inverseErfCoefficient[k] x^k, {k, 0, 20}] - (InverseErf[x] + O[x]^21)
   O[x]^21

As it turns out, the series coefficients for the Taylor expansion of the inverse error function at $0$ satisfy a simple nonlinear recurrence. The formulae used in the following are adapted from Carlitz and Dominici:

ck[0] = 0; ck[1] = 1;
ck[n_Integer?Positive] := ck[n] = 2 Sum[Binomial[n - 1, k] ck[k - 1] ck[n - k], {k, n - 1}]

SetAttributes[inverseErfCoefficient, Listable];
inverseErfCoefficient[m_Integer?NonNegative] := (Sqrt[π]/2)^m ck[m]/m!

This version should give the same results.

Improve this answer
$\endgroup$
7
  • $\begingroup$ An example of why your answers are among my favorites. +1 $\endgroup$ – ciao Nov 5 '17 at 6:32
  • $\begingroup$ Glad to know you like it! I try my best when it comes to special functions... $\endgroup$ – J. M.'s ennui♦ Nov 5 '17 at 7:25
  • $\begingroup$ Thanks @J.M., I tested the code . It answers my question. However for larger values of x, say x=.9, I have to keep k=61 which is high. For x=.99 value of k is really high and it takes more time to compute so I was looking at the alternative . I came across another approximate formula on wikipedia, nand tested , it works better. it goes like this: anotherErfInv[x_] := Sqrt[Sqrt[(2/(.147 Pi ) + Log[1 - x^2]/2)^2 - Log[1 - x^2]/.147] - (2/(.147 Pi) + Log[1 - x^2]/2)]. This works better. $\endgroup$ – nutan Nov 6 '17 at 15:33
  • $\begingroup$ @nutan, two things you need to remember: 1. a Taylor series is usually only accurate for values near its expansion point; and 2. no polynomial has ever had a vertical asymptote. Mull on those two statements, and then think about what you were doing again when you were evaluating at $x\approx 1$ $\endgroup$ – J. M.'s ennui♦ Nov 6 '17 at 19:58
  • 1
    $\begingroup$ This paper may be of interest: ams.org/journals/mcom/1976-30-136/S0025-5718-1976-0421040-7/… Though free access to the body of the paper is given at the AMS archive, the contents of the 62-page microfiche supplement (comprising 88 tables of rational approximants) are not present there. $\endgroup$ – user15994 Nov 6 '17 at 22:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.