HP Forums
HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Printable Version

+- HP Forums (https://www.hpmuseum.org/forum)
+-- Forum: HP Calculators (and very old HP Computers) (/forum-3.html)
+--- Forum: HP Prime (/forum-5.html)
+--- Thread: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) (/thread-20844.html)



HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Gil - 11-15-2023 04:00 PM

Try for instance on the HP PRIME approximate(erfc(5*i))

ie with no real part in the complex number input (argument).

It returns -82982738 79.76*i (a complex number with no real part).
According to Wolfram Alpha, the result should be
1 - 8.2982738 806 7680 ×10^9 i
(with the number 1 as the real part of the resulting complex number).

Fazit
1) As no real part of the HP PRIME result is shown, that part has to be considered being equal to 0, and the PRIME result seems therefore clearly to be wrong here (a bug?).
2) As a basic program on the HP50G allows to find (1.,-82982738 80.61), we see that the imaginary part of the HP PRIME is not so well "calibrated", letting suppose that the internal algorithm of the HP PRIME for erfc function could or should be improved.

Logically, on HP PRIME, approximate(erf(5*i)) returns
1 -82982738 79.76*i (a complex number with the number 1 as its real part, instead of 0).

Strange.

Or did I miss something?

Regards

Gil


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Albert Chan - 11-15-2023 04:53 PM

(11-15-2023 04:00 PM)Gil Wrote:  Logically, on HP PRIME, approximate(erf(5*i)) returns
1 -82982738 79.76*i (a complex number with the number 1 as its real part, instead of 0).

Strange.

erf(x) = 2/sqrt(pi) * x * hyp1f1(1/2, 3/2, -x*x)

erf(x)/x = 2/sqrt(pi) * hyp1f1(1/2, 3/2, -x*x)

If x is purely real or imaginery, RHS is real --> LHS is real too.

Without computation, erf(5i) is purely imaginery.


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Albert Chan - 11-15-2023 05:16 PM

(11-15-2023 04:53 PM)Albert Chan Wrote:  Without computation, erf(5i) is purely imaginery.

I would not consider this a bug though.
Difference is very tiny, especially compared in polar form.

p2> from mpmath import *
p2> x = erf(5j)
p2> x
mpc(real='0.0', imag='8298273880.6768036')
p2> polar(x)
(mpf('8298273880.6768036'), mpf('1.5707963267948966'))
p2> polar(x + 1)
(mpf('8298273880.6768036'), mpf('1.5707963266743896'))


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Brian Zilli - 11-15-2023 06:05 PM

The cutoff appears to be 3; that is, for any a>3, erf(a*i) will return a real part of 1. Near 3, the relative error is orders of magnitude larger.

It also leads to wonky results when computing purely real integrals in the CAS (and numerical Home in the 2020 firmware of the Mac emulator but not on a physical G2 with the latest firmware). e.g., int(e^(x^2),x,0,3.001) returns 1452.55199576+0.886226925453*I, when it should be approximately 1452.67.


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Brian Zilli - 11-15-2023 08:04 PM

Also interesting: try erf(3.16845945394*i) or int(e^(x^2),x,0,3.16845945394) in the CAS.


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - parisse - 11-15-2023 08:23 PM

I confirm a small bug for erfc if |z| is between 3 and 6.5 and re(z) is small (missed a 1 in real part). In this area, erf is computed with a continued fraction expansion
2*exp(z^2)*int(exp(-t^2),t=z..inf)=1/(z+1/2/(z+1/(z+3/2/(z+...))))
There is a warning of low accuracy if re(z) is small : I must indeed stop the fraction expansion and since we do not have multiprecision floats there are some wrong digits.


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Albert Chan - 11-15-2023 08:26 PM

Cas> a := 3.16845945394*i
Cas> int(exp(-x*x), x=0..a)       → 0.886226925453-51940029439.9*i
Cas> Ans * 2/sqrt(pi)                → 1.-58608047158.3*i
Cas> erf(a)                               → 1.-58608047158.3*i 

