+- 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: Beta and Airy Functions (/thread-4216.html)
Beta and Airy Functions - douganc - 06-24-2015 11:20 AM
Don't seem to be fully implemented on HP Prime.
XCAS:
Input: Beta(x,y)
Output: Gamma x * Gamma y / Gamma (x+y)
Prime: ??????
XCAS:
Input: Beta(n,2)
Output: Gamma n / Gamma (n+2)
Prime: ??????
XCAS:
Input: Airy_Ai(1)
Output: 0.135292416313
Prime: ????????
XCAS:
Input: Airy_Bi(1)
Output: 1.20742359495
Prime: ?????????
Latest Software Rev.
RE: Beta and Airy Functions - Tim Wessman - 06-24-2015 02:39 PM
All of these seem to be working fine for me and match xcas output.
Airy cant' be implemented on our end due to it using a numerical library we can't use at the moment - hence why it wasn't included.
RE: Beta and Airy Functions - Gerald H - 06-24-2015 03:27 PM
Anyone know an efficient way of caculating the Airy function?
RE: Beta and Airy Functions - DrD - 06-24-2015 03:33 PM
A burrito rich diet?
RE: Beta and Airy Functions - parisse - 06-24-2015 06:48 PM
http://people.math.sfu.ca/~cbm/aands/page_446.htm
RE: Beta and Airy Functions - roadrunner - 06-27-2015 11:19 AM
(06-24-2015 06:48 PM)parisse Wrote: http://people.math.sfu.ca/~cbm/aands/page_446.htm
That method becomes unstable outside a range for z of about -10 to 10. I believe that is because the calculator is trying to divide one huge number by another huge number to get a small number.
[attachment=2227]
How can one calculate the airy function outside that range?
Code:
local f(z)
begin
local t,t2,k,p,p2;
t:=1;
p2:=1;
for k from 1 to 80 do
p:=3*k;
p2:=p2*(3*k-2);
t:=t+p2*z^(p)/(p)!;
if t == t2 then return t; end;
t2:=t;
end;
return t;
end;
local g(z)
begin
local t,t2,k,p,p2;
t:=z;
t2:=0;
p2:=1;
for k from 1 to 80 do
p:=3*k+1;
p2:=p2*(3*k-1);
t:=t+p2*z^(p)/(p)!;
if t == t2 then return t; end;
t2:=t;
end;
return t;
end;
export ai(z)
begin
0.355028053887817*f(z)-
0.258819403792807*g(z);
end;
export bi(z)
begin
(0.355028053887817*f(z)+
0.258819403792807*g(z))*√3;
end;RE: Beta and Airy Functions - Gerald H - 06-27-2015 12:15 PM
Thank you, parisse & more so roadrunner, here my translation of your programme for the 48G.
Code:
Ai
« DUP f .355028053888
* SWAP g
.258819403793 * -
»
Bi
« DUP f .355028053888
* SWAP g
.258819403793 * + 3.
√ *
»
g
« DUP 0. DUP 1. Z T
T2 P P2
« 1. 80.
FOR K K 3. * 1. +
'P' STO P2 3. K * 1.
- * 'P2' STO T P2 Z P
^ * P ! / + 'T' STO T
T2 ==
IF
THEN 81. 'K'
STO
ELSE T 'T2' STO
END
NEXT T
»
»
f
« 1. 4. NDUPN DROP
Z T T2 P P2
« 1. 80.
FOR K K 3. * 'P'
STO P2 3. K * 2. - *
'P2' STO T P2 Z P ^ *
P ! / + 'T' STO T T2
==
IF
THEN 81. 'K'
STO
ELSE T 'T2' STO
END
NEXT T
»
»RE: Beta and Airy Functions - parisse - 06-27-2015 12:41 PM
For |z| large, you can use asymptotic series expansion near infinity.
http://people.math.sfu.ca/~cbm/aands/page_448.htm
RE: Beta and Airy Functions - Gerald H - 06-27-2015 05:25 PM
& here's a sys version for the 49G.
Code:
Ai
::
CK1&Dispatch
BINT1
::
DUP
ID f
% .355028053888
%*
SWAP
ID g
% .258819403793
%*
%-
;
;
Bi
::
CK1&Dispatch
BINT1
::
DUP
ID f
% .355028053888
%*
SWAP
ID g
% .258819403793
%*
%+
%3
%SQRT
%*
;
;
f
::
%1
BINT4
NDUPN
DROP
'
NULLLAM
BINT5
NDUPN
DOBIND
BINT81
ONE_DO
INDEX@
UNCOERCE
DUP
%3
%*
2PUTLAM
1GETLAMSWAP_
%3
%*
%2
%-
%*
1PUTLAM
4GETLAM
1GETLAM
5GETLAM
2GETLAM
%^
%*
2GETLAM
%NFACT
%/
%+
DUP4PUTLAM
3GETLAM
%=
ITE
ExitAtLOOP
::
4GETLAM
3PUTLAM
;
LOOP
4GETLAM
ABND
;
g
::
DUP
%0
DUP
%1
'
NULLLAM
BINT5
NDUPN
DOBIND
BINT81
ONE_DO
INDEX@
UNCOERCE
DUP
%3
%*
%1+
2PUTLAM
1GETLAMSWAP_
%3
%*
%1-
%*
1PUTLAM
4GETLAM
1GETLAM
5GETLAM
2GETLAM
%^
%*
2GETLAM
%NFACT
%/
%+
DUP4PUTLAM
3GETLAM
%=
ITE
ExitAtLOOP
::
4GETLAM
3PUTLAM
;
LOOP
4GETLAM
ABND
;RE: Beta and Airy Functions - roadrunner - 06-28-2015 12:47 AM
(06-27-2015 12:41 PM)parisse Wrote: For |z| large, you can use asymptotic series expansion near infinity.
http://people.math.sfu.ca/~cbm/aands/page_448.htm
That's much better.
[attachment=2229]
Now it matches this web site: http://keisan.casio.com/exec/system/1180573401 pretty close. I only checked out to +/-100.
Code:
export c(k)
begin
if k == 0 then return 1; end;
return Gamma(3*k+1/2)/54^k/k!/Gamma(k+1/2);
end;
export d(k)
begin
return c(k)*(6*k+1)/(1-6*k);
end;
local f(z)
begin
local t,t2,k,p,p2;
t:=1;
t2:=0;
p2:=1;
for k from 1 to 80 do
p:=3*k;
p2:=p2*(3*k-2);
t:=t+p2*z^(p)/(p)!;
if t == t2 then return t; end;
t2:=t;
end;
return t;
end;
local g(z)
begin
local t,t2,k,p,p2;
t:=z;
t2:=0;
p2:=1;
for k from 1 to 80 do
p:=3*k+1;
p2:=p2*(3*k-1);
t:=t+p2*z^(p)/(p)!;
if t == t2 then return t; end;
t2:=t;
end;
return t;
end;
local aismallz(z)
begin
return
0.355028053887817*f(z)-
0.258819403792807*g(z);
end;
local aibigz(z)
begin
local s;
if z ≥ 0 then
s:=2/3*z^(3/2);
return π^(−1/2)/2*z^(−1/4)*e^(−s)*
sum((−1)^K*c(K)*s^(−K),K,0,20);
end;
z:=ABS(z);
s:=2/3*z^(3/2);
return π^(−1/2)*z^(−1/4)*(SIN(s+π/4)*
sum((−1)^K*c(2*K)*s^(−2*K),K,0,20)-
COS(s+π/4)*
sum((−1)^K*c(2*K+1)*s^(−2*K-1),K,0,20));
end;
local bismallz(z)
begin
return (0.355028053887817*f(z)+
0.258819403792807*g(z))*√3;
end;
local bibigz(z)
begin
local s;
if z ≥ 0 then
s:=2/3*z^(3/2);
return π^(−1/2)*z^(−1/4)*e^(s)*
sum(c(K)*s^(−K),K,0,20);
end;
z:=ABS(z);
s:=2/3*z^(3/2);
return π^(−1/2)*z^(−1/4)*(COS(s+π/4)*
sum((−1)^K*c(2*K)*s^(−2*K),K,0,20)+
SIN(s+π/4)*
sum((−1)^K*c(2*K+1)*s^(−2*K-1),K,0,20));
end;
export ai(z)
begin
HAngle:=0;
if ABS(z) < 10 then
return aismallz(z);
end;
return aibigz(z);
end;
export bi(z)
begin
HAngle:=0;
if ABS(z) < 10 then
return bismallz(z);
end;
return bibigz(z);
end;RE: Beta and Airy Functions - douganc - 07-20-2015 08:10 PM
Thank You !!
A very useful program.
Request for you: Could you please put the program in the HP Prime Software library section.
So that other people can use the functions as well.
Thank You
Colm
RE: Beta and Airy Functions - roadrunner - 07-20-2015 08:55 PM
Colm,
I had some issues with the version that I posted in this thread so I posted a newer version that clears up those issues and runs faster to boot.
Road
PS. Thanks to Mr. Parisse for teaching me asymptotic series expansion near infinity.