QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - Han - 12-11-2013 05:24 PM
Jan. 31, 2018: This program now exists in the HP Prime firmware.
This program takes a decimal value and returns a one of the following expressions that is a "close" rational approximation of the specified decimal value:
\[ \frac{p}{q}, \quad \frac{a}{b}\cdot\sqrt{\frac{p}{q}},
\quad \frac{p}{q}\cdot \pi, \quad e^{\frac{p}{q}}, \quad
\text{ or } \quad \ln \left(\frac{p}{q}\right) \]
Also works for complex numbers and lists of real/complex numbers. The main algorithm basically finds the continued fraction representation of a decimal and the process either self-terminates (in the case of a rational value) or terminates due to reaching the accuracy limit. For example,
\[ \frac{47}{13} = 3 + \frac{1}{1+\frac{1}{1+\frac{1}{1+\frac{1}{1+\frac{1}{2}}}}} \]
To get the continued fraction representation, note that
\[ \frac{47}{13} \approx 3.61538461538 \]
This decimal is converted to a rational expression
\[ \frac{361538461538}{10^{11}} \]
which is converted into the continued fraction as follows:
\[ \frac{361538461538}{10^{11}} = 3 +
\frac{61538461538}{10^{11}} =
3 + \frac{1}{\frac{10^{11}}{61538461538}}
= 3 + \frac{1}{1+ \frac{38461538462}{61538461538}} = \dotsm \]
and so on. While computing the continued fraction, the algorithm simultaneously reduces the partial continued fraction into a rational value of the form \(\frac{p}{q}\). The other forms are merely variations in which the decimal value is either squared, divided by \(\pi\), etc.
Source code below, and also included in attached zip file below.
Code:
// QPI by Han Duong
// ported from QPI 4.3 for the HP48G/GX by Mika Heiskanen & Andre Schoorl
export qpiEXPLN:=100; // max denom for exp(p/q) or ln(p/q)
export qpiMAXINT:=2^20; // largest n allowed for sqrt(n)=a*sqrt(b)
export qpiDIGITS:=10; // controls accuracy (best results at 9 or 10)
qpi_approx();
qpi_asqrtb();
qpi_out();
qpi_outsqrt();
qpi_real();
qpi_complex();
qpi_list();
qpi_root();
qpi_pi();
qpi_ln();
qpi_exp();
EXPORT QPI(r)
BEGIN
case
if TYPE(r)==0 then qpi_real(r); end;
if TYPE(r)==1 then RETURN(r); end;
if TYPE(r)==3 then qpi_complex(r); end;
if TYPE(r)==6 then qpi_list(r); end;
if TYPE(r)==8 then QPI(approx(r)); end;
DEFAULT msgbox("Object type: " + TYPE(r) + " not supported!");
end;
END;
qpi_real(r)
BEGIN
local frac;
if r then
frac:=qpi_approx(r);
if frac(2)<100 then
qpi_out(frac);
else
qpi_root(r,frac);
end;
else
RETURN(0);
end;
END;
qpi_complex(c)
BEGIN
local rpart, ipart;
rpart:=STRING(qpi_real(RE(c)));
ipart:=STRING(qpi_real(abs(IM(c))));
if IM(c)>0 then
expr("'" + rpart + "+" + ipart + "*'"); // bold i symbol
else
expr("'" + rpart + "-" + ipart + "*'"); // bold i symbol
end;
END;
qpi_list(l)
BEGIN
local i,n;
n:=SIZE(l);
for i from 1 to n do
l(i):=QPI(l(i));
end;
RETURN(l);
END;
qpi_root(r,frac)
BEGIN
local frac1;
if r^2<500001 then
frac1:=qpi_approx(r^2);
if r<0 then frac1(1):=-frac1(1); end;
frac1(3):=1;
if (frac1(2)<1000) AND (frac1(2)<=frac(2)) then
if frac1(2)<10 then
qpi_out(frac1);
else
qpi_pi(r,frac1);
end;
else // sqrt denom not smaller
qpi_pi(r,frac);
end;
else // r^2>500000
qpi_pi(r,frac);
end; // end_if r^2<500000
END;
qpi_pi(r,frac)
BEGIN
local frac1;
if abs(r/pi)<101 then
frac1:=qpi_approx(r/pi);
frac1(3):=2;
if (frac1(2)<1000) AND (frac1(2)<=frac(2)) then
if frac1(2)<10 then
qpi_out(frac1);
else
qpi_ln(r,frac1);
end;
else // (r/pi) denom not smaller
qpi_ln(r,frac);
end;
else // abs(r/pi)>100
qpi_ln(r,frac);
end; // end_if abs(r/pi)<101
END;
qpi_ln(r,frac)
BEGIN
local frac1,tmp;
tmp:=e^(r);
if tmp<1001 then
// check for LN(0)
if tmp then
frac1:=qpi_approx(tmp);
else
frac1:=qpi_approx(MINREAL);
end;
frac1(3):=3;
if (frac1(1)*frac1(2)==1) OR (frac1(2)>qpiEXPLN) then
qpi_exp(r,frac);
else
if (frac1(2)<=frac(2)) then
if frac1(2)<10 then
qpi_out(frac1);
else
qpi_exp(r,frac1);
end;
else
qpi_exp(r,frac);
end;
end; // end_if p*q==1 or q>50
else // e^(r)>1000
qpi_exp(r,frac);
end; // end_if e^(r)<1001
END;
qpi_exp(r,frac)
BEGIN
local frac1;
if r<0 then
qpi_out(frac);
else
frac1:=qpi_approx(LN(r));
frac1(3):=4;
if frac1(2)>qpiEXPLN then
qpi_out(frac);
else
if frac1(2)<=frac(2) then
qpi_out(frac1);
else
qpi_out(frac);
end;
end;
end;
END;
// returns frac(t+1) where
// frac:={p/q, sqrt(p/q), p/q*pi, ln(p/q), e^(p/q)}
// and list:={p,q,t}
qpi_out(list)
BEGIN
local s0="(", s1=")'";
if list(3)==1 then
qpi_outsqrt(list);
else
if list(1)<0 then s0:="-" + s0; end;
if list(3) then s1:=")" + s1; end;
case
if list(3)==2 then s1:=")*π'"; end;
if list(3)==3 then s0:="LN(" + s0; end;
if list(3)==4 then s0:="e^(" + s0; end;
end;
s0:="'" + s0;
if list(2)==1 then
expr(s0 + abs(list(1)) + s1);
else
expr(s0 + abs(list(1)) + "/" + list(2) + s1);
end;
end;
END;
qpi_outsqrt(list)
BEGIN
local ab1, ab2;
local s0="'";
if list(1)<0 then s0:=s0+"-"; end;
ab1:=qpi_asqrtb(abs(list(1)));
ab2:=qpi_asqrtb(list(2));
if ab1(1)<>ab2(1) then
if ab2(1)==1 then
s0:=s0 + ab1(1) + "*";
else
s0:=s0 + "(" + ab1(1) + "/" + ab2(1) + ")*";
end;
end;
s0:=s0 + "√(";
if ab2(2)==1 then
s0:=s0 + ab1(2);
else
s0:=s0 + ab1(2) + "/" + ab2(2);
end;
expr(s0+")'");
END;
// returns {a,b} where n=a*sqrt(b)
qpi_asqrtb(n)
BEGIN
local div,quo,rem,num,den,nodd;
if n>qpiMAXINT then RETURN({1,n}); end;
div:=1;
num:=n;
den:=4;
nodd:=3;
repeat
quo:=IP(num/den);
rem:=num MOD den;
if rem==0 then
div:=div*(IP(nodd/2)+1);
num:=quo;
else
nodd:=nodd+2;
den:=den+nodd;
end;
until (quo==0) OR (nodd>qpiMAXINT) end;
RETURN({div,num});
END;
// returns {p,q,0} where r=p/q
qpi_approx(r)
BEGIN
local num,inum,den,iden;
local p0,q0,p1,q1,p2,q2;
local quo,rem;
if NOT(r) then RETURN({0,1,0}); end;
num:=abs(r);
inum:=IP(num);
den:=1;
while num-inum do
num:=num*10;
den:=den*10;
inum:=IP(num);
end;
iden:=den;
rem:=den; den:=num;
p1:=0; p2:=1;
q1:=1; q2:=0;
repeat
p0:=p1; p1:=p2;
q0:=q1; q1:=q2;
num:=den; den:=rem;
quo:=IP(num/den);
rem:=num MOD den;
p2:=quo*p1+p0;
q2:=quo*q1+q0;
until 10^qpiDIGITS*abs(inum*q2-iden*p2)<iden*p2 end;
if (r>0) then
RETURN({p2,q2,0});
else
RETURN({−p2,q2,0});
end;
END;
[attachment=9]
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - ramanvda - 05-31-2014 09:18 PM
Thank you very much!!!
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - salvomic - 02-06-2015 02:57 PM
Han, I like very much your program!
I'd use it almost everywhere...
It, however, should treat also matrices...
As I need this, I use QPI with this ancillary program, after your advice:
Code:
#cas
qpimat(m):=
BEGIN
local s:=dim(m);
m:=mat2list(m);
m:=QPI(m)
m:=list2mat(m,s(2));
return m;
END;
#end
Do you think to include some routine to handle matrices in your original program or it couldn't work in Home mode with "extension"?
Thank you a lot,
Salvo
***
P.S. a part of π, there is a way to rationalize also √π, π^2...?
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - Han - 02-06-2015 07:54 PM
(02-06-2015 02:57 PM)salvomic Wrote: Han, I like very much your program!
I'd use it almost everywhere...
It, however, should treat also matrices...
Do you think to include some routine to handle matrices in your original program or it couldn't work in Home mode with "extension"?
Thank you a lot,
Salvo
***
P.S. a part of π, there is a way to rationalize also √π, π^2...?
If you want it to work in Home view, then a workaround is to convert the matrices into lists of lists. There is no getting around the fact that Home view forces all matrices to have numerical (non-symbolic) values. The drawback is that lists of lists are not displayed to look like matrices.
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - salvomic - 02-06-2015 08:07 PM
(02-06-2015 07:54 PM)Han Wrote: If you want it to work in Home view, then a workaround is to convert the matrices into lists of lists. There is no getting around the fact that Home view forces all matrices to have numerical (non-symbolic) values. The drawback is that lists of lists are not displayed to look like matrices.
yes, for that I prefer works almost always in CAS (I come from HP50g, with no difference between Home and CAS)...
I thought to list of lists, but I prefer real matrices.
For now my little program let me to use your QPI also with matrices (in CAS) and works well. Also I'd prefer that QPI could do the job also with matrices...
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - salvomic - 02-06-2015 08:51 PM
please, Han, could you explain if with QPI would be possible to have rational approximation also for √π (1.77245385091), π^2 (9.86960440109), 2/π, and others "classic" irrational numbers?
Thank you!
Salvo
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - Eddie W. Shore - 02-12-2015 03:36 AM
This is awesome, Han. Great program!
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - compsystems - 11-05-2016 09:49 PM
Request#1: please HP-Prime team adhere this function to the firmware
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - ggauny@live.fr - 11-06-2016 10:46 AM
Hi,
Yes a great good program !
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - salvomic - 11-06-2016 03:02 PM
(11-05-2016 09:49 PM)compsystems Wrote: Why this function is not built-in in the CAS?
I quote!
hoping soon...
Salvo
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - compsystems - 11-06-2016 10:51 PM
Request#2: Please Han, QPI also for symbolic expressions
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - Han - 02-12-2017 04:30 AM
For symbolic expressions you can use this:
Code:
#cas
qpif(f):=
begin
local j,n,g,r;
if (type(f) <> DOM_SYMBOLIC) then return(f); end;
n:=dim(f)+1;
for j from 2 to n do
g:=f[j];
if (type(g) == DOM_FLOAT) then
r:=QPI(g);
f[j]:=r;
end;
if (type(g) == DOM_SYMBOLIC) then
r:=qpif(g);
f[j]:=r;
end;
end;
return(f);
end;
#end
When I have more time I may write a version that handles all known object types. But right now I am focusing on other projects.
[attachment=4476]
RE: QPI: convert decimal to p/q, ln(p/q), p/q*pi, e^(p/q), or sqrt(p/q) - compsystems - 11-19-2017 10:05 AM
Thanks Han
Version QPI_4.4 that includes function with symbolic expressions, available below
Observing the source code, I see that the output could return in various formats,
Request#3:
Please Han adhere a new function called QPIRLNE( EXPR, FORMAT) where FORMAT = 0/1/.../6
case 0: (Default):
qpirlne( expr, 0) -> Expression as pi or root or ln or e
case 1: output only as a expression of QUOTIENT 1:
qpirlne( expr, 1) if it does not find the equivalent to quotient (1), but without
pi, root, ln, e, it returns the same value
case 2: output only as a expression of QUOTIENT 2:
qpirlne( expr, 2) if it does not find the equivalent to quotient (2), but without
pi, root, ln, e, it returns the same value
case 3: output only as a expression of PI:
qpirlne( expr, 3) if it does not find the equivalent to PI, it returns the same value
case 4: output only as a expression of ROOT
qpirlne( expr, 4) if it does not find the equivalent to ROOT, it returns the same value
case 5: output only as a expression of LN
qpirlne( expr, 5) if it does not find the equivalent to LN, it returns the same value
case 6: output only as a expression of EXPR
qpirlne( expr, 6) if it does not find the equivalent to EXP, it returns the same value
PHP Code:
ex#0: qpirlne( (2*π/3)+(3*π/4) , 0) -> 17/12*π
qpirlne( (2*π/3)+(3*π/4) , 1) -> 1137949/255685 // Q1
qpirlne( (2*π/3)+(3*π/4) , 2) -> 4+(115209/255685) // Q2
qpirlne( (2*π/3)+(3*π/4) , 3) -> 17/12*π // PI
qpirlne( (2*π/3)+(3*π/4) , 4) -> (2*π/3)+(3*π/4) // ROOT
qpirlne( (2*π/3)+(3*π/4) , 5) -> (2*π/3)+(3*π/4) // LN
qpirlne( (2*π/3)+(3*π/4) , 6) -> (2*π/3)+(3*π/4) // e
ex#1: qpirlne( LN(3*π)-LN(√(5)), 0) -> LN( (3*π*√(5)/5) ) qpirlne( LN(3*π)-LN(√(5)), 1) -> 55715/38728 // Q1 qpirlne( LN(3*π)-LN(√(5)), 2) -> 1+(16987/38728) // Q2 qpirlne( LN(3*π)-LN(√(5)), 3) -> LN( (3*π*√(5)/5) ) // PI qpirlne( LN(3*π)-LN(√(5)), 4) -> LN( (3*π*√(5)/5) ) // ROOT qpirlne( LN(3*π)-LN(√(5)), 5) -> LN( (3*π*√(5)/5) ) // LN qpirlne( LN(3*π)-LN(√(5)), 6) -> LN(3*π)-LN(√(5)) // e
ex#2: qpirlne( LN((2/5))-LN(√(2)), 0) -> -LN((25/2))/2 qpirlne( LN((2/5))-LN(√(2)), 1) -> -116599/92329 // Q1 qpirlne( LN((2/5))-LN(√(2)), 2) -> -1+(-24270/92329) // Q2 qpirlne( LN((2/5))-LN(√(2)), 3) -> LN((2/5))-LN(√(2)) // PI qpirlne( LN((2/5))-LN(√(2)), 4) -> LN((2/5))-LN(√(2)) // ROOT qpirlne( LN((2/5))-LN(√(2)), 5) -> -LN((25/2))/2 // LN qpirlne( LN((2/5))-LN(√(2)), 6) -> LN((2/5))-LN(√(2)) // e
ex#3:
qpirlne( e^(2*π/(3*√(7))), 0) -> e^((2*π*√(7)/21)) qpirlne( e^(2*π/(3*√(7))), 1) -> 224192/101585 // Q1 qpirlne( e^(2*π/(3*√(7))), 2) -> 1+(21022/101585) // Q2 qpirlne( e^(2*π/(3*√(7))), 3) -> e^(2*π/(3*√(7))) // PI qpirlne( e^(2*π/(3*√(7))), 4) -> e^(2*π/(3*√(7))) // ROOT qpirlne( e^(2*π/(3*√(7))), 5) -> e^(2*π/(3*√(7))) // LN qpirlne( e^(2*π/(3*√(7))), 6) -> e^(2*π/(3*√(7))) // e
ex#4: qpirlne( 7*π/√(90), 0) -> 7*π*√(10)/30
qpirlne( 7*π/√(90), 1) -> 171470/73971
qpirlne( 7*π/√(90), 2) -> 260521/112387
qpirlne( 7*π/√(90), 3) -> 7*π/√(90)
ex#5: qpirlne( 1/(3+i*√(3)), 0) -> (1/4)-i*((√(3)/12)) qpirlne( 1/(3+i*√(3)), 1) -> (1/4)-(1/4)*i*√(1/3) // Q1 qpirlne( 1/(3+i*√(3)), 2) -> (1/4)-(i*37829/262087) // Q1 qpirlne( 1/(3+i*√(3)), 3) -> 1/(3+i*√(3)) // PI qpirlne( 1/(3+i*√(3)), 4) -> (1/4)-i*((√(3)/12)) // ROOT qpirlne( 1/(3+i*√(3)), 5) -> 1/(3+i*√(3)) // LN qpirlne( 1/(3+i*√(3)), 6) -> 1/(3+i*√(3)) // e
ex#6: qpirlne( ACOS((-1/2)), 0) -> 2/3*PI qpirlne( ACOS((-1/2)), 1) -> 138894/66317 // Q1 qpirlne( ACOS((-1/2)), 2) -> 2*(6260/66317) // Q2 qpirlne( ACOS((-1/2)), 3) -> 2/3*PI // PI qpirlne( ACOS((-1/2)), 4) -> ACOS((-1/2) // ROOT qpirlne( ACOS((-1/2)), 5) -> ACOS((-1/2) // LN qpirlne( ACOS((-1/2)), 6) -> ACOS((-1/2) // e
ex#7: qpirlne( COS((3*π/4)), 0) -> -√(-2)/2 qpirlne( COS((3*π/4)), 1) -> -195025/275807 // Q1 qpirlne( COS((3*π/4)), 2) -> -195025/275807 // Q2 qpirlne( COS((3*π/4)), 3) -> COS((3*π/4)) // PI qpirlne( COS((3*π/4)), 4) -> -√(1/2) // ROOT qpirlne( COS((3*π/4)), 5) -> COS((3*π/4)) // LN qpirlne( COS((3*π/4)), 6) -> COS((3*π/4)) // e
ex#8: qpirlne( COS(π/12), 0) -> (√(3)+1)*(√(2)/4) qpirlne( COS(π/12), 1) -> 129209/133767 // Q1 qpirlne( COS(π/12), 2) -> 272847/282472 // Q2 qpirlne( COS(π/12), 3) -> COS((3*π/4)) // PI qpirlne( COS(π/12), 4) -> (√(3)+1)*(√(2)/4) // ROOT qpirlne( COS(π/12), 5) -> COS(π/12) // LN qpirlne( COS(π/12), 6) -> COS(π/12) // e
ex#9: qpirlne( SIN(π/10), 0) -> (-1+√((5)))/4 qpirlne( SIN(π/10), 1) -> 98209/317811 // Q1 qpirlne( SIN(π/10), 2) -> 98209/317811 // Q2 qpirlne( SIN(π/10), 3) -> SIN(π/10) // PI qpirlne( SIN(π/10), 4) -> (-1+√((5)))/4 // ROOT qpirlne( SIN(π/10), 5) -> SIN(π/10) // LN qpirlne( SIN(π/10), 6) -> SIN(π/10) // e
ex#10: qpirlne( SIN(π/8), 0) -> √(2-√(2))/2 qpirlne( SIN(π/8), 1) -> 69237/180925 // Q1 qpirlne( SIN(π/8), 2) -> 69237/180925 // Q2 qpirlne( SIN(π/8), 3) -> SIN(π/8) // PI qpirlne( SIN(π/8), 4) -> √(2-√(2))/2 // ROOT qpirlne( SIN(π/8), 5) -> SIN(π/8) // LN qpirlne( SIN(π/8), 6) -> SIN(π/8) // e
ex#11: qpirlne( COS(π/5), 0) -> (1+(√(5)))/4 qpirlne( COS(π/5), 1) -> 98209/121393 // Q1 qpirlne( COS(π/5), 2) -> 317811/392836 // Q2 qpirlne( COS(π/5), 3) -> COS(π/5) // PI qpirlne( COS(π/5), 4) -> (1+(√(5)))/4 // ROOT qpirlne( COS(π/5), 5) -> COS(π/5) // LN qpirlne( COS(π/5),, 6) -> COS(π/5) // e
|