Try numerical integration instead (double exponential method, quad6)

Cas> quad(x -> exp(-x*x), 0, a)       → [3833.13077361*i, 1.16978024776e−12]
Cas> Ans[1] * 2/sqrt(pi)                   → 4325.2249097*i == erf(3.16845945394*i)


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - parisse - 11-15-2023 08:32 PM

The continued fraction expansion "low accuracy" is probably full wrong for some values of z near the imaginary axis...


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - parisse - 11-16-2023 01:22 PM

I have disabled the continued fraction, it should now work correctly (it will just be a little bit slower if |z| is between 4 and 6.5).


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Gil - 11-16-2023 01:33 PM

I have the Android version for mobile phone.

Should be I desinstall it and install it again?

Or wait for an automatic update?

Right now erf (5i)returns 1+b×i.


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Gil - 11-17-2023 12:50 AM

An interesting approach for erf(z) is the one I found in
https://math.stackexchange.com/questions/712434/erfaib-error-function-separate-into-real-and-imaginary-part

for which the resulting calculations (for a number of about 30 summations) seem to give better results (in the last digit) in my HP50G than the built-in function in HP PRIME.


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - StephenG1CMZ - 11-17-2023 09:11 AM

(11-16-2023 01:33 PM)Gil Wrote:  I have the Android version for mobile phone.

Should be I desinstall it and install it again?

Or wait for an automatic update?

Right now erf (5i)returns 1+b×i.

I would not recommend uninstalling and reinstalling because Google Play does not reliably list the HP Prime Pro app.
So at best, you might get back what you had. At worst, there is a risk that you might be unable to reinstall.
https://www.hpmuseum.org/forum/thread-20078.html?highlight=hit+and+miss
Eventually, a new version of the App might be published, and that might include the updates described.
KlassKuperus said in September that an update was being worked on.


RE: HP Prime bug for approximate(erfc(z)), z as (0+b*i), or erf(b*i) - Albert Chan - 04-07-2024 12:06 PM

(11-15-2023 08:23 PM)parisse Wrote:  There is a warning of low accuracy if re(z) is small : I must indeed stop the fraction expansion and
since we do not have multiprecision floats there are some wrong digits.

Throwing more precisions may not help.
Problem is more to do with how deep "fraction expansion" default we pick.

Let n = "last" numerator --> erfc(z, n=1) = exp(-z*z)/√(pi) / (z + (1/2)/(z + 1/z))

Let z = 3.16845945394*I, precisions = 1000 digits, results rounded to integer.

erfc(z, n=0.5) = -4292 I
erfc(z, n=1.0) = -4317 I
erfc(z, n=1.5) = -4323 I
erfc(z, n=2.0) = -4324 I
erfc(z, n=2.5) = -4325 I
erfc(z, n=3.0) = -4325 I
erfc(z, n=3.5) = -4326 I
erfc(z, n=4.0) = -4311 I
erfc(z, n=4.5) = -4324 I
erfc(z, n=5.0) = -4325 I
...
erfc(z, n=19.0) = -4324 I
erfc(z, n=19.5) = -4326 I
erfc(z, n=20.0) = 58957691655 I
erfc(z, n=20.5) = -4325 I
erfc(z, n=21.0) = -4326 I
...

For z=3.16845945394*I, n=20 is really bad.
But changing n will not help. It would just shift problematic z somewhere else.

---

More precisions can't help because CF formula denominator is close to 0
But, erfc(z) is continuous ... something is going to break.

1 / ±0 = ±∞

Cas> erfc(3.16845945394*i)      → +58608047158.3*i
Cas> erfc(3.16845945395*i)      → −14620291004.7*i

CF Denomaintor = x * p(x²)/q(x²), rational polynomial with all positive terms
--> erfc CF formula for non-zero real argument will not get divide-by-0 issue.

(04-06-2024 11:04 PM)Albert Chan Wrote:  denominator may get to 0 only if re(z)=0