48/49: Best Fraction Challenge
Joseph K. Horn
(An Extremely Difficult Challenge Disguised as a Mini-Challenge)
DEFINITION OF "BEST FRACTION": What fraction is equal to 0.5? There
are, of course, infinitely many correct answers to that question. If
we limit the answers to ratios of positive integers, there are still
infinitely many answers: 1/2, 2/4, 3/6 ... 1234/2468 ... all evaluate
to 0.5 on a calculator. But the *smallest* integers (more precisely,
the least non-negative integers), with a ratio of 0.5, are 1 and 2.
In this sense, the "best fraction" for 0.5 is 1/2.
CHALLENGE: Write a program (none such has ever been published!) that
inputs any positive decimal number, and outputs the "best fraction"
for it on the HP48 or HP49, that is, the smallest integers which, when
divided on the calculator, yield exactly the same number as the input,
in 5 seconds or less for any input, in 100% User RPL.
Hundreds of *similar* programs have been written over the years. In
fact, the ->Q function is similar. But NONE of these meet the "Best
Fraction Challenge". All of them yield some answers not exactly equal
to the input, or some answers that aren't the "best".
WARNING! I have posted this challenge here (in various forms) over the
past 13 years, and almost every time some people responded with
haughty disdain, crowing that it's a trivial and already-solved
problem... but of course, when taken to task, they all failed to post
a program that would have proved them right. If you really think that
the above is as easy as it sounds, then to save face *please* code it
up and test it on the following examples before bragging here about
how easily you could do it if you only had the time or some other
transparent excuse. Thank you.
EXAMPLES of inputs and CORRECT outputs (NOT ->Q outputs!) for the
HP48/49:
Input: 2 SIN (in RAD mode)
Output: 1029101 / 1131754
Input: 6.0030000004 (the strange number on the HP38G case)
Output: 14823406 / 2469333
Input: 0.999999999975 (the "worst case" ->Q number)
Output: 39215686274 / 39215686275
Input: 0.9498580222 (my telephone number!)
Output: 1306259 / 1375215
Input: 50000.00005
Output: 999050001 / 19981
Input: 0.500000000001
Output: 166666666668 / 333333333335
Input: 0.69999937
Output: 777777/1111111
Input: 0.903493239856
Output: 528118/584529
Of course, the program must also be able to handle ordinary inputs,
0.5 etc.
Yes, I've done it. I showed it off at HHC 2003... but I didn't give
copies of the program to anybody, for reasons to be revealed soon. If
ANYBODY can meet or beat my program, it will be of monumental
importance to me and to all users of calculator fractions.
Good luck! Happy Programming!
-Joe-
Joseph K. Horn
p...@holyjoe.us
James
> The "Best Fraction Challenge"
> (An Extremely Difficult Challenge Disguised as a Mini-Challenge)
>
> CHALLENGE: Write a program (none such has ever been published!) that
> inputs any positive decimal number, and outputs the "best fraction" for
> it on the HP48 or HP49, that is, the smallest integers which, when
> divided on the calculator, yield exactly the same number as the input,
> in 5 seconds or less for any input, in 100% User RPL.
Pardon my ignorance, but can't this already be accomplished using the XQ
function? I worked through some of your example problems with it, and they
seemed to be properly reduced.
James
FX-502P Geek
FX-502P back in the 80s (except I preferred the output in vulgar fractions),
but by simple inspection of your challenge you've got the wrong answer for
the 5th example. The fractional part of 50000.0005 is clearly 1/20000 and so
your correct answer should be 1000000001/20000 as given by ->Q of a HP48.
Incidentally the program that I wrote for the FX-502P is available at
www.fx502p.pwp.blueyonder.co.uk and gives the correct answer of
50000+1/20000. I would not pretend however that the program would be up to
the challenge in general. I would of course like to see your HP48 program...
especially when it is fixed.
Peter.
"Joseph K. Horn" <joe...@holyjoe.net> wrote in message
news:87233f9e.03101...@posting.google.com...
Tony Hutchins
Hi FX-502P Geek ! <use_th...@mysite.co.uk>
01h56m ago, on Sat, 11.10.03 at 9:52 p.m. +0100, you wrote
in message ID <H0_hb.1364$ny5...@news-binary.blueyonder.co.uk> :
> I do not deny that the problem is very difficult - I
> tried it on the Casio FX-502P back in the 80s (except I
> preferred the output in vulgar fractions), but by simple
> inspection of your challenge you've got the wrong answer
> for the 5th example. The fractional part of 50000.0005
> is clearly 1/20000 and so your correct answer should be
> 1000000001/20000 as given by ->Q of a HP48.
That is thw whole point of the challenge - here is Joe's answer:
Input: 50000.00005
Output: 999050001 / 19981
And 19981 is LESS than 20000.
Joe's answer is surprisingly correct :)
Not only "better" than the 20000 denominator, but he claims it
to be the "best" - quite sensational. I can't wait to see
where his program will be published. It must be a totally new
approach.
> Incidentally the program that I wrote for the FX-502P is
> available at www.fx502p.pwp.blueyonder.co.uk and gives
> the correct answer of 50000+1/20000. I would not pretend
> however that the program would be up to the challenge in
> general. I would of course like to see your HP48 program...
> especially when it is fixed.
I suspect that it won't need fixing<G>. Although a real
counter-example would no doubt be welcomed.
- Tony
BlackLotus
"Tony Hutchins" <jus...@csi.com> a écrit dans le message de
news:1078940...@news.cis.dfn.de...
> -=[ Sun, 12.10.03 11:48 a.m. +1300 (NZDT) ]=-
>
> Hi FX-502P Geek ! <use_th...@mysite.co.uk>
>
> 01h56m ago, on Sat, 11.10.03 at 9:52 p.m. +0100, you wrote
> in message ID <H0_hb.1364$ny5...@news-binary.blueyonder.co.uk> :
>
> > I do not deny that the problem is very difficult - I
> > tried it on the Casio FX-502P back in the 80s (except I
> > preferred the output in vulgar fractions), but by simple
> > inspection of your challenge you've got the wrong answer
> > for the 5th example. The fractional part of 50000.0005
> > is clearly 1/20000 and so your correct answer should be
> > 1000000001/20000 as given by ->Q of a HP48.
>
> That is thw whole point of the challenge - here is Joe's answer:
>
> Input: 50000.00005
> Output: 999050001 / 19981
>
> And 19981 is LESS than 20000.
> Joe's answer is surprisingly correct :)
It is because of a tiny approximation error. I checked it on a Math soft on
the PC, which performs "exact" calculation, and substracting the two
fractions shows an error of 19/399620000. So, no, Joe's answer's
mathmatically wrong, it's just a twist from the HP48/49 that makes it
"right" (4E-8 is neglected in regard to the number of digits and their order
in 50000.00005 as "too small to be significant")
I'm still curious as to how you achieved that, though.
BlackLotus
"Joseph K. Horn" <joe...@holyjoe.net> a écrit dans le message de
news:87233f9e.03101...@posting.google.com...
> The "Best Fraction Challenge"
Here's how I'd do it, assuming I had a previous prime algorithm: (I'll call
it PP, I suppose Joe has that somewhere for the 48, but I can't find it, and
I think Erable has it in built-in for the 49)
Input: Positive integer
Output:
<< -> X
<< X ->Q OBJ-> DROP DROP SWAP DROP
'Q' STO Q 'R' STO
WHILE Q X * 0 RND Q / X ==
REPEAT Q DUP 'R' STO PP 'Q' STO
END
"'"
Q R * 0 RND + "/" + R + "'" + OBJ->
>>
>>
Needless to say (and as I stated in another post), this should give (I can't
test it) the best fraction that will return a correct result in the HP48/49
(and I insist, it only works because of how the ints are in the HP48/49)
BlackLotus
Output: Fraction
Tony Hutchins
Hi BlackLotus ! <black...@nospam.me>
34m ago, on Sun, 12.10.03 at 02:02 a.m. +0200, you wrote
in message ID <bma5m5$76e$1...@news-reader4.wanadoo.fr> :
news client =Microsoft Outlook Express 6.00.2800.1158
previous messages in thread =3 + 1 for your message
tot ref chars=133
chars per id =44.3333333333333
Top Thread id=<87233f9e.03101...@posting.google.com>
Organisation =Wanadoo, l'internet avec France Telecom
--
Tony Hutchins
Wellington
New Zealand
#113 Fate often saves a warrior when his courage
endures. Beowolf
Joseph K. Horn
> Pardon my ignorance, but can't this already be
> accomplished using the XQ function?
No. XQ has the same faults as ->Q, namely, it gets answers that are
either (a) not the best, (b) not equal to the input, or (c) both of
the above.
> I worked through some of your example problems
> with it, and they seemed to be properly reduced.
Umm... No. Look again. XQ gets horribly wrong answers for all of
the examples except 50000.00005 which it can't even handle ("XQ Error:
Overflow"!!!). Compare the answers you get with the answers that I
listed with each example. None match, and some are WAY off. ->Q is
better than XQ, but even that fails to get the right answer to ANY of
my examples.
Peter wrote:
> ... by simple inspection of your challenge you've
> got the wrong answer for the 5th example. The
> fractional part of 50000.0005 is clearly 1/20000
> and so your correct answer should be
> 1000000001/20000 as given by ->Q of a HP48.
Ummm... No. Look again. I want the BEST FRACTION (not the obvious
one), with BEST FRACTION being defined as "the ratio of the smallest
two integers which, when divided ON AN HP48 or HP49, yield exactly the
input."
Yes, 50000.00005 can be represented by 1000000001/20000, but please
look closer at the answer I gave. It has a SMALLER numerator and
denominator, AND it evaluates to the same thing as the input ON AN
HP48 or HP49. Therefore it is a "better" answer. And since there are
no smaller integers that work, it's the "best fraction" for the input.
Recap:
1000000001 / 20000 = 50000.00005 on an HP48/49, but
999050001 / 19981 = 50000.00005 on an HP48/49 too.
So the latter is better.
Nothing is better than it, so it's the best fraction.
> Incidentally the program that I wrote for the FX-502P
> is available at www.fx502p.pwp.blueyonder.co.uk and
> gives the correct answer of 50000+1/20000.
The "correct" answer depends on the machine, just like "the correct
sine of 2" depends on the machine. If you tell me how many digits the
FX-502P calculates with, and whether it truncates or rounds the last
digit, then I can tell you what it should return for the examples.
BRIEF VERSION OF THE CHALLENGE: I have a decimal fraction. I want to
divide two integers ON AN HP48 or HP49 and get EXACTLY that decimal
fraction. What are the SMALLEST two integers that'll work?
> I would of course like to see your HP48 program...
> especially when it is fixed.
It is working perfectly... unlike any other "decimal to fraction"
program or function that I know of. It'll be posted here after this
challenge has had enough time to simmer.
HOMEWORK: What is the BEST fraction for sqrt(LOG(2)) on an HP48 or
HP49?
-Joe-
Virgil
joe...@holyjoe.net (Joseph K. Horn) wrote:
>
> HOMEWORK: What is the BEST fraction for sqrt(LOG(2)) on an HP48 or
> HP49?
>
> -Joe-
For LOG, base 10, I get 998995/1820784, which works but may not be
best.
For LN, base e, I get 462042/554969, which also works, but may not
be best.
James
> HOMEWORK: What is the BEST fraction for sqrt(LOG(2)) on an HP48 or
> HP49?
>
> -Joe-
On my 49g+, using XQ: 29422/53625
Tony Hutchins
Hi Virgil ! <vmh...@comcast.net>
05h14m ago, on Sat, 11.10.03 at 10:22 p.m. -0600, you wrote
in message ID <vmhjr2-016FFC....@news.newsguy.com> :
> For LOG, base 10, I get 998995/1820784, which works but
> may not be best.
FWIW that gets my vote as well. I haven't a clue how to get
there is less than 5 seconds though. The continued fraction
method leapfrogs way over that.
--
Tony Hutchins
New Zealand
#71 All computers wait at the same speed.
Eddy Hondo
joe...@holyjoe.net (Joseph K. Horn) wrote in message news:<87233f9e.03101...@posting.google.com>...
> James wrote:
>
> The "correct" answer depends on the machine, just like "the correct
> sine of 2" depends on the machine. If you tell me how many digits the
> FX-502P calculates with, and whether it truncates or rounds the last
> digit, then I can tell you what it should return for the examples.
How many digits does the 49g calculate with?? Does it truncate or
round the last digit??
Need that for a shot at the program. TY
Joseph K. Horn
>> HOMEWORK: What is the BEST fraction for
>> sqrt(LOG(2)) on an HP48 or HP49?
>
> For LOG, base 10, I get 998995/1820784,
> which works but may not be best.
>
> For LN, base e, I get 462042/554969,
> which also works, but may not be best.
Congrats! You got the best fractions for both cases!
PUH-LEEEEEZE tell me how you did it! Did you use an HP48 or 49? How
long did it take? I'm dying to know!!!
For what it's worth, my program converts sqrt(LOG(2)) into
998995/1820784 in 0.745 seconds on an hp49g+, and runs unmodified on
an HP48GX.
-Joe-
FX-502P Geek
OK Joseph I see what mean, although I'm not sure that it is useful. I
suspect that the System RPL of '->Q' calculates errors on the at a higher
accuracy than is returned to the user and so the denominator is closer to
the true answer of the input number. Your answer is also valid for your
interpretation on the FX-502P which also returns exactly 50000.00005 for
999050001/19981.
>
> > Incidentally the program that I wrote for the FX-502P
> > is available at www.fx502p.pwp.blueyonder.co.uk and
> > gives the correct answer of 50000+1/20000.
>
> The "correct" answer depends on the machine, just like "the correct
> sine of 2" depends on the machine. If you tell me how many digits the
> FX-502P calculates with, and whether it truncates or rounds the last
> digit, then I can tell you what it should return for the examples.
This is not straight forward. The FX-502P works to 12 places (but truncates
to 10 if you put a value in a memory!). It truncates to 12 places... unless
the last 2 digits are greater than 49 and the 8th, 9th & 10th digits are 999
in which case it rounds up the 9s to carry forward into the 7th digit. This
rather dubious logic prevents decimal errors leading to .99999... when the
next integer is the correct answer.
>
> BRIEF VERSION OF THE CHALLENGE: I have a decimal fraction. I want to
> divide two integers ON AN HP48 or HP49 and get EXACTLY that decimal
> fraction. What are the SMALLEST two integers that'll work?
>
> > I would of course like to see your HP48 program...
> > especially when it is fixed.
>
> It is working perfectly... unlike any other "decimal to fraction"
> program or function that I know of. It'll be posted here after this
> challenge has had enough time to simmer.
I look forward to it.
Joseph K. Horn
>> What is the BEST fraction for sqrt(LOG(2))
> On my 49g+, using XQ: 29422/53625
Close, but no cigar. 29422/53625 is not EQUAL to sqrt(LOG(2)) on the
HP48/49. The challenge, as clearly stated, requires that the fraction
evaluate to EXACTLY the same number as the input. Using your
solution, the input and output are, respectively:
.548662004939
.548662004662 (compare last three digits)
Even being off by 1 digit in the last place invalidates a solution;
this answer is off by 277 digits in the last place. Not acceptable.
XQ answers a diferent question than my challenge. XQ is meant for the
removal of "noise" in numbers that are known to be clean ratios. For
example, suppose you see a reference in a book to "0.8947368421" and
the context makes it clear that they're talking about some kind of
ratio of reasonably small integers. They obviously rounded it off to
10 digits, but HP uses 12 digits, so applying ->Q to that decimal
yields a ratio of two huge integers... clearly not what we're looking
for. Applying XQ yields the nice answer of 17/19.
WHICH ANSWER IS "RIGHT" DEPENDS ON THE QUESTION.
If the question is instead, "What are the smallest two integers which
divide on the HP48 or HP49 to get a result of exactly 0.8947368421?"
then the best answer is neither what ->Q yields nor what XQ yields,
but is instead 8171112721/9132420100 (obtained on my hp49g+ in 1.12
seconds).
-Joe-
"It is difficult to answer when one does not understand the question."
-- Sarek (Star Trek IV)
BlackLotus
I searched for other solutions, and it puzzles me how you can do that in
less than 5 seconds. Here's why:
Bruteforcing my way down a big problem: it can take a lot of time if ->Q
gives a wrong result (one that once ->NUMed is different from the input).
sqrt(ln(2)) is a good example. Another problem is a stopping condition, if I
start from a valid ->Q result, I can tell it to keep lowering Q until it
reaches a Q that makes Q X * 0 RND Q / different from X
A less calculating method seems to be, starting from the "right" fraction
(not the best):
Let's call the fraction Q/D (not best) and Qb/Db (best).
If you divide Q by d from D going down, the stop point is the result which
has the lowest number of digits, while still having RND(d*X)/d == X As you
can guess, the problem is to have Q and D to start with, since ->Q can gives
wrong results, and so does a continuous fraction method (tried it, it fails
because of the approximation, I might retry it with different parameters
than the one I found (it was a VB algorithm) to see if I can adjust it to
the HP's precision).
In both cases, what really kills me is how to find a starting point that's
close enough and/or a method to converge faster to a solution.
"BlackLotus" <black...@nospam.me> a écrit dans le message de
news:bma85n$rgf$1...@news-reader2.wanadoo.fr...
Rafael Mill?n
> It is because of a tiny approximation error. I checked it on a Math soft on
> the PC, which performs "exact" calculation, and substracting the two
> fractions shows an error of 19/399620000. So, no, Joe's answer's
> mathmatically wrong, it's just a twist from the HP48/49 that makes it
> "right" (4E-8 is neglected in regard to the number of digits and their order
> in 50000.00005 as "too small to be significant")
The 49 CAS itself (in exact mode) will do:
'1000000001/20000 - 999050001/19981' EVAL
'-19/399620000'
But the challenge does not ask for a mathematically exact result, but
an exact result when the fraction is approximated to a 48 or 49 real.
RaM.
Virgil
James <nos...@ham.org> wrote:
On my 49, in approximate mode,
'SQRT(LOG(2,)' evaluates to .548662004939
but
'29422/53625' evaluates to .548662004662.
Virgil
joe...@holyjoe.net (Joseph K. Horn) wrote:
> Virgil wrote:
>
> >> HOMEWORK: What is the BEST fraction for
> >> sqrt(LOG(2)) on an HP48 or HP49?
> >
> > For LOG, base 10, I get 998995/1820784,
> > which works but may not be best.
> >
> > For LN, base e, I get 462042/554969,
> > which also works, but may not be best.
>
> Congrats! You got the best fractions for both cases!
>
> PUH-LEEEEEZE tell me how you did it! Did you use an HP48 or 49? How
> long did it take? I'm dying to know!!!
I used a very primitive and very slow method based on Farey
fractions. as outlined below.
Let x be the decimal number to be aproximated, and assume it is not
already an integer.
Find integers n1 < x and n2 > x, and start with fractions n1/d1 and
n2/d2 where d1 = d2 = 1.
Then decrease the interval by taking n/d in place of the appropriate
endpoint, where n = n1 + n2 and d = d1 + d2.
Eventually, but very slowly, n/d gives the correct decimal on the
49, at least for all the decimals I have tried. I also did several
tests using reals for which XQ gives the wrong values without
getting any wrong answers.
A variant on continued fractions should give the same results in
considerably less time, but I haven't got the bugs worked out of it
yet.
>
> For what it's worth, my program converts sqrt(LOG(2)) into
> 998995/1820784 in 0.745 seconds on an hp49g+, and runs unmodified on
> an HP48GX.
>
> -Joe-
My program on my 49 takes about 53 seconds for the same result.
BlackLotus
"Virgil" <vmh...@comcast.net> a écrit dans le message de
news:vmhjr2-C09590....@news.newsguy.com...
> I used a very primitive and very slow method based on Farey
> fractions. as outlined below.
>
> Let x be the decimal number to be aproximated, and assume it is not
> already an integer.
> Find integers n1 < x and n2 > x, and start with fractions n1/d1 and
> n2/d2 where d1 = d2 = 1.
>
> Then decrease the interval by taking n/d in place of the appropriate
> endpoint, where n = n1 + n2 and d = d1 + d2.
>
> Eventually, but very slowly, n/d gives the correct decimal on the
> 49, at least for all the decimals I have tried. I also did several
> tests using reals for which XQ gives the wrong values without
> getting any wrong answers.
>
> A variant on continued fractions should give the same results in
> considerably less time, but I haven't got the bugs worked out of it
> yet.
>
Thank you so much. I finally managed to get something that works.
I thought about that problem all day long, and I never came to think about
that. Maybe I should have taken maths and not economics. :p
> My program on my 49 takes about 53 seconds for the same result.
It takes less than 3 seconds on my HP48 using the method you described; my
program takes 357 bytes.
Joseph K. Horn
> How many digits does the 49g calculate with??
> Does it truncate or round the last digit??
User RPL calculates internally with a truncated 15-digit mantissa,
then rounds that down to 12 digits and returns that to the user as the
final answer.
-Joe-
Sam Hughes
Try it with one of his examples...
0.69999937 XQ produces '111106/158723'. Then pressing ->NUM produces
the decimal .699999369972, which means that XQ's result is not the
best possible fraction. The best fraction, of course, is
'777777/1111111'.
reth
Veli-Pekka Nousiainen
news:3f8a9c77$1...@dnews.tpgi.com.au...
> what about 3.14159265359?
pi is a special case.
It is taken into consideration as a 15-nibble version internally
and the 12 number version si actually 13 digits precise:
3.14159265358979
15 => round to 13
3.141592653590
BUT what is the answer, Joe?
VPN
BlackLotus
Tony Hutchins
Hi reth ! <re...@nospan.com.at>
03h13m ago, on Mon, 13.10.03 at 10:37 p.m. +1000, you wrote
in message ID <3f8a9c77$1...@dnews.tpgi.com.au> :
> what about 3.14159265359?
I the ->Q answer for that (1146408/364913) satisfies the
challenge.
Like others, I lost a few hours yesterday, but the challenge
did come with a warning<G>. My result works well on most
examples, but fails *miserably* on .50000000001 and also the
HP38G number. I start with a continued fraction method and so
lose the last "1" as .500000000001 1/X 1/X = .5 only :(
Maybe I need a MANT somewhere. But, where?<G>. I look forward
to seeing Joe's code :)
--
Tony Hutchins
New Zealand
#2 Save whales - collect the whole set.
BlackLotus
part of the problem.
I tried Virgil's method, and it seemed to work fine until I tried a number
like 50000.00005. There it kind of killed my hopes, because it takes
horribly long to get the answer. I looked a little, and it seems it chockes
on the limit cases such as 0.001 (it requires 999 iterations), especially
when it uses all the digits and a lot of them are similar. I think my
iterations can go up to (1/X)-1 (or 1/(1-X)-1), which means close 1E11
iterations to reach a result for 0.999999999975. Way too long.
I spent a good while trying to figure how to get the solution with as little
iterations as possible, and while some CF algorithms looked nice, they often
return the exact result, not the best result. I tried changing the way the
interval is shrunk, but it doesn't cut much time.
I looked deeper and tried to find what is special with the result we want,
and I found a few things.
First, if X is my number, and Q/D the fraction I'm looking for, then it
appears that D*X <> Q, but Q/X = D, due to rounding. We know X, we'll search
Q, then. Also, the Q/X we're looking for has to particualrity of being an
integer, which means that CEIL(Q/X)=FLOOR(Q/X)=IP(Q/X)=Q/X. How does that
help? Well, the Q that we're looking for is the smallest Q so that
IP(Q)/CEIL(Q/X)=X (or FLOOR(Q/X)/CEIL(Q/X)=1, with Q an integer. The problem
is, the 48's built-in solver stops on extremums or sign changes, or leaps
over the solution.
So, you guessed it, my problem is to find the solution to that fairly simple
equation, given that the built-in solver chokes (not to say "is useless"),
and that derivation won't work. So, the possibilities I see are these:
Either I can find a solving algorithm and chances are it will go faster than
the convergence I used so far, or I can find a derivable mathematical
formula that says exactly the same as my IP(Q)/CEIL(Q/X)-X=0, except without
IP and CEIL, since the HP won't derivate them.
I don't know if they're the way Joe chose (I doubt it, I usually find the
most complicated solution to a problem), but I'd like to explore it, and I
don't have the knowledge to solve that equation. I figured it could be
interesting, since ideally, finding a solution would take pretty much the
same time regardless of the starting value.
Joseph K. Horn
> what about 3.14159265359?
The best fraction for pi on an HP 12-digit calculator is what ->Q
returns for it in STD mode: 1146408/364913.
Same is true for pi^e, which for all I know is a rational number,
though I doubt it.
-Joe-
Tony Hutchins
Hi BlackLotus ! <black...@nospam.me>
35m ago, on Mon, 13.10.03 at 7:37 p.m. +0200, you wrote
in message ID <bmentb$84o$1...@news-reader5.wanadoo.fr> :
> I have a question here, because I lack the mathematical
> knowledge to solve a part of the problem.
I lack the knowledge of how to differentiate those "IP"
functions as well - you have step functions there. Possibly
there is no classical way to handle them.
> I tried Virgil's method, and it seemed to work fine until
> I tried a number like 50000.00005. There it kind of killed
> my hopes, because it takes horribly long to get the answer.
The CF algorithm really helps on that one - it gives the 20000
denominator instantly. AFAIK the only way to test if smaller
denominators work is to decrement the 20000 by 1 until it
doesn't work. It doesn't take long. I suspect some process
like this is what makes Joe's algorithm sometimes take 5
seconds, instead of under 1 second.
> I looked a little, and it seems it chockes on the limit cases
> such as 0.001 (it requires 999 iterations), especially when
> it uses all the digits and a lot of them are similar. I
> think my iterations can go up to (1/X)-1 (or 1/(1-X)-1),
> which means close 1E11 iterations to reach a result for
> 0.999999999975. Way too long. I spent a good while trying
> to figure how to get the solution with as little iterations
> as possible, and while some CF algorithms looked nice, they
> often return the exact result, not the best result.
But CF is great to get near the best answer - I think
successive results always bound the answer. The trick is not to
let the CF go too far. A good place to stop is when
Denominator*X=Numerator, which often happens before the exact
result.
[...]
BlackLotus
"Tony Hutchins" <jus...@csi.com> a écrit dans le message de
news:1049541...@news.cis.dfn.de...
> -=[ Tue, 14.10.03 07:13 a.m. +1300 (NZDT) ]=-
>
> denominator instantly. AFAIK the only way to test if smaller
> denominators work is to decrement the 20000 by 1 until it
> doesn't work. It doesn't take long. I suspect some process
> like this is what makes Joe's algorithm sometimes take 5
> seconds, instead of under 1 second.
>
418507862138/762779012161, which is far from the best solution (and lower
than the exact solution), and going your way down from here is going to take
long. On the same idea, I realized that sometimes the exact fraction and the
best fraction are such that there's a set of values where the fraction is
true, then a set of value where it's false (I used that as a stopping
condition in the first try I made) and then back to a set where it's true
again like with sqrt(log(2)).
Virgil
"reth" <re...@nospan.com.at> wrote:
> what about 3.14159265359?
>
>
1146408/364913?
Tony Hutchins
Hi BlackLotus ! <black...@nospam.me>
29m ago, on Mon, 13.10.03 at 9:01 p.m. +0200, you wrote
in message ID <bmesqe$846$1...@news-reader1.wanadoo.fr> :
> But for example, on sqrt(log(2)), the CF algorithm I still have returns
> 418507862138/762779012161,
Calling that N/D, and the target X, I stop CF when I have D*X=N
This happens when D=470069 I think. This really helps in most
of the examples, but not all ... so it is not the "key" by any
means<G>
FX-502P Geek
The FX-502P actually has it built in pi wrong at 3.14159265358 i.e.
truncated to 12 places rather than rounded to 3.14159265359. Casio even got
it successor the FX-602P wrong with 3.14159265360!
Peter.
BlackLotus
> -=[ Tue, 14.10.03 08:31 a.m. +1300 (NZDT) ]=-
>
> This happens when D=470069 I think. This really helps in most
> of the examples, but not all ... so it is not the "key" by any
> means<G>
257909/470069 is different (by 1E-12, the last digit) from sqrt(log(2)).
And I know it's not the perfect solution, which puzzles me as how Joe did
it.
Tony Hutchins
Hi BlackLotus ! <black...@nospam.me>
21m ago, on Mon, 13.10.03 at 10:22 p.m. +0200, you wrote
in message ID <bmf1ie$ovm$1...@news-reader4.wanadoo.fr> :
> Just for the record, 470069*sqrt(log(2)) is an integer (257909), but
> 257909/470069 is different (by 1E-12, the last digit) from sqrt(log(2)).
Yup the CF method is pretty good at getting close, quickly.
What eludes me is a general rule :(
I tried lots of variations hoping that one would suddenly work
for all examples, and then I could see if I could understand
it by general reasoning, after playing about. I learned to
program quickly, and it was fun - I like working with the
49g+. It seems I was lucky there as my keyboard is fine. No
keys stand out as non-functional, and I don't use a hammer<G>.
> And I know it's not the perfect solution, which puzzles me
> as how Joe did it.
It will be very interesting to see the answer-the logic, and
the implementation.
I can't even do that example of .5 + E-12 - I keep losing the
E-12 <G>... and CF jumps to D=250,000,000,001 on the second
term!! If I divide that by 15 I get close. I almost have a
different program for each example<G>. I reckon if I could see
how to do that one ... I might see how to find 2 denominators
D1,D2 such that I know the answer is D1,D2 or something bigger
like D1+D2. Then if an answer is found it always needs
checking in case there is a lower equivalent, separated by
1,2,3 etc. Also I may be totally wrong about everything <G>.
What we need is probably a number theory expert. But it sounds
like Joe himself has broken new ground here!
Joseph K. Horn
> ... AFAIK the only way to test if smaller
> denominators work is to decrement the 20000
> by 1 until it doesn't work. It doesn't take
> long. I suspect some process like this is
> what makes Joe's algorithm sometimes take 5
> seconds, instead of under 1 second.
You're hot on the trail, but FYI I'm pretty sure that the slowest
input for my program is the "worst case" number 0.999999999975, which
takes 1.25 seconds (on a 49g+).
> But CF is great to get near the best answer -
> I think successive results always bound the
> answer. The trick is not to let the CF go too
> far.
That's one of the charms of my program: it never introduces ANY
roundoff error at any point, not even with 1/x. It is partly based on
continued fractions, but it does NOT use FP(1/x) [like everybody else
does] because that loses information at every iteration. My method
loses NO information, and in fact marches inexorably towards the best
fraction every time, with at most a handful of iterations, even with
inputs that are pure evil, such as the ones with large partial
quotients such as .500000000001, whose continued fraction expansion is
exactly [0; 1, 1, 249999999999, 2].
MASSIVE HINT: Write a program that converts any input number to its
*exact* continued fraction expansion, and vice versa. From there to
the solution to this Best Fraction Challenge is just a hop, skip and a
jump. Literally.
In case it helps, here are the *EXACT* continued fraction expansions
of the examples that have appeared in this thread so far. Instead of
the usual notation [integer part; partial quotient, pq, pq, ...], I'll
just show 'em in a list, HP style. Remember, the operative word here
is EXACT! For example, pi is not exactly equal to
314159265359/100000000000, but HP's version of pi *is* exactly equal
to that, mathematically. Let's call it "HP pi".
HP 2 SIN (in RAD mode) =
{ 0 1 10 39 1 12 1 2 1 45 1 25 1 22 14 1 1 24 }
6.0030000004 =
{ 6 333 3 2499 1 2 333 }
0.999999999975 =
( 0 1 3999999999 }
0.9498580222 =
{ 0 1 18 1 16 1 1 1 12 1 3 24 3 6 5 1 14 2 }
50000.00005 =
{ 50000 20000 }
0.500000000001 =
{ 0 1 1 249999999999 2 }
0.69999937 =
{ 0 1 2 3 15872 1 2 1 1 12 3 2 }
0.903493239856 =
{ 0 1 9 2 1 3 4 1 2 6 1 5 1 6 7 1 11 15 18 4 }
HP sqrt(LOG(2)) = 0.548662004939 =
{ 0 1 1 4 1 1 1 3 7 7 9 1 6 1 5 1 1 2 3 10 2 20 1 7 2 2 }
HP sqrt(LN(2)) = 0.832554611158 =
{ 0 1 4 1 34 1 5 6 3 1 1 9 6 1 14 1 6 1 1 178 1 2 }
0.8947368421 (almost 17/19) =
{ 0 1 8 2 526315789 }
17/19 (*not* HP 17/19, but the real thing) =
{0 1 8 2 }
HP pi =
{ 3 7 15 1 292 1 1 1 2 1 4 1 1 1 1 1 3 1 68 2 4 2 }
HP pi^e = 22.4591577184 =
{ 22 2 5 1 1 1 1 1 3 2 1 1 3 8 1 8 1 11 1 3 4 2 }
-Joe-
Tony Hutchins
Hi Joseph K. Horn ! <joe...@holyjoe.net>
55m ago, on Mon, 13.10.03 at 4:16 p.m. -0700, you wrote
in message ID <87233f9e.03101...@posting.google.com> :
[...]
> even with inputs that are pure evil, such as the ones
> with large partial quotients such as .500000000001,
> whose continued fraction expansion is exactly [0; 1, 1,
> 249999999999, 2].
I get [0; 1,2], which converts back to .666', so indeed it is
pure evil, evil extended <G>.[1]
Many many thanks for your tips Joe! I see the problem I now
have to deal with. I seem to have to do reciprocals - I can
see how to recover lost digits (DUP 1/x 1/x - for example),
but I can't see what to do with them yet - especially with the
above example I want a reciprocal of the .500000000001 that is
at least under 2. A nice puzzle - how to adjust a reciprocal.
> MASSIVE HINT: Write a program that converts any input number
> to its *exact* continued fraction expansion, and vice versa.
> From there to the solution to this Best Fraction Challenge
> is just a hop, skip and a jump. Literally.
MANY THANKS!!!!
Literally? I can't find these hop,skip & jump functions in user
RPL <G>
That should keep me on the straight and narrow over the next
month ;-)
Today, for the first time, I really followed the continued
fraction notation - in the past all those divisions made me
dizzy <G>
I already tried to reverse the process but failed there too -
with the evil one above - going backwards in RPN style from
the end of the fraction, I lost significance at the first
+ below:
2 1/x 249999999999 + 1/x 1 + 1/x 1 + 1/x 0 +
... and I just end up with .5 again. Still that's better than
.666<G>[1]. Once I crack this I'll be on the way, I hope.
So, I have to find a way to avoid reciprocals - but this looks
rather tricky, given the usual simple continued fraction
notation :)
No, no, I just need accurate reciprocals.
Or somehow turn everything on its head.
> In case it helps, here are the *EXACT* continued fraction expansions
> of the examples that have appeared in this thread so far.
[...]
Thanks for the detailed notation and description - I can now
follow what you mean!
I'll be working on this on the train and hope to soon see the
49g+ march ineluctably to the correct answers - after I figure
out the hop and the skip and the jump :) I can almost imagine
I will have to hop to the right denominator interval, and
maybe skip to one end and jump about a bit in the interval -
*that* part I have been doing already<G>.
This is a terrific challenge - although I consider it
recreational I feel it also might have uses in computer
science - hope it does!
[1] Oops I see I made a mistake above - for the .500000000001 I
get [0;2], not [0;1,2] which converts back to .5, not .666'. Unintentional.
I was learning the notation as I wrote.
--
Tony Hutchins
Wellington
New Zealand
#248 I praise loudly; I blame softly. Queen
Catherine II
Tony Hutchins
Hi Joseph K. Horn ! <joe...@holyjoe.net>
09h16m ago, on Mon, 13.10.03 at 4:16 p.m. -0700, you wrote
in message ID <87233f9e.03101...@posting.google.com> :
[...]
Joe, me again. I made progress :)
> inputs that are pure evil, such as the ones with large partial
> quotients such as .500000000001, whose continued fraction expansion is
> exactly [0; 1, 1, 249999999999, 2].
I can do this now. Took me a while to un-confuse myself<G>
> MASSIVE HINT: Write a program that converts any input number to its
> *exact* continued fraction expansion, and vice versa. From there to
> the solution to this Best Fraction Challenge is just a hop, skip and a
> jump. Literally.
I seem to have a problem doing the reverse - getting back to
the exact number..see below
[...]
> For example, pi is not exactly equal to
> 314159265359/100000000000,
*That* is what gave me the clue to doing the expansion - I
don't start with .500000000001 but that*E12 /E12 - I start
with an integral numerator and a denominator. Thanks for that
hidden tip!!
Than I saw a worked example of old Euclid's algorithm for GCD
and that gave me another hint. This is stuff I never studied
at school.
[...]
The problem:
Input: 0.500000000001
Output: 166666666668 / 333333333335
The expansion in partial quotients:
[0; 1, 1, 249999999999, 2].
a0 a1,a2, a3 , a4 (notation)
To reverse the process I assume I calculate successive
"convergents":p1/q1,p2/q2,p3/q3 and p4/q4.
I thought p4/q4 would give me the original .500000000001
p0=a0=0 p(-1)=1 q0=1 q(-1)=0
p1=a1*p0+p(-1)=1*0+1 =1
q1=a1*q0+q(-1)=1*1+0 =1 p1/q1=1/1=1
p2=a2*p1+p0=1*1+0 =1
q2=a2*q1+q0=1*1+1 =2 p2/q2=1/2=.5
p3=a3*p2+p1=249999999999*1+1 =250000000000
q3=a3*q2+q1=249999999999*2+1 =499999999998
Oh!!! I see the error!! q3 sb=499999999999
p3/q3=.500000000001
Yeah!
p4=a4*p3+p2=2*250000000000+1 =500000000001
q4=a4*q3+q2=2*499999999999+2 =E12
p4/q4=.500000000001 - exactly.
OK, I think I can now do the reversal.
I seem to be stumped as to how to "hop, skip and jump" to
near the 33333333335 denominator.
Is there a different way to do the convergents?
Maybe backwards?
You can tell I haven't the foggiest idea how to be sure I am
homing in on the Best Fraction for the HP48/49. For all I know
it could be lurking almost anyware<G>.
However I do have p3/q3. So I must check backwards with
smaller denominators until I find the earliest example where
the "p/q" equals our value, on the HP48/49.
Ahha .. maybe it is as I once suspected that there is a whole
series of *consecutive* results that give the answer. So,
maybe I:
1. hop to the first endpoint that gives the answer.
2. That means the previous convergent doesn't give the
answer<G>. So I just skip to the denominator in the middle,
and check that out.
3. This way I establish a new interval to look in, and I jump
to that one :)
Yup, successive bisection - might be a general method.
Joe, in writing this I think I found my way! I started to
write because I was stumped.
I'm still a bit worried that somehow an earlier value may
exist - but although the results go up and down in a sawtooth
like way, overall they are probably well-behaved.
Thanks so much for the tips!!!!!
--
Tony Hutchins
Wellington
New Zealand
#77 Allen's Axiom: When all else fails, follow instructions.
Sam Hughes
Horn says...
> For what it's worth, my program converts sqrt(LOG(2)) into
> 998995/1820784 in 0.745 seconds on an hp49g+, and runs unmodified on
> an HP48GX.
How long does it take to run on a 48GX, or on a 49G?
Joseph K. Horn
> > For example, pi is not exactly equal to
> > 314159265359/100000000000,
>
> *That* is what gave me the clue to doing
> the expansion - I don't start with .500000000001
> but that*E12 /E12 - I start with an integral
> numerator and a denominator. Thanks for that
> hidden tip!!
That's the trick for avoiding the 1/FP(x) roundoff errors that have
plagued decimal-to-fraction algorithms forever.
"You have taken your first steps into a larger world."
-- Obi-Wan Kenobi
> I seem to be stumped as to how to "hop, skip
> and jump" to near the 33333333335 denominator.
MASSIVE HINT (spoiler, actually) #2: With the list of partial
quotients handy, examine the *differences* between successive
intermediate convergents' numerators. Now do the same for their
denominators. When the pattern hits you, it'll hit you hard, so be
sure to be sitting down!
> Yup, successive bisection - might be a general method.
Put that together with the previous paragraph, and you have my
program, sans bells & whistles.
> although the results go up and down in a sawtooth
> like way, overall they are probably well-behaved.
Indeed they are. Have fun!
-Joe-
Tony Hutchins
Hi Tony Hutchins ! <jus...@csi.com>
17h17m ago, on Tue, 14.10.03 at 10:30 p.m. +1300 (NZDT), you wrote
in message ID <1147579...@news.cis.dfn.de> :
> Joe, me again. I made progress :)
Joe, I can now make all the examples work. The final hopping
skipping and jumping *did* stump me - took about 5 attempts.
I have it working, based on sort of heuristic theory ;-)
eg .500000000001 BF leaves this in the stack:
166666666668
333333333335
Here is BF, with some comments where the line begins with '
<<
' don't allow negative input
ABS 'Z' STO
' create partial quotients, and the convergents
Z MAN MPQN
'could use PN QN / and a counter
'now, HOP to the first convergent which produces
'an "exact-HP48/49 result"
' This is called N2/D2 & the previous one is NJ/DJ
' DJ is actually the JUMPING unit
1 'D2' STO
PN 1 GET 'N2' STO
1 'K' STO
WHILE N2 D2 / Z NOTEQUAL (i.e. "=/=" )
REPEAT
K 1 + 'K' STO
N2 'NJ' STO
PN K GET 'N2' STO
D2 'DJ' STO
QN K GET 'D2' STO
END
' now we do the jumping backwards from N2/D2
' we start about the middle of the whole interval 0->D2
' There, we find D3 and N3
' The trick is that we hop back from D2
' in steps of DJ
' If N3/D3 is indeed a BETTER FRACTION candidate we use D3
' as D2 - otherwise we shift D1 up to D3, and work
' from the middle of the new D1->D2
' J is the distance from 'D3 to D2' measured in steps of
' DJ
0 'D1' STO
D2 DJ / 2 / IP 'J' STO
WHILE J 0 >
REPEAT
D2 DJ J * - 'D3' STO
N2 NJ J * - 'N3' STO
N3 D3 / Z ==
<<N3 'N2' STO D3 'D2' STO>>
<<D3 'D1' STO>>
IFTE
D2 D1 - DJ / 2 / IP 'J' STO
END
N2 D2
>>
MAN - this makes the AN list
of partial quotients
for a positive fraction Z
I make the list end with a 0.
<< ->Z
<< {} Z IP + 'AN' STO
E12 DUP Z FP *
WHILE DUP 0 >
REPEAT
DUP 'D1' STO NDQR
SWAP AN SWAP + 'AN' STO
D1 SWAP
END
AN 0 + 'AN' STO
>>
>>
NDQR (used by MAN)
takes positive integers N D in stack and outputs integers Q R
where N=Q*D + R and R is less than D, and not negative.
The bit before the WHILE is to quickly get a good first guess
for Q.
This puts Q and R on the stack as well as making global
variables for them - I'm in a hurry<G>
<< -> N D
<< N D / IP DUP 0 > <<1 ->>IFT 'Q' STO
N D Q * - 'R' STO
WHILE R D >=
REPEAT
R D - 'R' STO
Q 1 + 'Q' STO
END
Q R
>>
>>
Here is the last piece -
MPQN: uses AN (list of partial quotients)
and produces PN and QN - lists of pk and qk
such that pk/qk is kth convergent.
<<
1 'P0' STO 0 'Q0' STO
AN 1 GET 'P1' STO 1 'Q1' STO
{} P1 + 'PN' STO
{} Q1 + 'QN' STO
AN 2 GET 'AK' STO 2 'K' STO
WHILE AK 0 >
REPEAT
AK P1 * P0 + 'P2' STO
AK Q1 * Q0 + 'Q2' STO
P1 'P0' STO P2 'P1' STO
Q1 'Q0' STO Q2 'Q1' STO
PN P2 + 'PN' STO
QN Q2 + 'QN' STO
K 1 + 'K' STO
AN K GET 'AK' STO
END
>>
None of the examples take more than 3-4 seconds, on a 49g+,so it
is 3-4 times slower than your program Joe.
I only just got this functional, with what I am sure is a very
limited subset of RPL<G>
I see now I could calculate pk/qk as I go, calculating the ak
and can stop as soon as I have the first "exact HP48/49"
result. Also, I make far too many global variables. But, I'm
just so pleased to see it working :)
Really great challenge - thanks!
I'm not really sure about my jumping unit, but, as it does
'hit' on the likely candidates, there must be something to
it<G>. I soon realised I couldn't just look at all
possibilities, unless of course the prior convergent has
denominator 1. The likely candidates are more or less
"quantized" by the prior convergent.
- Tony
Tony Hutchins
Hi Joseph K. Horn ! <joe...@holyjoe.net>
1 day 06h24m ago, on Mon, 13.10.03 at 4:16 p.m. -0700, you wrote
in message ID <87233f9e.03101...@posting.google.com> :
> In case it helps, here are the *EXACT* continued fraction
> expansions of the examples that have appeared in this thread
> so far.
It did help :)
Now I have it all at least functional (in case anyone keys in
the programs from my last post, the MAN can do with a DROP
DROP at the end), I've been experimenting.
The square root of 3 is fascinating - on the HP48/49 it
produces nearly 30 partial quotients!!!! Starting with 9 pairs
of 1 2. I get 716035/413403 for the best fraction.
Joe, you must tell us the story as to how you came to
investigate this - .. if there is a story :-)
Is it part of a book/paper you intend to publish?
--
Tony Hutchins
Wellington
New Zealand
#252 What you will do matters. All you need is to
do it. Judy Grahn
Joseph K. Horn
>> For what it's worth, my program converts
>> sqrt(LOG(2)) into 998995/1820784 in 0.745
>> seconds on an hp49g+, and runs unmodified on
>> an HP48GX.
>
> How long does it take to run on a 48GX,
> or on a 49G?
It takes approximately 3 times longer, for almost any input.
-Joe-
Joseph K. Horn
Challenge! In fact, he's the first person to post a program that
finds the two smallest integers which divide to any given decimal in
less than 5 seconds! Well done!
> Joe, you must tell us the story as to
> how you came to investigate this -
> .. if there is a story :-) Is it part
> of a book/paper you intend to publish?
Sorry, no story. Just one of the things I've stumbled into while
playing around with numbers and HP calculators for the past 27 years.
However, I'm currently working with HP to see that this algorithm,
which I call the PDQ algorithm ("P Divided by Q"), will be used in
place of the current ->Q algorithm in future models. However, my
proposal goes beyond that, allowing the accuracy to be specified
several ways: by a maximum tolerable error, or maximum tolerable
denominator, or minimum exact-match decimal places, or real-world
tolerance by a user-specified scalar. Finally, some sort of Fraction
Hunter application (Fraction Ranger? Fraction Scout? Best Fraction
Finder? PDQ?) would allow the user to see simultaneously the input,
the output, the previous and next fraction (if any), their respective
errors, and real-world scalar error, with options to single-step and
back-step through the continued-fraction convergents and/or
intermediate convergents. This would be helpful in many real-world
applications that rely on precise tolerances, and it would also be
ideal for teaching fractions, decimals, continued fractions, rationals
and irrationals, and many other concepts.
If HP doesn't bite the bait, I'll just "publish" it all on the web for
all to see and use as they wish. Meanwhile, have fun with it!
Legal Fine Print for TI, Casio, Sharp et al: The PDQ Algorithm is the
intellectual property of Joseph K. Horn. Ownership was demonstrated
at HHC 2003 on September 21, 2003, in front of dozens of witnesses,
and the entire proceedings are available on videotape. Any use of
this algorithm for commercial purposes without proof of prior
ownership or written permission of Joseph K. Horn would be a Bad
Thing. You Have Been Warned.
-Joe-
Joseph K. Horn
p...@holyjoe.us
Gerald Squelart
[...]
> The PDQ Algorithm is the
> intellectual property of Joseph K. Horn. Ownership was demonstrated
> at HHC 2003 on September 21, 2003, in front of dozens of witnesses,
> and the entire proceedings are available on videotape. Any use of
> this algorithm for commercial purposes without proof of prior
> ownership or written permission of Joseph K. Horn would be a Bad
> Thing. You Have Been Warned.
> -Joe-
> Joseph K. Horn
> p...@holyjoe.us
Hi Joe!
AFAIK, unless you patent your algorithm for every possible way it may be
used, anybody may use it. All you own at the moment is the program you wrote
(that uses this algorithm) and I guess it could be named 'Horn Algorithm' if
it's truly new. So if you don't want others to freely use your algorithm,
patent it or keep it for yourself; Both being Bad Things anyway in my
opinion.
IANAL, of course, so prove me wrong if I am ;-)
HAGD,
G.
Tony Hutchins
Hi Joseph K. Horn ! <joe...@holyjoe.net>
1 day 30m ago, on Tue, 14.10.03 at 4:08 p.m. -0700, you wrote
in message ID <87233f9e.03101...@posting.google.com> :
> > *That* is what gave me the clue to doing
> > the expansion - I don't start with .500000000001
> > but that*E12 /E12 - I start with an integral
> > numerator and a denominator. Thanks for that
> > hidden tip!!
>
> That's the trick for avoiding the 1/FP(x) roundoff errors that have
> plagued decimal-to-fraction algorithms forever.
They have indeed - I looked into decimal->fraction using
continued fractions earlier this year. And everything I saw
used the 1/x. I started with the PPC ROM version - always a
good place to start :)
Joe, I don't know how, but I missed this post yesterday -
despite eagerly looking out for another hint!!
> "You have taken your first steps into a larger world."
> -- Obi-Wan Kenobi
Indeed, the wholesome world of integers :)
> > I seem to be stumped as to how to "hop, skip
> > and jump" to near the 33333333335 denominator.
>
> MASSIVE HINT (spoiler, actually) #2: With the list of
> partial quotients handy, examine the *differences* between
> successive intermediate convergents' numerators. Now do
> the same for their denominators. When the pattern hits
> you, it'll hit you hard, so be sure to be sitting down!
Ahha! I missed this spoiler, and I only started looking at
this out of desperation<G>. I noticed one example by chance -
that 2 LOG SQRT case. I got there by trial and error with
little understanding. But, what a thrill to see it work for
all cases! Actually successive convergents are very closely
related - multiplying numerators by the others denominators
and taking the difference one seems to get 1!!!
[...]
> > although the results go up and down in a sawtooth
> > like way, overall they are probably well-behaved.
>
> Indeed they are. Have fun!
I did - and now I might try and improve the RPL by using INCR
and STO+ etc and see how many globals I can just keep on the
stack instead. And I should prove to myself that there are no
other better candidates lurking.
Thanks for your more recent post too - telling about "PDQ" :)
Fascinating. It has a certain beauty and simplicity, once one
sees it :) Thanks for sharing it with us!
--
Tony Hutchins
Wellington
New Zealand
#229 It's not over till it's over. Yogi Berra
Rodger Rosenbaum
"the smallest integers which, when
divided on the calculator, yield exactly the same number as the input."
Let us call this a Best(JH) fraction. There might be other "best" fractions which are not
Best(JH). To use one of Joe's examples:
"Input: 6.0030000004 (the strange number on the HP38G case)"
"Output: 14823406 / 2469333"
Two other fractions which both evaluate to 6.0030000004 on the HP48 are:
15009499/2500333 and 15195592/2531333
They have special properties; who can characterize them?
I think one might say that "best" is in the eye of the user.
For example, suppose one were designing a synthesiser to generate a 6.0030000004 Hertz
frequency. One way to do it would be to buy a crystal oscillator running at a frequency
of 14823406 Hertz and divide that output by 2469333; this would be the Best(JH) way to do
it, but would it be the "best" way? What does "best" mean for you?
Joe's challenge is fascinating and difficult, but I'm not sure that the Best(JH) fraction
is the "best" except in curiosity and educational value.
Can anybody think of a situation where 14823406/2469333 would be preferred to
15009499/2500333 as an approximation to 6.0030000004? When would the slightly smaller
numerator and denominator of 14823406/2469333 compared to 15009499/2500333 be of value,
given the special property of 15009499/2500333 (said property to be determined by the
audience). The ->Q function on the HP48 gives 15003496/2499333 for 6.00300000004 (in STD
mode), which lacks the special property that 15009499/2500333 has.
This posting is not meant to belittle or minimize this latest excellent challenge of Joe
Horn (who is a friend of mine). :-)
Tony Hutchins
Hi Joseph K. Horn ! <joe...@holyjoe.net>
09h42m ago, on Wed, 15.10.03 at 12:25 p.m. -0700, you wrote
in message ID <87233f9e.03101...@posting.google.com> :
> Congratulations to Tony Hutchins! He wins the Best Fraction
> Challenge! In fact, he's the first person to post a program that
> finds the two smallest integers which divide to any given decimal in
> less than 5 seconds! Well done!
Thanks Joe!
I'll still be looking at this till Christmas I think.
My old paper HP49G Advanced User's Guide (the one with the
black cover - must be rare<G>), just told me that IDIV2 does
exactly what my 'NDQR' code does. And it works great too - eg
E12 500000000001 IDIV2 gives 1 and 499999999999 - even in
"approx" mode (contrary to the manual).
> > Joe, you must tell us the story as to
> > how you came to investigate this -
> > .. if there is a story :-) Is it part
> > of a book/paper you intend to publish?
>
> Sorry, no story. Just one of the things I've stumbled into
> while playing around with numbers and HP calculators for
> the past 27 years. However, I'm currently working with HP
> to see that this algorithm, which I call the PDQ algorithm
> ("P Divided by Q"), will be used in place of the current
> ->Q algorithm in future models.
Great idea!!! I look forward to that. It might even be in a
future 49g+ ROM :)
> However, my proposal goes beyond that, allowing the
> accuracy to be specified several ways: by a maximum
> tolerable error, or maximum tolerable denominator,
> or minimum exact-match decimal places, or real-world
> tolerance by a user-specified scalar.
Terrific. I see the old ->Q in there somewhere ;-)
> Finally, some sort of Fraction Hunter application (Fraction
> Ranger? Fraction Scout? Best Fraction Finder? PDQ?)
HP can invent a name :)
Harvester ?
Picker?
Your PDQ is so thorough in finding excellent and useful
fractions. And surprising fractions! It is almost forensic.
It finds new fractions. And it works right up to the limits of
the HP48/49 precision, because of that 15 digit thing you
mentioned, and the internal rounding.
> would allow the user to see simultaneously the input,
> the output, the previous and next fraction (if any), their
> respective errors, and real-world scalar error, with options
> to single-step and back-step through the continued-fraction
> convergents and/or intermediate convergents.
Will be very educational.
I hope there is also a built-in partial quotient generator.
They are so interesting. And the reverse. ->PQ and PQ-> :)
GCD must already do something like that internally.
> This would be helpful in many real-world applications that
> rely on precise tolerances, and it would also be ideal for
> teaching fractions, decimals, continued fractions, rationals
> and irrationals, and many other concepts.
Indeed. Now I know why your domain is holyjoe - you like whole
numbers :)
> If HP doesn't bite the bait, I'll just "publish" it all on
> the web for all to see and use as they wish. Meanwhile, have
> fun with it!
*I* am sure HP will use it :)
And soon, I hope.
I looked up my old 1965 Encyclopedia Britannica - it has a
great article on Partial Fractions. But everything I see even
there jumps to more complex features.
> Legal Fine Print
[...]
Hehe I have done you a favour then - anybody studying my code
would just get very confused - no chance of translating it :)
It's amazing the piracy that happens in business circles
though - recently I was pointed to an *non hp* financial
calculator which seems to be a 12C clone!!!!! (I ordered
one<G>).
Tony Hutchins
Hi Tony Hutchins ! <jus...@csi.com>
05h06m ago, on Thu, 16.10.03 at 1:45 p.m. +1300 (NZDT), you wrote
in message ID <1090147...@news.cis.dfn.de> :
> Joe, I don't know how, but I missed this post yesterday -
> despite eagerly looking out for another hint!!
I figured it out now - I missed 11 consecutive articles -
about the time I did a chkdsk /f and found cross-linked files.
My poor old DOS machine was the problem, not cis.dfn.de, nor
my home-grown newsreader. First time I have seen this -
everything looked like it was working fine, but ... the
articles came in and I'm not sure where they went<G>.
- Tony
Gjermund Skailand
"Joseph K. Horn" <joe...@holyjoe.net> wrote in message
news:87233f9e.03101...@posting.google.com...
> Congratulations to Tony Hutchins! He wins the Best Fraction
> Challenge! In fact, he's the first person to post a program that
> finds the two smallest integers which divide to any given decimal in
> less than 5 seconds! Well done!
I still havent got it right with the "Evil" ones,
But for what it is worth, in my longfloat (arbitrary precision) library
there are
3 simple commands
- for converting a longreal to a simple fraction
- for converting a simple fraction to a continued fraction
- for converting a continued fraction to a simple fraction
The continued fractions use lists. But no "minimum fractions":-(
Best regards
Gjermund Skailand
Werner Huysegoms
@ IDIV2 - Quotient and Remainder
@
@ bytes: 30.0
@ check: #C142h (48)
@ not needed on 49
@
@ In: 2: a
@ 1: b
@
@ Out: 2: q
@ 1: r
@
@ so that b*q + r = a
@
\<< DUP2 MOD ROT OVER - ROT / SWAP \>>
@ \->CF - Create Simple Continued Fraction
@
@ bytes: 98.0 (48)
@ 96.0 (49)
@ check: #E52Bh (48)
@ #5800h (49)
@
@ In: 1: x
@
@ Out: 1: { a0 a1 .. an }, ai integer (any i), and ai>0 for i>0
@
@ x = a0 + 1/(a1 + 1/(.. + 1/an))
@
\<<
DEPTH
\-> D
\<<
@ turn x into p/q
DUP MANT 1E11 * OVER / SWAP OVER *
@ stack holds 2:q 1:p, with ri = p/q
WHILE
@ calculate ai
OVER IDIV2 DUP
REPEAT @ stop when remainder is zero
@ calculate ri+1
ROT
END
ROT DROP2
D DEPTH SWAP - \->LIST
\>>
\>>
Cheers,
Werner
Tony Hutchins
Hi Werner Huysegoms ! <werner-h...@freegates.be>
02h16m ago, on Thu, 16.10.03 at 09:09 a.m. -0700, you wrote
in message ID <44ec85ff.03101...@posting.google.com> :
> Here's the first subprograms needed (working on the rest):
Thanks Werner - they are so ** elegant ** !
> @ IDIV2 - Quotient and Remainder
> @
> @ bytes: 30.0
> @ check: #C142h (48)
> @ not needed on 49
> @
> @ In: 2: a
> @ 1: b
> @
> @ Out: 2: q
> @ 1: r
> @
> @ so that b*q + r = a
> @
> \<< DUP2 MOD ROT OVER - ROT / SWAP \>>
I *never* realized that MOD would do the job there.
I did wonder, but never tested it.
Great - you have made an IDIV2 for the 48 :)
This is much better than my current version:
<< DUP IP 1 ->LIST SWAP E12 SWAP FP OVER *
WHILE DUP 0 >
REPEAT DUP 4 ROLLD IDIV2 3 ROLLD + 3 ROLLD
END * +
>>
This returns a zero on the last position (that last * + does
that), but I can see you will use DEPTH with your string to
figure out how far to go with it.
Interesting how you use the stack - I must expand my
imagination - it only extends to a 4 levels <G>
I really look forward to your next instalment :)
I will try your ->CF today.
It is a work of ART!!! I wondered about MANT as well.
Your MANT E11 * is probably better than FP E12 *
- then your OVER / creates a perfect E12 :)
Very interesting.
Is it true that I could take your text above, store it in a
file on the desktop and then transfer it to a HP48/49 and it
would work?
Could you point me to something at www.hpcalc.org that
explains this format, especially the "@" for comments ?
Joseph K. Horn
> I *never* realized that MOD would do the job there.
Is there any other way? Any other method (e.g. / IP first) gets
roundoff errors.
> Your MANT E11 * is probably better than FP E12 *
> - then your OVER / creates a perfect E12 :)
Because of MOD's accuracy, you don't need to turn the input into a
ratio of huge integers; just start with the unmodified input as the
numerator, 1 as the denominator, and start looping. Cool, huh? It's
*mathematically* identical to what you were doing, but faster and just
as accurate, thanks to HP's MOD code.
-Joe-
Joseph K. Horn
> AFAIK, unless you patent your algorithm for
> every possible way it may be used, anybody
> may use it. All you own at the moment is the
> program you wrote (that uses this algorithm)
> and I guess it could be named 'Horn Algorithm'
> if it's truly new. So if you don't want
> others to freely use your algorithm, patent it
> or keep it for yourself; Both being Bad Things
> anyway in my opinion.
I agree that keeping it to myself would be a Bad Thing, which is why I
(a) presented it at HHC 2003, (b) indirectly presented in this thread
and will post it in full shortly, and (c) took a vow never to do that
with anything. ;-)
As for patenting it, wouldn't it be a Good Thing IF I did that AND (a)
either gave or licensed it to HP for them do put it into their
calculators, and (b) released it for free use by everybody EXCEPT as a
built-in feature of any handheld device? That way, everybody would
benefit, even TI and their ilk could offer it as add-on software.
Only HP would be able to advertise "the best fraction calculators"
where "best" refers to the fractions, not to the calculators, of
course. >:-O
Since all this is currently being appraised by Powers That Be,
*please* advise if I'm blissfully strolling down a primrose path,
since I am not familiar with the legalities involved.
If it DOES get patended, it will be known as the PDQ Algorithm, not
the Horn Algorithm. As I told the HHC 2003 assemblage, "No way!
Students would see a button that says HORN, and they'd push it, and
complain that it doesn't work." The preceding is a printable
paraphrase of what was actually said. ;-)
I'll bet anybody here up to $1.53 that it is in fact a new algorithm.
This offer will stand until the first time I lose the bet, at which
time all bets are off. :-b
-Joe-
Joseph K. Horn
> I think one might say that "best" is in the eye of the user.
Hurwitz and most other mathematicians agree with you. In fact, my
full PDQ package outputs the "Hurwitz Accuracy" of each convergent
(calculated as the actual error times the denominator squared times
the square root of five) so that the user can base his or her decision
on that, if they wish.
Loosely, the Hurwitz Accuracy of a convergent is a kind of ratio of
the number of digits of input over the number of digits of output. If
making gears, for example, it would make no sense to go from a 22/7
gear ratio to a 179/57 gear ratio in an attempt to approximate pi,
because even though 179/57 is a "better fraction" for pi, it offers
only 0.0000227 increase in accuracy, roughly 2 centimeters per
kilometer! This can be seen easily by looking at the Hurwitz Accuracy
of 22/7 for pi (namely, -0.138546713972020) versus the Hurwitz
Accuracy of 179/57 (namely 9.021486720965779). Hurwitz proved that
any convergent with an Hurwitz Accuracy less than 1 (absolute value)
must be a continued fraction convergent, and not "merely" one of the
intermediate convergents.
In brief, the continued fraction convergents yield the "best"
fractions IF AND ONLY IF the "bang per buck" is considered, that is,
the number of digits you get out compared to the number of digits you
put in. In this sense, there is no "better" fraction than 355/113 for
pi until you get out to something 440 digits. The calculation of that
fraction is left as an exercise for the reader. ;-)
Hurwitz Accuracy is only one kind of accuracy. There are others, each
appropriate for some particular tasks. I discussed several in my
conference paper.
> What does "best" mean for you?
As I hope I explained during this thread already, it means different
things for different *purposes*. If one's task demands a complete and
accurate examination of ALL the fractions that approach any given
decimal number, then every fraction in that list is a "best" fraction,
and the PDQ Algorithm is the only way to find all of them.
On the other hand, if one has no need of precision or completeness,
then the age-old continued-fraction-only algorithm based on 1/FP(x)
will suffice, such as HP's ->Q.
> Joe's challenge is fascinating and difficult,
> but I'm not sure that the Best(JH) fraction
> is the "best" except in curiosity and
> educational value.
Thank you. I can think of no loftier objectives than those. ;-)
"Fun and learning. Everything comes down to those two."
-- Donald Shimoda, reluctant messiah.
PDQ solves a previously unsolved problem. Why? Because it was there,
and nobody else was trying to solve it. Amen.
> Can anybody think of a situation where
> 14823406/2469333 would be preferred to
> 15009499/2500333 as an approximation to
> 6.0030000004?
Yes! HP *claimed* to offer "the smallest fraction equal to the
displayed decimal number" with the ->Q function. It cannot read your
mind to know what YOU mean by "best", so it simply churned out the
smallest fraction equal to the input (in STD mode)... or at least it
was SUPPOSED to. PDQ does exactly what ->Q claimed to do but never
did.
> When would the slightly smaller numerator and
> denominator of 14823406/2469333 compared to
> 15009499/2500333 be of value, given the special
> property of 15009499/2500333 (said property to
> be determined by the audience).
Is "said property" that the latter has a much better Hurwitz Accuracy,
namely 0.7447095671868303496, as compared to -67.732319542344457 for
the former? If that's not the one you had in mind, it's still true.
> The ->Q function on the HP48 gives
> 15003496/2499333 for 6.00300000004
> (in STD mode), which lacks the special
> property that 15009499/2500333 has.
The former has a Hurwitz Accuracy of -1.491059671430066483, which also
not as good as the latter's (given above).
Something that Hurwitz seemed not to notice (or care about) is that
his definition of "best" omits not only all of the intermediate
convergents (as it seems you too would prefer to do) but also all the
regular ones based on a partial quotient that is smaller than the
largest one so far. This would mean that pi only has 3/1, 22/7, and
355/113 as "best" fractions, with the next "best" fraction having
something around 440 digits in it. An algorithm that spits out THOSE
fractions is trivial, and pretty useless considering the paucity of
fractions it generates. It doesn't find ANY for the square root of 3
(other than 2/1, which is not very satisfying!).
> This posting is not meant to belittle
> or minimize this latest excellent challenge
> of Joe Horn (who is a friend of mine). :-)
You, Bob Wilson, and Bob Hall have given me no end of delightful food
for thought for many years. This thread is my way of saying "Thank
You" to y'all!
"The best payment for a good idea is another good idea."
-- Richard J. Nelson
-Joe-
Tony Hutchins
Hi Joseph K. Horn ! <joe...@holyjoe.net>
02h31m ago, on Thu, 16.10.03 at 2:27 p.m. -0700, you wrote
in message ID <87233f9e.03101...@posting.google.com> :
> Tony Hutchins wrote:
>
> > I *never* realized that MOD would do the job there.
>
> Is there any other way? Any other method (e.g. / IP first) gets
> roundoff errors.
There was another way <VBG>
I did it my way with IP first.
This was my "IDIV2":
<< N D / IP DUP 0 > <<1 ->>IFT 'Q' STO
N D Q * - 'R' STO
WHILE R D >=
REPEAT
R D - 'R' STO
Q 1 + 'Q' STO
END
Q R
>>
If my IP was positive I took 1 away from it, and then, I gave
it back if it were needed. Hehe. Rather pointless I admit,
given that we have MOD ;-)
LOL! You should have seen my first version - it built up the
'Q', one by one, from zero, but it was a bit slow<G>.
> > Your MANT E11 * is probably better than FP E12 *
> > - then your OVER / creates a perfect E12 :)
>
> Because of MOD's accuracy, you don't need to turn the input into a
> ratio of huge integers; just start with the unmodified input as the
> numerator, 1 as the denominator, and start looping. Cool, huh? It's
> *mathematically* identical to what you were doing, but faster and just
> as accurate, thanks to HP's MOD code.
WOW is all I can say. The code shrinks<G>. IDIV2 is not needed
(well, can't be used then). MOD gives the remainder to carry
forward - meanwhile we need to reconstruct the integral quotient that
we capture: N D MOD gives R, and Q=(N/D -R)/D which will be an
exact integer, thanks to the MOD :)
Thanks! BTW I received my DataFile today and it has Richard's
HHC2003 conference report in, with a fascinating paragraph or
two about PDQ :)
--
Tony Hutchins
Wellington
New Zealand
#191 Whom God wishes to destroy, he first makes
mad. Proverb
Tony Hutchins
Hi Tony Hutchins ! <jus...@csi.com>
31m ago, on Fri, 17.10.03 at 1:33 p.m. +1300 (NZDT), you wrote
in message ID <1088903...@news.cis.dfn.de> :
[...]
> WOW is all I can say. The code shrinks<G>. IDIV2 is not needed
> (well, can't be used then). MOD gives the remainder to carry
> forward - meanwhile we need to reconstruct the integral quotient that
> we capture: N D MOD gives R, and Q=(N/D -R)/D which will be an
> exact integer, thanks to the MOD :)
Whoops! Q=(N-R)/D
--
Tony Hutchins
Wellington
New Zealand
#302 When you win, nothing hurts. Joe Namath
Werner Huysegoms
>
> > > Your MANT E11 * is probably better than FP E12 *
> > > - then your OVER / creates a perfect E12 :)
> >
> > Because of MOD's accuracy, you don't need to turn the input into a
> > ratio of huge integers; just start with the unmodified input as the
> > numerator, 1 as the denominator, and start looping. Cool, huh? It's
> > *mathematically* identical to what you were doing, but faster and just
> > as accurate, thanks to HP's MOD code.
>
[To Joe:] Ah yes of course!
> WOW is all I can say. The code shrinks<G>. IDIV2 is not needed
> (well, can't be used then). MOD gives the remainder to carry
> forward - meanwhile we need to reconstruct the integral quotient that
> we capture: N D MOD gives R, and Q=(N/D -R)/D which will be an
> exact integer, thanks to the MOD :)
>
[To Tony:] well, both of these combined is what IDIV2 does... so I keep it in.
Cheers,
Werner
Werner Huysegoms
.. but of course, I cannot call it IDIV2 any more as on the 49, IDIV2
will accept only integers. let's call it QUOTREM, then.
Werner
Rodger Rosenbaum
>Yes! HP *claimed* to offer "the smallest fraction equal to the
^^^^^!!
>displayed decimal number" with the ->Q function. It cannot read your
>mind to know what YOU mean by "best", so it simply churned out the
>smallest fraction equal to the input (in STD mode)... or at least it
>was SUPPOSED to. PDQ does exactly what ->Q claimed to do but never
>did.
>
Did HP actually say just as you have quoted? This a pretty bold statement; I can't find
it in any of my manuals. Where did they say that?
Gjermund Skailand
Finally I got the "evil" tests right.
Here is my entry, all test examples are correct.
The program builds the fraction one step at a time with cont. fractions, and
check lower and upper end if that is better. Test of upper end may be
unnecessary as it was not used in any of the test examples, and maybe I can
remove it.
The strange letter Z is the left-arrow, and the S is the >= .
On my old hp49g it runs in about 6 seconds even for the "evil" ones.
Best regards
Gjemund Skailand
main program
« 1. CF 2. CF DUP MANT
100000000000. * OVER XPON 11.
SWAP - ALOG DUP2 GCD ROT OVER
/ UNROT / ROT 1. 0. 0. 1. Q
R ZB ZQ1 ZR1 ZQ0 ZR0
«
DO Q R IDIV2 R 'Q' STO 'R'
STO DUP
IF 1. - DUP FTRUE
THEN OVER F * OVER F *
IF S
THEN 1. SF -1. BRKT F
ROT DUP
END
END DROP DUP
IF 1. + DUP FTRUE
THEN OVER F * OVER F *
IF S
THEN 2. SF 1. BRKT F
ROT DUP
END
END DROP F DUP2 ZR1
'ZR0' STO 'ZR1' STO ZQ1 'ZQ0'
STO 'ZQ1' STO
UNTIL / ZB == 1. FS? OR 2.
FS? OR
END ZQ1 ZR1
IF DUP2 / ZB == NOT
THEN DROP2 1.
9.99999999999E499
END
IF 2. FS?C
THEN BMIN
END
IF 1. FS?C
THEN BMIN
END
»
»
'R2F' STO
«
WHILE DUP2 + FTRUE
REPEAT 2. *
END OVER +
DO OVER 2. / OVER 2. / + IP
IF DUP FTRUE
THEN SWAP ROT
ELSE SWAP
END DROP
UNTIL DUP2 - ABS 1. ?
END DROP
»
'BRKT' STO
« DUP ZQ1 * ZQ0 + SWAP ZR1 *
ZR0 +
IFERR /
THEN DROP
END ZB ==
»
'FTRUE' STO
« DUP ZQ1 * ZQ0 + SWAP ZR1 *
ZR0 +
»
'F' STO
« DUP2 * 5. ROLL 5. ROLL DUP2
* 4. ROLL
IF <
THEN 4. ROLL 4. ROLL
END DROP2
»
'BMIN' STO
"Werner Huysegoms" <werner-h...@freegates.be> wrote in message
news:44ec85ff.03101...@posting.google.com...
Veli-Pekka Nousiainen
news:ka3vov815dehkde5g...@4ax.com...
Date, Time and Fraction Arithmetic 16-5
"The accuracy of the fractional approximation is dependent on the display
mode.
If the display mode is Dtd, the approximation is accurate to 11 significant
digits.
If the display mode is n Fix, the approximation is accurate to n significant
digits."
***Note word used: approximation
HP 48G Series Advanced User's Reference Manual
Command Reference 3-251
->Q To Quontient Command: Returns a rational from of the argument.
"The rational result is a "best guess", since there might be more than one
rational
expression consistent with the argument. ->Q finds a quontient of integers
that
agrees with the argument to within the number of decimal places specified by
the display format mode."
***Note the "quotas" by HP: "best guess"
VPN
Gjermund Skailand
> fractions IF AND ONLY IF the "bang per buck" is considered, that is,
> the number of digits you get out compared to the number of digits you
> put in. In this sense, there is no "better" fraction than 355/113 for
> pi until you get out to something 440 digits. The calculation of that
> fraction is left as an exercise for the reader. ;-)
1901870728566923076090143944714770339621590768313546337192526115562704339680
9635643200078081079293702997523451876888357413870030368533612856711580598677
02399073227994426905220194699766118756059055619036488502928002591
6053842551464203261023610232159405317163914781503450207392312531721347406882
3247694600005871377454979656144746826774641287402271754410094658714414873962
6803435133473281606663121381125761746030151344353855924025288111
With a total of 217+217 digits you get an accuary of pi of about 437 decimal
digits.
Best regards
Gjermund Skailand
( sorry couldn't resist)
Tony Hutchins
Hi Werner Huysegoms ! <werner-h...@freegates.be>
02h45m ago, on Thu, 16.10.03 at 10:50 p.m. -0700, you wrote
in message ID <44ec85ff.0310...@posting.google.com> :
> (follow up to the previous message, posted too quickly)
[I am always doing that]
> .. but of course, I cannot call it IDIV2 any more as on
> the 49, IDIV2 will accept only integers. let's call it
> QUOTREM, then.
Yep. Now your ->CF will be very short. I tried to incorporate
your QUOTREM in it:
<< DEPTH ->D
<< 1
WHILE DUP 0 >
REPEAT
SWAP OVER DUP2 MOD
ROT OVER - ROT / ROT ROT
END
DROP2
D DEPTH SWAP - ->LIST >> >>
I did this to generate either convergent denominators or
numerators depending on the little string input
CF-string little-string ->CONV
with little string = {0 1} we get the denominators
with it equal to {1 a0} we get the numerators.
<< DUP SIZE -> p a k
<< 2 k FOR i
a i GET p i GET * p i 1 - GET +
p SWAP + 'p' STO
NEXT
p 2 k OVER + SUB >> >>
I look forward to your next instalment Werner!
--
Tony Hutchins
Wellington
New Zealand
#177 A silly remark can be made in Latin as well as
in Spanish. Miguel de Cervantes
Tony Hutchins
Hi Rodger Rosenbaum ! <rodg...@siteconnect.com>
02h21m ago, on Thu, 16.10.03 at 11:34 p.m. -0700, you wrote
in message ID <ka3vov815dehkde5g...@4ax.com> :
> On 16 Oct 2003 16:23:33 -0700, joe...@holyjoe.net (Joseph K. Horn) wrote:
>
> >Yes! HP *claimed* to offer "the smallest fraction equal to the
> ^^^^^!!
> >displayed decimal number" with the ->Q function.
[...]
> Did HP actually say just as you have quoted? This a pretty bold statement; I can't find
> it in any of my manuals. Where did they say that?
eg HP48SX Programmer's Reference Manual, " ->Q finds a
quotient of integers that agrees with the given number to the
number of decimal places specified by the display format mode"
"agrees with" probably has to imply "equals" here
--
Tony Hutchins
Wellington
New Zealand
#56 A clean desk is a sign of a cluttered desk drawer.
Joseph K. Horn
> Did HP actually say just as you have quoted?
Bill Wickes (the Father of RPL) wrote in his book "HP48 Insights, Part
II: Problem Solving Resources" on page 526:
------- begin Wickes quotation -------
->Q finds the fraction with the *smallest denominator* that matches
the argument within the error limit set by the display format. For
example,
STD 1.61803398875 ->Q --> '514229/317811'
Notice that the denominator does not have 11 digits. The fraction
139583862445/86267571272 is another representation of the original
argument that evaluates to the same decimal value as 514229/317811;
->Q chooses the latter because the denominator is smaller.
------- end of Wickes quotation -------
(emphasis is Wickes')
Also, the author of the ->Q function himself (not Bill Wickes) wrote,
"->Q's aim is to produce the fraction with the *smallest denominator*
and a given error size" [emphasis is his] in the documentation for
NEW2Q on Goodies Disk #3, which was his response to my 1990 DEC2FRAC
program, my first attempt to improve ->Q.
If both he and Bill Wickes clearly stated that ->Q was supposed to get
the smallest possible denominator, you might ask why the official HP
documentation only go so far as to call it a "best guess." [Cf. the
HP 48SX Owner's Manual (page 134), and the HP 48 Owner's Manual (page
9-5), and the HP 48 Programmer's Reference Manual (page 283), and the
HP 48G Series Advanced User's Reference Manual (page 3-251)]. The
reason can also be found in the NEW2Q documentation:
"The HP48 manuals did little to clarify this issue, especially since
we wanted to make no specific claim with only a conjecture of
correctness."
BINGO! There it is! What HP *intended* to do is clear, and the above
sentence reveals that they never really wrapped their brains around
the problem, and lived in hope that it worked well enough to minimize
complaints, but couldn't guarantee the *smallest* denominator. ->Q
never lived up to their original stated goal, like PDQ does.
Gjermund Skailand
All test examples are now calulated within 5 sec. on the old hp49g.
Directory size is now 632 bytes.
Best regards
Gjermund Skailand
------------------------------------
@ Main program: continous fraction + check if accuracy is exact equal
/<< DUP MANT 100000000000. *
OVER XPON 11. SWAP - ALOG ROT
1. 0. 0. 1. /-> /<-B /<-Q1 /<-R1 /<-Q0
/<-R0
/<<
DO DUP UNROT IDIV2 SWAP
IF DUP FTRUE
THEN TTTTT
ELSE F /<-R1 '/<-R0' STO
'/<-R1' STO /<-Q1 '/<-Q0' STO '/<-Q1'
STO 0.
END
UNTIL
END UNROT DROP2 F
/>>
/>>
'R2F' STO
@ subprogram: Pick the fraction with smallest size
/<< DUP F * OVER
IF 1. - DUP FTRUE
THEN DUP F * ROT
IF /<=
THEN -1. BRKT SWAP DROP
ELSE DROP
END
ELSE DROP2
END 1
/>>
'TTTTT' STO
@ subprogram FTRUE, return 1 if fraction is exact
/<< DUP /<-Q1 * /<-Q0 + SWAP /<-R1 *
/<-R0 + / /<-B ==
/>>
'FTRUE' STO
@ subprogram F: calculate fraction
/<< DUP /<-Q1 * /<-Q0 + SWAP /<-R1 *
/<-R0 +
/>>
'F' STO
@ subprogram BRKT: search for limit of continous fraction interval resulting
in
@ same (exact) result (bisesction?)
/<<
WHILE DUP2 + FTRUE
REPEAT
@ while still exact, expand interval until it includes oundary
2. *
END OVER +
@ determine boundary, using midpoint at each iter.
DO OVER 2. / OVER 2. / + IP
IF DUP FTRUE
THEN SWAP ROT
ELSE SWAP
END DROP
UNTIL DUP2 - ABS 1. /<=
END DROP
/>>
'BRKT' STO
Wayne Brown
> Hurwitz and most other mathematicians agree with you. In fact, my
> full PDQ package outputs the "Hurwitz Accuracy" of each convergent
> (calculated as the actual error times the denominator squared times
> the square root of five) so that the user can base his or her decision
> on that, if they wish.
Will you be making the "full PDQ package" available somewhere (like
maybe hpcalc.org)? I poked around a bit on your web site last night
but couldn't find it.
--
Wayne Brown (HPCC #1104) | "When your tail's in a crack, you improvise
fwb...@bellsouth.net | if you're good enough. Otherwise you give
| your pelt to the trapper."
"e^(i*pi) = -1" -- Euler | -- John Myers Myers, "Silverlock"
Werner Huysegoms
@ QUOTREM - Quotient and Remainder
@
@ bytes: 30.0
@ check: #C142h (48)
@ #14A6h (49)
@
@ In: 2: a
@ 1: b
@
@ Out: 2: q
@ 1: r
@
@ so that b*q + r = a
@
\<< DUP2 MOD ROT OVER - ROT / SWAP \>>
@ \->CF - Create Simple Continued Fraction
@
@ bytes: 74.5
@ check: #25D3h (48)
@ #6FABh (49)
@
@ In: 1: x
@
@ Out: 1: { a0 a1 .. an }, ai integer (any i), and ai>0 for i>0
@
@ x = a0 + 1/(a1 + 1/(.. + 1/an))
@
\<<
DEPTH
\-> D
\<<
@ turn x into p/q - start with p=x and q=1
1. SWAP
@ stack holds 2:q 1:p, with ri = p/q
WHILE
@ calculate ai
OVER QUOTREM DUP
REPEAT @ stop when remainder is zero
@ calculate ri+1
ROT
END
ROT DROP2
D DEPTH SWAP - \->LIST
\>>
\>>
@ SAF - determine Smallest Approximating Fraction
@
@ bytes: 237.0
@ check: #7A2Ah (48)
@ #AE9Fh (49)
@
@ In: 1: x
@ Out: 2: p
@ 1: q
@
@ p,q are integers
@ a simple / will yield the exact argument back.
@
@ timing:
@ .500000000001 48SX: 3.6084_s
@ 48GX: 2.4226_s
@ 49G: 2.649_s
@ 3.14159265359 48SX: 1.3682_s
@ 48GX: 0.8965_s
@ 49G: 1.0372_s
@ timings obtained with EMU48
\<<
@ determine Simple Continued Fraction
DUP \->CF 0.
\-> z cf i
\<<
@ *forward* evaluate the CF till the quotient matches the argument
@
@ P[i] = ai*P[i-1] + P[i-2]
@ Q[i] = ai*Q[i-1] + Q[i-2]
@
@ where
@ | -2 -1 0
@ -------------
@ P | 0 1 a0
@ Q | 1 0 1
@
@ we'll use a complex number to hold a (p,q) pair
(0. 1.) (1. 0.)
DO
DUP cf 'i' INCR GET * ROT +
UNTIL DUP C\->R / z SAME
END
SWAP cf i GET 0.
@ bisection to determine largest N for which
@
@ P[i+1] - N*P[i]
@ EVAL(---------------) = z
@ Q[i+1] - N*Q[i]
@
@ stack: R[i+1] R[i] a[i+1] 0.
@ Ri = (Pi Qi)
WHILE DUP2 - 1. >
REPEAT
4. DUPN + 2. / IP
DUP 5 ROLLD
* - C\->R / z \=/
{ ROT } IFT DROP
END
SWAP DROP
* - C\->R
\>>
\>>
Werner
Werner Huysegoms
> Sam Hughes writes:
>
> >> For what it's worth, my program converts
> >> sqrt(LOG(2)) into 998995/1820784 in 0.745
> >> seconds on an hp49g+, and runs unmodified on
> >> an HP48GX.
> >
> > How long does it take to run on a 48GX,
> > or on a 49G?
>
> It takes approximately 3 times longer, for almost any input.
>
> -Joe-
My SAP program (posted somewhere in this thread a few minutes ago)
takes 1.31_sec on a 48GX for this problem.
That certainly beats 3*0.745, even when calculated with a TI ;-)
And I forgot to mention that it converts all your examples flawlessly,
which doesn't mean the program is without errors ;-O
Cheers, Werner
Gilles Carpentier
Cool ! A FX502P web site... :-) My first programmable calculator ! I've not
yet finished but it's already on line :
Of course, i will put a link toward your 502 site.
Gilles
"FX-502P Geek" <use_th...@mysite.co.uk> a écrit dans le message news:
H0_hb.1364$ny5...@news-binary.blueyonder.co.uk...
> I do not deny that the problem is very difficult - I tried it on the Casio
> FX-502P back in the 80s (except I preferred the output in vulgar
fractions),
> but by simple inspection of your challenge you've got the wrong answer for
> the 5th example. The fractional part of 50000.0005 is clearly 1/20000 and
so
> your correct answer should be 1000000001/20000 as given by ->Q of a HP48.
> Incidentally the program that I wrote for the FX-502P is available at
> www.fx502p.pwp.blueyonder.co.uk and gives the correct answer of
> 50000+1/20000. I would not pretend however that the program would be up to
> the challenge in general. I would of course like to see your HP48
program...
> especially when it is fixed.
>
> Peter.
Tony Hutchins
Hi Werner Huysegoms ! <werner-h...@freegates.be>
05h27m ago, on Fri, 17.10.03 at 07:01 a.m. -0700, you wrote
in message ID <44ec85ff.0310...@posting.google.com> :
> second 'instalment'. I'm not too happy with 'SAF' yet, though.
But you *must* be happy with this first part!!!!!
> @ determine Simple Continued Fraction
> DUP \->CF 0.
> \-> z cf i
> \<<
> (0. 1.) (1. 0.)
> DO
> DUP cf 'i' INCR GET * ROT +
> UNTIL DUP C\->R / z SAME
> END
I racked my brain trying to figure out an elegant way to do
that - this complex number usage is so ingenious. Thanks!
BTW your timings *do* indicate you have approached Joe's speed
on this!! I think he got about .8 seconds for the
.500000000001 on the 49g+? That would probably be about 3
seconds on the 49G, and you have 2.649!!
Also, FWIW, I think your second part, that does the bisection,
is certainly very efficient, as it uses the stack only, and it
has to be the way it is to reflect the comlexity of the
calculation itself. For my money, you have won the challenge
Werner!!!
Overall I suppose one could possibly achieve something faster
by not formally computing all of the partial quotients, as
they are not needed past a certain point. But, that would
entail doing too many things at once, for my brain anyway<G>
I wonder if Joe uses your complex number method?
Now I see it, it seems so natural for doing 2 things at once
:)
--
Tony Hutchins
Wellington
New Zealand
#8 Quantized Revision of Murphy's Law: Everything goes wrong all at once.
Werner Huysegoms
Tony Hutchins <jus...@csi.com> wrote in message
news:1060035...@news.cis.dfn.de...
> -=[ Sat, 18.10.03 08:28 a.m. +1300 (NZDT) ]=-
>
> Hi Werner Huysegoms ! <werner-h...@freegates.be>
Hey, Werner will do!
> BTW your timings *do* indicate you have approached Joe's speed
> on this!! I think he got about .8 seconds for the
> .500000000001 on the 49g+? That would probably be about 3
> seconds on the 49G, and you have 2.649!!
>
IIRC, that was for SQRT(LOG(2.)), for which I get 1.5 secs (on a real 49G -
shall
we call them 49G- from now on?)
> Also, FWIW, I think your second part, that does the bisection,
> is certainly very efficient, as it uses the stack only, and it
> has to be the way it is to reflect the comlexity of the
> calculation itself. For my money, you have won the challenge
> Werner!!!
Nono -*you* have won it.
All I did was copycat your work and write it the way I am used to, and I do
have a
fair bit of 'mini-challenge experience' ;-)
>
> Overall I suppose one could possibly achieve something faster
> by not formally computing all of the partial quotients, as
> they are not needed past a certain point. But, that would
> entail doing too many things at once, for my brain anyway<G>
For Pi (well its 12-digit approximation) that would be a bit faster.
But the bisecting part is the most time-consuming, and I'm too lazy to
combine
the 'continued fraction' part and the 'forward evaluation part'.
They are pretty elegant as they are.
>
> I wonder if Joe uses your complex number method?
First, I used a list - but that doesn't work an a 48S.
Then I used a vector, and only thought about complex numbers by accident..
It *halves* the time.
> Now I see it, it seems so natural for doing 2 things at once
My feeling exactly - when Joe pointed out the superfluous x -> p/q
conversion at the beginning
BTW:
In SAF, replace (0,1) (1,0) by -1. \v/ 1. (-1. SQRT 1.), 29.5 bytes less,
and it
works the same way even if the second number is real i.o. complex
Werner
Rodger Rosenbaum
>Rodger Rosenbaum writes:
>
>> Did HP actually say just as you have quoted?
>
>Bill Wickes (the Father of RPL) wrote in his book "HP48 Insights, Part
-------SNIP
>BINGO! There it is! What HP *intended* to do is clear, and the above
>sentence reveals that they never really wrapped their brains around
>the problem, and lived in hope that it worked well enough to minimize
>complaints, but couldn't guarantee the *smallest* denominator. ->Q
>never lived up to their original stated goal, like PDQ does.
>
So I take it that the answer to my question:
-----------
>Yes! HP *claimed* to offer "the smallest fraction equal to the
^^^^^!!
>displayed decimal number" with the ->Q function. It cannot read your
>mind to know what YOU mean by "best", so it simply churned out the
>smallest fraction equal to the input (in STD mode)... or at least it
>was SUPPOSED to. PDQ does exactly what ->Q claimed to do but never
>did.
>
Did HP actually say just as you have quoted? This a pretty bold statement; I can't find
it in any of my manuals. Where did they say that?
-----------
is no?
The examples you gave don't use the unqualified word "equal". They use phrases like
"given error size", "within the error limit", "best guess", etc.
Tony Hutchins
Hi Werner
01h16m ago, on Fri, 17.10.03 at 10:37 p.m. +0200, you wrote
in message ID <bmpjm0$cf$1...@news.worldonline.be> :
> Hi, Tony.
>
> > Hi Werner Huysegoms ! <werner-h...@freegates.be>
>
> Hey, Werner will do!
That's my automated header program - I should change the
template so it doesn't reproduce the date as well - too lazy!
> > BTW your timings *do* indicate you have approached Joe's speed
> > on this!! I think he got about .8 seconds for the
> > .500000000001 on the 49g+? That would probably be about 3
> > seconds on the 49G, and you have 2.649!!
> >
>
> IIRC, that was for SQRT(LOG(2.)), for which I get 1.5 secs
Ah, yes, that rings a bell.
Wow, you may even have a faster version than Joe's then!!
Certainly very close.
> (on a real 49G - shall we call them 49G- from now on?)
Hehe. Maybe they are just a complex pair, the real 49G and the
imaginary 49g+ :)
> > Also, FWIW, I think your second part, that does the bisection,
> > is certainly very efficient, as it uses the stack only, and it
> > has to be the way it is to reflect the comlexity of the
> > calculation itself. For my money, you have won the challenge
> > Werner!!!
>
> Nono -*you* have won it.
OKok, I accept :)
> All I did was copycat your work and write it the way I
> am used to, and I do have a fair bit of 'mini-challenge
> experience' ;-)
Ahha, I am not surprised. I have learned a lot!
[...]
> For Pi (well its 12-digit approximation) that would be a bit
> faster. But the bisecting part is the most time-consuming,
> and I'm too lazy to combine the 'continued fraction' part
> and the 'forward evaluation part'. They are pretty elegant
> as they are.
Yes they are, and it's nice to have the full partial quotient
string as a by product. I can't think of any other way to do
the final part.
> > I wonder if Joe uses your complex number method?
>
> First, I used a list - but that doesn't work an a 48S.
Why's that? Does a 48G have extra LIST features?
Couldn't you use LIST-> just like you use C->R?
> Then I used a vector, and only thought about complex numbers by accident..
> It *halves* the time.
WOW!!
> > Now I see it, it seems so natural for doing 2 things at once
>
> My feeling exactly - when Joe pointed out the superfluous x -> p/q
> conversion at the beginning
My first reaction is to avoid MOD since bad experiences with
it years ago with Basic. Old habits die hard, especially when
a race is on :)
> BTW:
> In SAF, replace (0,1) (1,0) by -1. \v/ 1. (-1. SQRT 1.), 29.5 bytes less,
Hehe. You have me there. You see I am way behind the times.
I have never coded on an emulator and don't know what \v/ is
for a start. I haven't a clue as to how 2 complex numbers
spring out here. I know that "\" prefixes a command not easily
displayable in regular ascii, and I just guess most of them,
when I see them here, like \->. I do see how the -1 SQRT gives
(0,1), but where does the (1,0) come from?
> and it works the same way even if the second number is real
> i.o. complex
Which second number?
--
Tony Hutchins
Wellington
New Zealand
#103 Life is tons of discipline. Robert Frost
James Horn
> If it DOES get patended, it will be known as the PDQ Algorithm, not
> the Horn Algorithm. As I told the HHC 2003 assemblage, "No way!
> Students would see a button that says HORN, and they'd push it, and
> complain that it doesn't work." The preceding is a printable
> paraphrase of what was actually said. ;-)
> I'll bet anybody here up to $1.53 that it is in fact a new algorithm.
> This offer will stand until the first time I lose the bet, at which
> time all bets are off. :-b
What a trooper! Thanks for a truely fascinating and fun thread, Joe and
friends!!!
JimH (older'n Joe)
BlackLotus
"Tony Hutchins" <jus...@csi.com> a écrit dans le message de
news:1075661...@news.cis.dfn.de...
> -=[ Sat, 18.10.03 10:53 a.m. +1300 (NZDT) ]=-
>
*snip*
>
> Hehe. You have me there. You see I am way behind the times.
> I have never coded on an emulator and don't know what \v/ is
> for a start. I haven't a clue as to how 2 complex numbers
> spring out here. I know that "\" prefixes a command not easily
> displayable in regular ascii, and I just guess most of them,
> when I see them here, like \->. I do see how the -1 SQRT gives
> (0,1), but where does the (1,0) come from?
I think that's what Werner was saying, it works even if the second number
(1) is a real and not its complex form ((1,0)): in terms of how the calc
handles it (in Werner's prog), having (0,1) 1 on the stack works the same as
(0,1) (1,0)
>
> > and it works the same way even if the second number is real
> > i.o. complex
>
> Which second number?
See above.
Tony Hutchins
Hi BlackLotus ! <black...@nospam.me>
27m ago, on Sat, 18.10.03 at 07:55 a.m. +0200, you wrote
in message ID <bmqkkq$8ml$1...@news-reader1.wanadoo.fr> :
> Heya
[...]
> I think that's what Werner was saying, it works even if the second number
> (1) is a real and not its complex form ((1,0)): in terms of how the calc
> handles it (in Werner's prog), having (0,1) 1 on the stack works the same as
> (0,1) (1,0)
Thanks Blacklotus :)
So, Werner explained *everything* to me in the parentheses!
Now I see that \v/ is SQRT!!
I thought the (-1 SQRT 1) was part of the code.(I even tried
to execute it<G>).
Thanks for putting me straight :)
--
Tony Hutchins
Wellington
New Zealand
#207 We are all worms, but I do believe that I am a
glowworm. Winston Churchhill
Rodger Rosenbaum
It looks to me like there is an earlier convergent with small Hurwitz accuracy.
22/7 gives a Hurwitz accuracy of .0076167
The continued fraction of 307 partial quotients has a Hurwitz accuracy of .00511777
and gives a convergent of 474490.../151034... (sorry; I'm not gonna type in all those
digits) for an absolute accuracy of about 311 digits of PI.
Your example is a continued fraction of 431 partial quotients with a Hurwitz accuracy of
.00010762
Joseph K. Horn
> Will you be making the "full PDQ package"
> available somewhere (like maybe hpcalc.org)?
> I poked around a bit on your web site last
> night but couldn't find it.
Real Soon Now; it's up to the Powers That Be. It'll first appear
here. However, now that the fundamental concepts have been revealed
in this thread, it'll be kinda anticlimactic. Yes, I'll make it
available on my website at that time.
-Joe-
Joseph K. Horn
> My SAP program takes 1.31_sec on a 48GX
> for this problem. That certainly beats
> 3*0.745, even when calculated with a TI ;-)
Your program is undoubtedly faster than mine at this stage. Mine
still uses variables instead of the stack (locals though, not globals;
they're too slow), and it also inputs the desired accuracy and tests
against that, and it keeps track of the Hurwitz ratio along the way,
all of which slows it down. But for the core task of this Challenge,
you're WAY ahead of me! :-)
-Joe-
Tony Hutchins
Hi Werner
12h ago, on Fri, 17.10.03 at 10:37 p.m. +0200, you wrote
in message ID <bmpjm0$cf$1...@news.worldonline.be> :
> Nono -*you* have won it.
> All I did was copycat your work and write it the way I
> am used to, and I do have a fair bit of 'mini-challenge
> experience' ;-)
There should be a library of mini-challenges - must be some
great coding tips in those answers.
And you did much more than copycat - you gave me many puzzles
eg - I wondered about
@ stack: R[i+1] R[i] a[i+1] 0.
where you use a[i+1] as the maximum integral N - it took me a
while to confirm this. IMHO your bisection code is the BEST :)
I never quite realised how to use the stack this way.
[...]
> BTW:
> In SAF, replace (0,1) (1,0) by -1. \v/ 1. (-1. SQRT 1.), 29.5 bytes less,
> and it
> works the same way even if the second number is real i.o. complex
Big byte saving! Plus, I understand this now too :)
--
Tony Hutchins
Wellington
New Zealand
#338 If a man has a talent and cannot use it,
he has failed. Thomas Wolfe *
Joseph K. Horn
> So I take it that the answer to my question is no?
> The examples you gave don't use the unqualified word
> "equal". They use phrases like "given error size",
> "within the error limit", "best guess", etc.
(A) SEMANTIC ARGUMENT: Hey, you left out Wickes' use of the word
"matches"! He wrote: "->Q finds the fraction with the *smallest
denominator* that MATCHES the argument within the error limit set by
the display format." (emphasis mine)
If "matches within the error limit", when the error limit is zero,
does not mean "equals" to you, then the answer to your question is
"no", and I apologize for misquoting from memory. However, as far as
I can tell, it means "equals". Perhaps we should ask him. Anybody
know where he is these days?
'A+B+C' { 'B' 'X' } ^MATCH --> 'A+B+X' 0.5 :-b
(B) PHILOSOPHICAL ARGUMENT: If ->Q in STD mode is NOT supposed to be
"EQUAL" to the input, then what IS it supposed to be? "Close"? HOW
close?!? Within precisely what range? If that range is unclear, then
of what use is the doggone thing? "Take this number and give me
something more or less close to it"?!? Come on now! How could one
who scrutinizes the precise difference between "equals" and "matches
with zero error" tolerate such vague thinking? ;-)
I am merely suggesting that ->Q be rewritten to have a user-specified
error range, from exact to Whatever The User Needs, in such a way that
it can be guaranteed to give a result in precisely that range, and not
just ANY result in that range, but EXACTLY the CLOSEST POSSIBLE one.
If that's too threatening or restrictive or something, we could just
*add* a new PDQ function, leaving the current ->Q as-is in the
machine, for those who prefer fuzzy mathematics. That way everybody
is happy.
(C) ARGUMENT FROM AUTHORITY: Because God said so. Oops, that only
works in theology. Sorry. ;-)
-Joe-
Joseph K. Horn
> > there is no "better" fraction than 355/113 for
> > pi until you get out to something 440 digits.
> > The calculation of that fraction is left as
> > an exercise for the reader. ;-)
> That would be (using my longfloat library): [snip]
Very good! That's the one I was thinking of. Unfortunately, I missed
one in between 355/113 and yours, namely:
474490238234227174515948361693811988401340891649950114502451~
856183877768471768408541740671737665969442224795157401432738~
79768480709913455315809719544718248
------------------------------------------------------------
151034933727656572741457537643244123466736157365649904806452~
295026919400951667087834668197436704093923364340930317605883~
86010342390699617474201564040715683
This uses 310 digits (total) to generate pi to 311 digits of
accuracy... slightly better than 355/113, and not quite as good as
yours.
It is worthy of note that HP's version of pi (3.14159265359) has no
better fraction than 355/113 if Hurwitz accuracy is the only
consideration:
Input: 3.14159265359
P.Q. Places Hurwitz Fraction
------ ------------------- ------------------- --------
7 0.848959279048027 0.499474276880018 3//1
15 2.898084852551616 0.858403770355093 22//7
292 6.573872807978566 2.118230918843718 355//113
(Do space-tabbed tables still look ok in newsgroups like they used
to?)
-Joe-
Rodger Rosenbaum
>Rodger Rosenbaum writes:
>
>> So I take it that the answer to my question is no?
>> The examples you gave don't use the unqualified word
>> "equal". They use phrases like "given error size",
>> "within the error limit", "best guess", etc.
>
>(A) SEMANTIC ARGUMENT: Hey, you left out Wickes' use of the word
>"matches"! He wrote: "->Q finds the fraction with the *smallest
>denominator* that MATCHES the argument within the error limit set by
>the display format." (emphasis mine)
>
>If "matches within the error limit", when the error limit is zero,
Where does he say that the error limit is zero in std mode? Even you aren't able to
find a small fraction whose value is "equal" to the input argument except in a relatively
small number of cases.
The HP48 AUR says, "..agrees (synonymous with "matches") with the argument to within
the number of decimal places specified by the display format mode", so we might expect 12
digits of agreement in Std mode, except that: The HP48 User's Guide says on page 16-5 near
the bottom of the page, "If the display mode is Std, the approximation is accurate to 11
significant digits."
>does not mean "equals" to you, then the answer to your question is
>"no", and I apologize for misquoting from memory. However, as far as
>I can tell, it means "equals". Perhaps we should ask him. Anybody
>know where he is these days?
>
>'A+B+C' { 'B' 'X' } ^MATCH --> 'A+B+X' 0.5 :-b
>
>(B) PHILOSOPHICAL ARGUMENT: If ->Q in STD mode is NOT supposed to be
>"EQUAL" to the input, then what IS it supposed to be?
The AUR says "best guess"
"Close"? HOW
>close?!?
The AUR says "..within the number of decimal places specified by the display format
mode"
Within precisely what range? If that range is unclear, then
>of what use is the doggone thing? "Take this number and give me
>something more or less close to it"?!? Come on now! How could one
>who scrutinizes the precise difference between "equals" and "matches
>with zero error"]
Wickes didn't say "matches with zero error" according to your quote above; he said
"within the error limit set by the display format", which as I said above would be, at
best, 12 digits of agreement. But the User's Guide cops out at 11 digits in Std mode.
tolerate such vague thinking? ;-)
Is it vague thinking to put "matches with zero error" in quotes, thereby suggesting
that this is what Wickes actually said, when the quotation we see above is "matches the
argument within the error limit set by the display format"? :-)
>
>I am merely suggesting that ->Q be rewritten to have a user-specified
>error range, from exact to Whatever The User Needs,
Hasn't HP provided this, in that the accuracy is controlled by the display mode?
Perhaps not as fine control as might be, but user-specified nonetheless.
Veli-Pekka Nousiainen
X
> (Do space-tabbed tables still look ok in newsgroups like they used
> to?)
BUT
is it not time to bring out documented programs
so that we all can admire the different approaches here
???
VPN
Tony Hutchins
Hi Rodger Rosenbaum ! <rodg...@siteconnect.com>
06h49m ago, on Sat, 18.10.03 at 06:16 a.m. -0700, you wrote
in message ID <m2d2pvor8kupc24oe...@4ax.com> :
[...]
> Even you aren't able to find a small fraction whose value is
> "equal" to the input argument except in a relatively small
> number of cases.
Joe's method as presented here does find a rational fraction
that equals the input in the sense that if you subtract output
from the input you get zero. In *every* case.
I am obviously not sure what you mean by '"equals"', or by
"small fraction".
Can you give an example of one input from the large number of
cases you have in mind, and approximately what proportion of
say 100 random inputs would fall into this class?
--
Tony Hutchins
New Zealand
Tony Hutchins
Hi Tony Hutchins ! <jus...@csi.com>
01h06m ago, on Sun, 19.10.03 at 09:32 a.m. +1300 (NZDT), you wrote
in message ID <1062568...@news.cis.dfn.de> :
> I am obviously not sure what you mean by '"equals"', or by
> "small fraction".
>
> Can you give an example of one input from the large number of
> cases you have in mind, and approximately what proportion of
> say 100 random inputs would fall into this class?
Oh! You must mean inputs like 1,2,3 etc which can only result
in 1/1,2/1,3/1... or say 3.14159 which only results in
314159/100000 - these are not on the "small fraction" class?
And, by "equals" you mean that if you subtract the input from
the output you get zero.
So, is your point not purely semantic?
One must allow Joe some poetic license ;-)
And, in the examples for Joe's challenge, all inputs and final
outputs both look equal, and are also indeed equal by your
stronger test. This is because they al have 12 significant
digits :)
--
Tony Hutchins
Wellington
New Zealand
#290 Nothing can be loved or hated unless it is
first known. Leonardo Da Vinci
Tony Hutchins
Hi Werner
1 day 08h04m ago, on Fri,
17.10.03 at 07:01 a.m. -0700, you wrote
in message ID <44ec85ff.0310...@posting.google.com> :
> @ SAF - determine Smallest Approximating Fraction
> @
> @ bytes: 237.0
[...]
> WHILE DUP2 - 1. >
> REPEAT
> 4. DUPN + 2. / IP
> DUP 5 ROLLD
> * - C\->R / z \=/
> { ROT } IFT DROP
> END
I've never seen an IFT with the {} used like that!
Saves 7.5 bytes over the << >>. Is there a rule for using {}
instead of << >> ?
I confirm on my old HP48SX, with 3/4 rubber feet<G>, that the
SAF of only 207.5 bytes (with the -1 SQRT 1 instead of (0 1)
(1 0) ) works like a charm!
--
Tony Hutchins
Wellington
New Zealand
#18 PATH=C:\WIN\HELP\NEEDHELP\MOREHELP\OH_NO\ARRRGGGH\KILL\DIE_WIN\;
Veli-Pekka Nousiainen
> -=[ Sun, 19.10.03 11:06 a.m. +1300 (NZDT) ]=-
>
> Hi Werner
>
> 1 day 08h04m ago, on Fri,
> 17.10.03 at 07:01 a.m. -0700, you wrote
> in message ID <44ec85ff.0310...@posting.google.com> :
>
> > @ SAF - determine Smallest Approximating Fraction
> > @
> > @ bytes: 237.0
> [...]
>
> > WHILE DUP2 - 1. >
> > REPEAT
> > 4. DUPN + 2. / IP
> > DUP 5 ROLLD
> > * - C\->R / z \=/
> > { ROT } IFT DROP
> > END
>
> I've never seen an IFT with the {} used like that!
> Saves 7.5 bytes over the << >>. Is there a rule for using {}
> instead of << >> ?
::ROT IFT
how many more nibbles do you save now?
BTW: I think that standardr IF is faster...
VPN
Tony Hutchins
Hi Veli-Pekka Nousiainen ! <DROP...@welho.com>
01h09m ago, on Sun, 19.10.03 at 01:19 a.m. +0300, you wrote
in message ID <bmse97$7db$1...@nyytiset.pp.htv.fi> :
[...]
> > > * - C\->R / z \=/
> > > { ROT } IFT DROP
> > > END
> >
> > I've never seen an IFT with the {} used like that!
> > Saves 7.5 bytes over the << >>. Is there a rule for using {}
> > instead of << >> ?
> How about testing with
> ::ROT IFT
> how many more nibbles do you save now?
another 1.5 :)
> BTW: I think that standardr IF is faster...
Ah, yes IF THEN ROT END takes the same size as { ROT } IFT
too!! << ROT >> IFT takes 7.5 bytes more.
Thanks VPN!
Rodger Rosenbaum
> -=[ Sun, 19.10.03 09:05 a.m. +1300 (NZDT) ]=-
>
>Hi Rodger Rosenbaum ! <rodg...@siteconnect.com>
>
>06h49m ago, on Sat, 18.10.03 at 06:16 a.m. -0700, you wrote
>in message ID <m2d2pvor8kupc24oe...@4ax.com> :
>
>[...]
>> Even you aren't able to find a small fraction whose value is
>> "equal" to the input argument except in a relatively small
>> number of cases.
>
>Joe's method as presented here does find a rational fraction
>that equals the input in the sense that if you subtract output
>from the input you get zero.
Only if you do the subtraction with finite precision arithmetic, as on a HP calculator.
But this doesn't mean that different fractions that seem to have the same value when
converted to a decimal on the HP48 aren't possessed of different relative errors, computed
using exact arithmetic, of course. The Hurwitz accuracy seems to be from an era when
gears were used to divide quantities; it was good to have gears without too many teeth.
Now, we do these things electronically and we don't need a small number of teeth on a
gear. Now, we might prefer the lowest relative error.
In *every* case.
>
>I am obviously not sure what you mean by '"equals"'
Let me give some definitions.
By "equal" I mean identical in numeric value. I do not mean equal to within 1 part in
10^12, which is a limitation imposed by the HP48's floating point system.
By "trivial case" I mean that given a decimal number on the HP48 such as x.xxxxxxxxxxx it
can always be represented exactly by a fraction xxxxxxxxxxxx/100000000000 (or its reduced
version), but this is not very interesting.
By "relative error" I mean, given a number N and an approximation to it M, the relative
error of M is given by (M-N)/N. This must be computed without error induced by a finite
precision floating point number system.
, or by
>"small fraction".
By "small fraction" I mean a fraction with relatively small numerator and denominator,
such as 355/113 for an approximation to PI. A large fraction approximating PI would be
something like 141360325243/44996389039 or the trivial case of 314159265359/10000000000,
this last one being "equal" to 3.14159265359, the HP48's representation of PI.
Unfortunately "small fraction" is not precisely defined, but I think my illustration shows
what Joe Horn and Bill Wickes meant by it.
>
>Can you give an example of one input from the large number of
>cases you have in mind, and approximately what proportion of
>say 100 random inputs would fall into this class?
Let me say some things about decimal numbers and fractions. All the decimal numbers
representable on the HP48 are of the form .xxxxxxxxxxxx and are representable exactly as
fractions in at least one way (the trivial case). For example, .1=1/10; .12=12/100;
.123=123/1000; .123456789=123456789/1000000000, etc.
About fractions.
Fractions are rational numbers. When fractions are converted to decimal numbers,
sometimes the decimal representation goes on indefinitely (a repeating decimal), such as
1/7=.142857142857142857... Other fractions terminate after a finite number of decimal
places, such as 1/8=.125; 1/10=.1; 777/800=.97125; etc. Because we use the base 10 number
system, the only fractions that have a terminating decimal representation are those
(reduced to lowest terms) whose denominator has only prime factors of 2 and/or 5. If we
consider all the fractions that can be represented on the HP48, they must have
denominators between 1 and 999999999999. In that range of integers, there are only 291
that have only prime factors of 2 and/or 5 (I hope I haven't made an error in the
computation). To get all of the fractions that can be had on the HP48, we need to let the
numerator and denominator both range from 1 to 999999999999, independently of one another,
and that's a LOT of fractions. Fractions with larger or smaller exponents can be reduced
to one of these. And some of these are not reduced to lowest terms. So, at most, 291 of
the possible denominators form fractions (reduced to lowest terms) that are terminating.
So, if you pick a numerator and denominator at random from the range 1 to 999999999999,
the fraction you have made is VERY unlikely to have a terminating decimal representation.
But decimal numbers that can be represented on the HP48 are all inherently terminating
(terminated?), since they can only be 12 digits at most. Therefore the decimal numbers
that can be represented on the HP48 are VERY unlikely to be equal (and remember, I mean by
this EXACTLY equal in numeric value) to any fraction except the trivial one.
So, we see that a fraction randomly chosen from the all possible fractions representable
on the HP48 is VERY unlikely to also be representable as a decimal fraction on the HP48.
Now, going from the other direction, since all decimal fractions on the HP48 are
terminating (at most 12 digits, remember), they all have a terminating continued fraction
representation. If the CF is used to form the sequence of convergents, we will see that
the LAST convergent is EXACTLY equal to the input decimal. And, NONE of the earlier
convergents is exactly equal to the decimal, or the CF sequence would have terminated
there. The only way we can get a small fraction exactly equal to the input decimal is if
the input decimal has only a few significant digits to start with, such as .125 = 1/8.
So, all this is what I meant when I said to Joe:
"Even you aren't able to find a small fraction whose value is
"equal" to the input argument except in a relatively small
number of cases."
To answer your question, the large number of cases I referred to is the set of all
possible decimal numbers representable on the HP48. Some members of the relatively small
number of cases I referred to are .125; .5; .222, etc. These have small fraction
representations with relative errors of zero; they are equal to the input decimal.
Maybe somewhat less than 291/999999999999 = .0000000291% of the cases will have a small
fraction representation that is exact. If this is all there are, then no algorithm can
find more.
Rodger Rosenbaum
> -=[ Sun, 19.10.03 10:38 a.m. +1300 (NZDT) ]=-
>
>Hi Tony Hutchins ! <jus...@csi.com>
>
>01h06m ago, on Sun, 19.10.03 at 09:32 a.m. +1300 (NZDT), you wrote
>in message ID <1062568...@news.cis.dfn.de> :
>
>> I am obviously not sure what you mean by '"equals"', or by
>> "small fraction".
>>
>> Can you give an example of one input from the large number of
>> cases you have in mind, and approximately what proportion of
>> say 100 random inputs would fall into this class?
>
>Oh! You must mean inputs like 1,2,3 etc which can only result
>in 1/1,2/1,3/1... or say 3.14159 which only results in
>314159/100000 - these are not on the "small fraction" class?
>
>And, by "equals" you mean that if you subtract the input from
>the output you get zero.
>
>So, is your point not purely semantic?
>
>One must allow Joe some poetic license ;-)
This is mathematics, not poetry; he gets only mathematical license.
But I'll try not to let ambition mock his useful toil.
>
>And, in the examples for Joe's challenge, all inputs and final
>outputs both look equal, and are also indeed equal by your
>stronger test.
This was obviously posted before another rather long posting I made on the topic. My
stronger test is to compute the value of the fractions in exact arithmetic, and then
compare them to the inputs, computing the relative error as I explained in the other
posting. You can't do this on the HP48, because it is just the limitations of the
floating point system that falsely make the fraction and input decimal appear to be equal
(equal as I defined it in the other posting; I mean EXACTLY equal).
Tony Hutchins
Hi Rodger Rosenbaum ! <rodg...@siteconnect.com>
11m ago, on Sat, 18.10.03 at 5:05 p.m. -0700, you wrote
in message ID <1uk3pv4bcfher54og...@4ax.com> :
[...]
> This was obviously posted before another rather long
> posting I made on the topic.
It was. Thanks for your detailed posting! My mindset was
focussed on display/representation on the HP48, not on the
exact arithmetic. It all makes sense now - "equal" certainly
needs qualification :)
--
Tony Hutchins
Wellington
New Zealand
#1 Speed Kills. Use Microsoft Windows.
Eddy Hondo
Whats the use for these best fractions????
I know that this has been asked and answered to a lot but i STILL DONT
GET IT.
If i use this to convert dec to fraction, it is at least for the
calculator exact and best, but in a math test for example, doing this
conversion would be wrong, because of the certain mispresicion in
play. IMO the best thing to do is just *1e12/1e12 and simplyfy, its
exact and it works.
MAybe the real use for this is internal for the calculator or somethin
like that.
Also, would these best fractions benefit in acurracy?? i mean,
assuming the way that the 48/49g series work.
Also, is there any guarantee that the best fractions found are truly
best fractions?
BTW, interesting thing you found Joe Horn.
Tony Hutchins
Hi Joe,
1 day 16h24m ago, on Fri,
17.10.03 at 03:35 a.m. -0700, you wrote
in message ID <87233f9e.03101...@posting.google.com> :
> ------- begin Wickes quotation -------
>
> ->Q finds the fraction with the *smallest denominator*
> that matches the argument within the error limit set by
> the display format. For example,
I now wonder what is meant by "the error limit set by
the display format". For example:
DO: 5 FIX 3 SQRT ->Q EVAL
SEE: 1.73205 '362/209' 1.73206
I was expecting the input/output displays to be the same.
In fact root 3 is 1.7320508..
and 326/209 is 1.7320574...
The difference is clearly more than .5*E-6, which is what I
was expecting.
Quite surprising?!
Probably PDQ, in emulating what ->Q is intended to do, would
produce 989/571 = 1.7320490...
Ironic that we would need a bigger denominator, in this case.
Possibly ->Q just produces "the smallest convergent that is
within +/- 'x' of the number_as_displayed", but they tried to
find a way to describe it without using the term "convergent".
I haven't a clue what 'x' is. Maybe if the next convergent is
suitable for "FIX 6", and the previous for "FIX 4", then it
uses the current one for "FIX 5", no matter how far away it
is?
By number_as_displayed I just refer to what we see, completely
ignoring how it might be represented *in* the HP48, or how we
might think of it *outside* the HP48, in it's full
arithmetical accuracy.
--
Tony Hutchins
Wellington
New Zealand
#312 ... for 'tis better to be alone
than in bad company. George Washington *
James M. Prange
<snip>
> The Hurwitz accuracy seems to be from an era when
> gears were used to divide quantities; it was good to have gears without too many teeth.
> Now, we do these things electronically and we don't need a small number of teeth on a
> gear. Now, we might prefer the lowest relative error.
Well, using gears to divide or multiply for calculating purposes is
probably obsolete, but gears are not obsolete by any stretch of the
imagination. As you noted, having too large a number of teeth on a gear
may not be feasible, and of course, too small a number won't work
either. And of course, the number of teeth has to be a whole number;
fractional teeth on a gear just don't work well at all. ;-) When
designing a gear set, the gear ratio is a very important consideration.
When an exact ratio is required, the math is rather trivial, although
some exact ratios just aren't feasible. But often the designer is
looking for an approximate ratio. It seems to me that Joe's method,
modified to allow the user to set minimum and maximum values for the
numerator and denominator, would be ideal for this purpose.
<snip>
--
Regards,
James
Wayne Brown
Thanks; I'm looking forward to it!
--
Wayne Brown (HPCC #1104) | "When your tail's in a crack, you improvise
fwb...@bellsouth.net | if you're good enough. Otherwise you give
| your pelt to the trapper."
"e^(i*pi) = -1" -- Euler | -- John Myers Myers, "Silverlock"
Werner Huysegoms
> > First, I used a list - but that doesn't work an a 48S.
>
> Why's that? Does a 48G have extra LIST features?
> Couldn't you use LIST-> just like you use C->R?
>
with a list, I would have used EVAL instead of LIST\-> DROP,
but on a 48S you cannot multiply a list with a number, it does
not have 'Automatic List Processing' as the 48G and 49 do.
(nor the command ADD to add up two lists of numbers)
Cheers,
Werner