ASKSAGE: Sage Q&A Forum - Individual question feedhttp://ask.sagemath.org/questions/Q&A Forum for SageenCopyright Sage, 2010. Some rights reserved under creative commons license.Wed, 18 May 2011 17:22:31 -0500evaluating inverse erfhttp://ask.sagemath.org/question/8118/evaluating-inverse-erf/Hi - The code below results in an expression, y, involving the inverse erf function i.e. exp(sqrt(2) inverse_erf(1/2)) which is approximately 1.96 according to another mathematical package. I tried to get a numerical value for y using y.n() but that crashes. Can anyone please advise how to evaluate such symbolic expressions?
var('z,t')
PDF = exp(-x^2 /2)/sqrt(2*pi)
integralExpr = integrate(PDF,x,z,oo).subs(z==log(t))
y = solve(integralExpr==z,t)[0].rhs().subs(z==1/4)
yWed, 18 May 2011 03:22:39 -0500http://ask.sagemath.org/question/8118/evaluating-inverse-erf/Answer by kcrisman for <p>Hi - The code below results in an expression, y, involving the inverse erf function i.e. exp(sqrt(2) inverse_erf(1/2)) which is approximately 1.96 according to another mathematical package. I tried to get a numerical value for y using y.n() but that crashes. Can anyone please advise how to evaluate such symbolic expressions? </p>
<pre><code>var('z,t')
PDF = exp(-x^2 /2)/sqrt(2*pi)
integralExpr = integrate(PDF,x,z,oo).subs(z==log(t))
y = solve(integralExpr==z,t)[0].rhs().subs(z==1/4)
y
</code></pre>
http://ask.sagemath.org/question/8118/evaluating-inverse-erf/?answer=12359#post-id-12359This is something relatively new to Sage, via Maxima. See [this Maxima devel post](http://www.ma.utexas.edu/pipermail/maxima/2009/017196.html) for details.
What this means is that we have yet to implement this conversion in Sage - not surprising, since it doesn't seem to be in Ginac. (In fact, erf doesn't even seem to be in Ginac!)
However, it *is* in mpmath.
sage: e^(sqrt(2)*RR(mpmath.erfinv(1/2))).n()
1.96303108415826
Unfortunately, the typical beginning user might not know the `RR` trick. There are two things we need to do for this.
- Have a symbolic inverse erf function, which we do not yet have.
- Make a `._sage_()` method for mpmath numbers, or make it easier to figure out how to get them symbolic. There is also no `.n()` method. I shouldn't have to to this if I forget `RR`!
sage: e^(sqrt(2)*sage.libs.mpmath.all.mpmath_to_sage(mpmath.erfinv(1/2),53)).n()
1.96303108415826
----------
This reminds me that askbot has still not fixed the issue of formatting after bullets (which, to be fair, was upstream from them...)Wed, 18 May 2011 07:52:35 -0500http://ask.sagemath.org/question/8118/evaluating-inverse-erf/?answer=12359#post-id-12359Comment by kcrisman for <p>This is something relatively new to Sage, via Maxima. See <a href="http://www.ma.utexas.edu/pipermail/maxima/2009/017196.html">this Maxima devel post</a> for details.</p>
<p>What this means is that we have yet to implement this conversion in Sage - not surprising, since it doesn't seem to be in Ginac. (In fact, erf doesn't even seem to be in Ginac!)</p>
<p>However, it <em>is</em> in mpmath.</p>
<pre><code>sage: e^(sqrt(2)*RR(mpmath.erfinv(1/2))).n()
1.96303108415826
</code></pre>
<p>Unfortunately, the typical beginning user might not know the <code>RR</code> trick. There are two things we need to do for this.</p>
<ul>
<li>Have a symbolic inverse erf function, which we do not yet have.</li>
<li><p>Make a <code>._sage_()</code> method for mpmath numbers, or make it easier to figure out how to get them symbolic. There is also no <code>.n()</code> method. I shouldn't have to to this if I forget <code>RR</code>!</p>
<p>sage: e^(sqrt(2)*sage.libs.mpmath.all.mpmath_to_sage(mpmath.erfinv(1/2),53)).n()</p>
<p>1.96303108415826</p></li>
</ul>
<hr/>
<p>This reminds me that askbot has still not fixed the issue of formatting after bullets (which, to be fair, was upstream from them...)</p>
http://ask.sagemath.org/question/8118/evaluating-inverse-erf/?comment=21719#post-id-21719Implementing this natively in Sage is now http://trac.sagemath.org/sage_trac/ticket/11349.Wed, 18 May 2011 08:00:18 -0500http://ask.sagemath.org/question/8118/evaluating-inverse-erf/?comment=21719#post-id-21719Comment by RossK for <p>This is something relatively new to Sage, via Maxima. See <a href="http://www.ma.utexas.edu/pipermail/maxima/2009/017196.html">this Maxima devel post</a> for details.</p>
<p>What this means is that we have yet to implement this conversion in Sage - not surprising, since it doesn't seem to be in Ginac. (In fact, erf doesn't even seem to be in Ginac!)</p>
<p>However, it <em>is</em> in mpmath.</p>
<pre><code>sage: e^(sqrt(2)*RR(mpmath.erfinv(1/2))).n()
1.96303108415826
</code></pre>
<p>Unfortunately, the typical beginning user might not know the <code>RR</code> trick. There are two things we need to do for this.</p>
<ul>
<li>Have a symbolic inverse erf function, which we do not yet have.</li>
<li><p>Make a <code>._sage_()</code> method for mpmath numbers, or make it easier to figure out how to get them symbolic. There is also no <code>.n()</code> method. I shouldn't have to to this if I forget <code>RR</code>!</p>
<p>sage: e^(sqrt(2)*sage.libs.mpmath.all.mpmath_to_sage(mpmath.erfinv(1/2),53)).n()</p>
<p>1.96303108415826</p></li>
</ul>
<hr/>
<p>This reminds me that askbot has still not fixed the issue of formatting after bullets (which, to be fair, was upstream from them...)</p>
http://ask.sagemath.org/question/8118/evaluating-inverse-erf/?comment=21716#post-id-21716Thanks for the answer (I have something to work with now) and for creating the ticket!Wed, 18 May 2011 17:22:31 -0500http://ask.sagemath.org/question/8118/evaluating-inverse-erf/?comment=21716#post-id-21716Answer by RossK for <p>Hi - The code below results in an expression, y, involving the inverse erf function i.e. exp(sqrt(2) inverse_erf(1/2)) which is approximately 1.96 according to another mathematical package. I tried to get a numerical value for y using y.n() but that crashes. Can anyone please advise how to evaluate such symbolic expressions? </p>
<pre><code>var('z,t')
PDF = exp(-x^2 /2)/sqrt(2*pi)
integralExpr = integrate(PDF,x,z,oo).subs(z==log(t))
y = solve(integralExpr==z,t)[0].rhs().subs(z==1/4)
y
</code></pre>
http://ask.sagemath.org/question/8118/evaluating-inverse-erf/?answer=12362#post-id-12362Thanks for the answer (I have something to work with now) and creating the ticket!Wed, 18 May 2011 17:21:42 -0500http://ask.sagemath.org/question/8118/evaluating-inverse-erf/?answer=12362#post-id-12362