# PDQ Unleashed: No More Limitations

51 views

### Joe Horn

Oct 7, 2004, 5:58:09 AM10/7/04
to
PDQ2.TXT: documentation for the PDQ2 Directory & Program

The PDQ Algorithm, Version 2, for the HP49G and hp49g+.
Unlimited Precision Version of the PDQ Fraction Finder Algorithm.
Dedicated to Rodger Rosenbaum, who said that it couldn't be done.
Code and documentation: Copyright (C) 2004 by Joseph K. Horn.
Not to be used in ANY commercial product without written permission.

"Somebody said that it couldn't be done,
But he with a chuckle replied
That MAYBE it couldn't, but he would be one
Who wouldn't say so till he tried.

So he buckled right in with a trace of a grin
On his face. If he worried, he hid it.
He started to sing as he tackled the thing
That couldn't be done, and he did it."
-- Edgar Guest

INTRODUCTION: Three Problems of Increasing Difficulty

(a) What is the simplest fraction (smallest possible numerator and
denominator) between 0.4 and 0.6? The answer is obviously 1/2. No
software is needed for such a simple problem. Even though there are an
infinite number of fractions between 0.4 and 0.6, it's obvious that
they all have numerators and denominators greater than 1/2 does, so 1/2
is clearly the "winner" of the "simplest fraction award" in that
interval.

(b) What is the simplest fraction between 0.671 and 0.672? The answer
is not quite as obvious this time! Answer: 43/64 (=0.671875). PDQ2
finds the answer, as do PDQ1, DEC2FRAC, the HP-32SII, the hp33s, the
Sparcom Math-Pro cards for the HP48, and many other tools.

(c) What is the simplest fraction equal to pi +/- 261 parts per
non-obvious, it's impossible to calculate in one's head. In fact, it is
not attainable by *anything* other than PDQ2! [Challenge: A prize of
pi^3 dollars +/- 1 part per hundred will be awarded to the first person
who demonstrates to me any software or hardware that gets the right
899125804609/286200632530.

*ALL* questions like these, regardless of complexity, can be answered
by PDQ2 swiftly and exactly. The inputs and outputs can be of any
desired length, and are always exact (no round-off errors ever occur).
Decimals can be entered either as a string (e.g. "3.14159265358979") or
as a ratio of integers (e.g. '314159265358979/10^14'), thus allowing
for inputs of any desired length. The outputs are always ratios of
integers.

CONTENTS of the PDQ2 directory:

(1) PDQ2 - finds simplest fraction within any interval.
(2) CF~I2 - toggles between Continued Fraction & two integers.
(3) ISQR - infinite-precision square root of integer.
(4) PI64 - Pi to 64 decimal places.
(5) TOF - decimal-to-fraction by definition of decimals.
(6) \$TOF - same as TOF but uses strings for unlimited precision.
(7) FRAC? - tests if object is a ratio of two integers.

INSTRUCTIONS
___________________________________________________________

(1) PDQ2, "Simplest Ratio of Integers = Target +/- Tolerance"

Input:
Level 2: Exact TARGET (real, string, or ratio of two integers)
Level 1: Exact TOLERANCE (real, or ratio of two integers)

If the input tolerance is a ratio of two integers, or a decimal between
0 and 1, it is used as-is.

If the input tolerance T is >1, it is taken to mean that many TRUNCATED
decimal places of accuracy are desired, and is automatically replaced
by 1/10^T. For exact results, use only integers.

If the input tolerance T is <0, it is taken to mean that many ROUNDED
decimal places of accuracy are desired, and is automatically replaced
by 5/10^(|T|+1).

Output:
Level 2: Simplest Fraction (ratio of two integers)
Level 1: N/X:(sign)Exact Error (ratio of two integers)

The "N" or "X" tag on the output indicates whether the output is a
"normal" or "extra" result. "Normal" results are convergents of the
continued fraction expansion of the input target. "Extra" results are
"intermediate convergents" that cannot be found by most other software.

The sign of the error indicates whether the output is greater than
(positive) or less than (negative) the target. For example, an error
of -1/10 would mean that the output is exactly 0.1 less than the
target.

Pressing [-] [EVAL] converts the output back into the original target.

Pressing ABS LOG ABS on the level-1 output changes it from the actual
error into the number of decimal places of accuracy.

If the specified tolerance exceeds the precision of the target, then
its invisible trailing zeros will be assumed to be significant digits.
For example, 3. 1/x 12 PDQ2 --> 1/3 as expected, but 3. 1/x 13 PDQ2 -->
256410256410/769230769231, because that's the simplest fraction that
rounds to .3333333333330 (notice the trailing zero in the 13th decimal
place!). The PDQ2 Algorithm follows the HP Calculator Philosophy of
assuming that the user knows what he's doing and that his inputs are
exactly what he intended to input. Whenever PDQ2 gets an answer that
you don't expect, raise your expectations accordingly.

Examples (the three from the Introduction above, plus one more):

(a) What is the simplest fraction (smallest possible numerator and
denominator) between 0.4 and 0.6, inclusive? In other words, what's
the simplest fraction equal to .5 +/- a tolerance of .1?
.5 .1 PDQ2 --> '1/2' N:0
This means that the answer is 1/2, with no error, and it's a "normal"
result.

(b) What is the simplest fraction between 0.671 and 0.672, inclusive?
In other words, what's the simplest fraction that's equal to 0.6715 +/-
0.0005?
.6715 .0005 PDQ2 --> '43/64' X:3/8000
This means that the answer is 43/64, which is 3/8000 above the target,
and is an answer that cannot be discovered by the simple continued
fraction algorithm that everybody else uses.

(c) What is the simplest fraction equal to pi with a tolerance of +/-
261 parts per septillion (US) or per quadrillion (UK)?
[PI64] 261 24 10^x / PDQ2 --> 899125804609/286200632530 X:big/bigger
That fraction (in level 2) is greater than the target by exactly the
amount in level 1... a fraction so long that you'll probably need to
press ->NUM to get an idea of its magnitude.

(d) What's the simplest fraction for sqrt(13) on an 80-bit calculator
(e.g. the hp30s)?
If the calculator truncates its results, use this:
13 30 ISQR 2 -80 y^x PDQ2 --> 6536561700241/1812916028881, slightly
less than sqrt(13).
If the calculator rounds its results, use this instead:
13 30 ISQR 2 -81 y^x PDQ2 --> 7955840589842/2206553168161, slightly
less than sqrt(13).

___________________________________________________________

(2) CF~I2, "Toggle Between Continued Fraction and Two-Integer Ratio"

Input --> Output options:

(a) { a0 a1 a2 a3 ... } --> P Q
(b) 'P/Q' --> { a0 a1 a2 a3 ... }
(c) P Q --> { a0 a1 a2 a3 ... }
(d) Decimal number --> { a0 a1 a2 a3 ... }

Converts numbers to and from "continued fraction" form, represented as
a list containing the integer part of the input (a0) followed by the
"partial quotients" (a1, a2, a3, etc). P and Q are the numerator and
denominator, respectively, of the equivalent regular fracton. All
outputs are exactly equal to the input; no roundoff errors occur.

Examples:

.123 CF~I2 --> { 0 8 7 1 2 5 }
This means that .123 = 0+1/(8+1/(7+1/(1+1/(2+1/5))))

{ 0 8 7 1 2 5 } CF~I2 --> 123 1000
This means that 0+1/(8+1/(7+1/(1+1/(2+1/5)))) = 123/1000

'15/37' CF~I2 --> { 0 2 2 7 }
This is the same as:
15 37 CF~I2 --> { 0 2 2 7 }

Some irrational numbers have continued fractions containing discernable
patterns. For example, the constant 'e' (approx. 2.71828182846) is
exactly equal to the continued fraction { 2 1 2 1 1 4 1 1 6 1 1 8 1 1
10 1 1 12 1 1 14 ... }. This fact allows CF~I2 to be used to generate
fractions that are better than the 12-digit limit of the hp49g+.

Example:
{ 2 1 2 1 1 4 1 1 6 1 1 8 1 1 10 1 1 12 1 1 14 2 } CF~I2
--> 848456353 312129649
This ratio (848456353/312129649) is the same as 'e' to 18 decimal
places. Use more partial quotients to get as close to 'e' as you want.

___________________________________________________________

(3) ISQR, "Infinite-Precision Square Root of Integer"

Input:
Level 2: The integer to be square-rooted.
Level 1: The number of desired ROUNDED decimal places.

Output: 'p/q', exactly equal to the decimal expansion of the square
root rounded to the specified number of digits.

Example:

What is the exact sqrt(50) rounded to 25 decimal places?
50 25 ISQR --> '17677669529663688110021109/2500000000000000000000000'
This fraction can then be used as an input for PDQ2, keeping in mind
that the tolerance should not exceed 24 decimal places. If more places
are desired, use a larger input for ISQR.

The ISQR program is included in the PDQ2 directory only to help you
generate super-accurate square roots for the PDQ2 program to chew on.

___________________________________________________________

(4) PI64, "Pi to 64 decimal places"

To experiment with pi as an input to PDQ2, you can use PI64 to put pi
on the stack in string form, accurate to 64 decimal places. Therefore,
when using PI64 with PDQ2, do not specify a tolerance greater than 64,
or trailing zeros will be assumed as significant digits, diverging away
from pi. If more than 64 digits of pi are needed (for what?!?), you'll
have to type them in.

Example:

PI64 49 PDQ2--> 17579871157983393066594923/5595846787423398162594219
This is the simplest fraction in the interval [pi-1/10^49, pi+1/10^49].
Furthermore, only PDQ2 can find this answer, and it finds it very
quickly.

___________________________________________________________

(5) TOF, "To Fraction" using the simple definition of decimals.

Input: any real number or integer or ratio of integers.
Output: P Q (where P/Q is exactly equal to the input).

TOF is only included in the PDQ2 directory because it is called by the
PDQ2 program. It is rather useless by itself, since its results are
rather brainless. It merely puts the input over 1 and then multiplies
by 10/10 until the numerator is an integer. If a fraction is input, it
merely performs a FXND.

Examples:
'2/3' TOF --> 2 3
2. 3. / TOF --> 666666666667 1000000000000
pi ->NUM TOF --> 314159265359 100000000000

___________________________________________________________

(6) \$TOF, "To Fraction" for long inputs in string form.

Input: any decimal number as a string.
Output: 'P/Q' (where P/Q is exactly equal to the input).

See TOF for details.

Examples:
".123456789012345" \$TOF --> '24691357802469/200000000000000'
which is *exactly* equal to .123456789012345, even though HP
calculators cannot directly verify that fact, since they are limited to
12-digit precision. The \$TOF program (and the PDQ2 program) do not
have that limitation.

___________________________________________________________

(7) FRAC?, "Is this a fraction?"

Input: any object
Output: same object; 0. or 1. (1. if ratio of two integers, 0.
otherwise)

Examples:
3.1416 FRAC? --> 3.1416 0. (because the object 3.1416 is not a ratio of
two integers)
'3927/1250' FRAC? --> '3927/1250' 1. (because that object IS a ratio of
two integers)
'1.2/3.4' FRAC? --> '1.2/3.4' 0. (it's NOT a ratio of two integers)
"Foo" FRAC? --> "Foo" 0.
153 FRAC? --> 153 0.

Included here only because it gets called by the TOF program.

___________________________________________________________

SOURCE CODE

The programs in this 'PDQ2' directory are written entirely in
un-optimized User RPL for the sake of those who wish to study how the
programs work. Optimizing the code, or rewriting it in System RPL, or
in Saturn Assembly Language, or in ARM Assembly Language, is left as an
exercise for the student.

Google apparently removes the indenting that HP so carefully puts into
program listings. Sorry 'bout that.

%%HP: T(3)F(.);
@ PDQ2 directory by Joe Horn (C) Copyright 2004
@ See PDQ2.TXT for details.
DIR
PDQ2
@ PDQ Algorithm demo, version 2
@ Input: Exact target (%|\$|'n/d'), Exact tolerance (%|'a/b')
@ Output: Exact rational approximation ('p/q'), Exact error ('x/y')
@ Example: "3.14159265358979", 14
@ --> '53047102/16885417', '-15109243/1688541700000000000000'
\<<
IF OVER TYPE 2. SAME
THEN SWAP \$TOF SWAP
END
IF OVER FXND MOD
THEN
IF DUP 0 <
THEN NEG 2. LOG +
END
IF DUP 1 \>=
THEN NEG ALOG
END TOF ROT TOF DUP2 1 1 0 0 0 1 0 1 1 0 0
\-> a b n d n0 d0 cn x pn cd y pd lo hi mid q r
\<<
DO n d IDIV2 'r' STO 'q' STO
q cn * pn + 'x' STO
q cd * pd + 'y' STO
cn 'pn' STO
x 'cn' STO
cd 'pd' STO
y 'cd' STO
d 'n' STO
r 'd' STO
UNTIL b n0 y * d0 x * - ABS * a d0 * y * \<=
END
IF q 1 >
THEN q 'hi' STO
DO
lo hi + 2 IQUOT 'mid' STO
cn pn mid * - 'x' STO
cd pd mid * - 'y' STO
IF b n0 y * d0 x * - ABS * a d0 * y * \<=
THEN mid 'lo' STO
ELSE mid 'hi' STO
END
UNTIL hi lo - 1 \<=
END
cn pn lo * - 'x' STO
cd pd lo * - 'y' STO
END x y / DUP n0 d0 / - EVAL
IF lo
THEN "X"
ELSE "N"
END \->TAG
\>>
ELSE DROP :N: 0
END
\>>

CF~I2
@ Toggles between Continued Fraction (list) and two integers (p/q).
@ {} --> p q; e.g. { 3 7 16 } --> 355 113
@ p q --> {}; e.g. 355 113 --> { 3 7 16 }
@ %.% --> {}; e.g. 1.2345 --> { 1 4 3 1 3 1 1 2 5 }
@ 'p/q' --> {}; e.g. '355/113' --> { 3 7 16 }
\<<
IF DUP TYPE 5. SAME
THEN LIST\-> R\->B 1 0 # 1h 4. ROLL
START OVER 4. ROLL * + SWAP
NEXT
ELSE
IF DUP FXND MOD
THEN FXND
END DUP2
IF 12 ALOG < SWAP 12 ALOG < AND
THEN I\->R SWAP I\->R SWAP
END DEPTH DUP 1. + ROLLD SWAP
WHILE OVER DUP2 MOD ROT OVER - ROT / SWAP DUP
REPEAT ROT
END ROT DROP2 DEPTH ROLL DEPTH SWAP - 1. + \->LIST R\->I
END
\>>

ISQR
@ "Infinite-Precision Square Root"
@ Input: x (integer), p (integer)
@ Output: 'p/q' exactly equal to sqrt(x) rounded to p decimal places.
@ Example: 3 14 --> '21650635094611/12500000000000'
\<< SWAP OVER 1 + 2 * ALOG * DUPDUP \->STR
1. OVER SIZE 2. / SUB OBJ\->
WHILE DUP2 IQUOT + 2 IQUOT ROT OVER - SIZE
REPEAT SWAP OVER
END NIP DUP 10 IREMAINDER + 10 IQUOT SWAP ALOG /
\>>

PI64
@ Pi to 64 decimal places, just for playing with.
"3.1415926535897932384626433832795028841971693993751058209749445923"

TOF
@ "To Fraction" using the definition of decimal numbers.
@ e.g. 1.2 --> 12/10 & .333333333333 --> 333333333333/1000000000000
@ Note: the output is always exactly equal to the input.
\<<
IF FRAC?
THEN FXND
ELSE
IF DUP FXND MOD
THEN 11. OVER XPON - ALOG SWAP OVER * R\->I SWAP R\->I SIMP2
ELSE 1
END
END
\>>

\$TOF
@ Same as 'TOF' but inputs strings (of unlimited length).
@ e.g. "1.2345678901234567" --> 1234568901234567/10000000000000000'
@ Note: the output is always exactly equal to the input.
\<< DUP "." POS DUP2 1. - 1. SWAP SUB PICK3 PICK3 1 +
OVER SIZE SUB + STR\-> ROT SIZE ROT - R\->I ALOG /
\>>

FRAC?
@ Input: obj
@ Output: obj, 0/1 (1 if obj is 'integer/integer', else 0)
\<<
IF DUP TYPE 9. -
THEN 0.
ELSE DUP \->LST DUP 1.
\<< TYPE
\>> DOSUBS R\->I { 28 28 18 } SAME SWAP DUP SIZE GET
{ / } 1. GET SAME AND
END
\>>
END

### Veky (hp49G)

Oct 7, 2004, 11:38:43 AM10/7/04
to
//* Challenge: A prize of

pi^3 dollars +/- 1 part per hundred will be awarded to the first person
who demonstrates to me any software or hardware that gets the right
answer to this and all similar problems. -jkh- *//

WRI Mathematica 5.0 :
Rationalize[Pi,261/10^(4*6)]
returns exactly what's needed.

I want my 31bucks. :-]

### Bhuvanesh

Oct 7, 2004, 2:44:13 PM10/7/04
to
"Joe Horn" <joe...@holyjoe.net> wrote:
> (c) What is the simplest fraction equal to pi +/- 261 parts per
> septillion (US) or per quadrillion (UK)? The answer is not only
> non-obvious, it's impossible to calculate in one's head. In fact, it is
> not attainable by *anything* other than PDQ2! [Challenge: A prize of
> pi^3 dollars +/- 1 part per hundred will be awarded to the first person
> who demonstrates to me any software or hardware that gets the right
> 899125804609/286200632530.

This one is easy :-) This has been a part of Mathematica since at least version 3.0.

\$Version = 5.0.1 for Microsoft Windows (November 18, 2003)
\$TopDirectory = C:\Program Files\Wolfram Research\Mathematica\5.0.1

In[1]:= Rationalize[Pi, 261/10^24]

899125804609
Out[1]= ------------
286200632530

Cheers,
Bhuvanesh.

### Rodger Rosenbaum

Oct 7, 2004, 6:44:24 PM10/7/04
to
On 7 Oct 2004 02:58:09 -0700, "Joe Horn" <joe...@holyjoe.net> wrote:

>PDQ2.TXT: documentation for the PDQ2 Directory & Program
>
>The PDQ Algorithm, Version 2, for the HP49G and hp49g+.
>Unlimited Precision Version of the PDQ Fraction Finder Algorithm.
>Dedicated to Rodger Rosenbaum, who said that it couldn't be done.

*Exactly* what was it I said that you believe you have provided a counterexample for
in PDQ2?

I say *exactly* because I want to avoid a long back and forth like we had in the

### Tony Hutchins

Oct 7, 2004, 10:53:05 PM10/7/04
to
-=[ Fri, 8.10.04 3:12 p.m. +1300 (NZDT) ]=-

Hi Rodger Rosenbaum ! <nos...@aol.com>

03h28m ago, on Thu, 7.10.04 at 3:44 p.m. -0700, you wrote
in message ID <uhhbm050o0dddbdi3...@4ax.com> :

> On 7 Oct 2004 02:58:09 -0700, "Joe Horn" <joe...@holyjoe.net> wrote:
>
> >PDQ2.TXT: documentation for the PDQ2 Directory & Program
> >
> >The PDQ Algorithm, Version 2, for the HP49G and hp49g+.
> >Unlimited Precision Version of the PDQ Fraction Finder
> >Algorithm. Dedicated to Rodger Rosenbaum, who said that
> >it couldn't be done.
>
> *Exactly* what was it I said that you believe you have
> provided a counterexample for>

> I say *exactly* because I want to avoid a long back and

I remember this sentence well:-

> Subject: Re: 48/49: Best Fraction Challenge
> Date: Sat, 18 Oct 2003 06:16:23 -0700
> 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.

You meant "equal" in the strictest mathematical sense. Joe's
PDQ1 was limited to calculator precision. His PDQ2 handles
unlimited precision, which approaches your "equal"
assymptotically at least. I don't know whether he has
expanded on your "relatively small number of cases" though.

Perhaps Joe took your statement above as a challenge.

--
Tony Hutchins
New Zealand

### The Phantom

Oct 7, 2004, 11:23:37 PM10/7/04
to
On 7 Oct 2004 02:58:09 -0700, "Joe Horn" <joe...@holyjoe.net> wrote:

Not if one is an Idiot Savant; in fact an Idiot Savant friend of mine wants to know
what is the "simplest fraction" equal to Pi +/- 13131 parts in 10^440 and he wants to know
how long it takes PDQ2 to find it. It took him less than a second! Way less!

In fact, it is
>not attainable by *anything* other than PDQ2! [Challenge: A prize of
>pi^3 dollars +/- 1 part per hundred will be awarded to the first person
>who demonstrates to me any software or hardware that gets the right
>899125804609/286200632530.
>
>*ALL* questions like these, regardless of complexity, can be answered
>by PDQ2 swiftly and exactly. The inputs and outputs can be of any
>desired length,

I don't think you really meant to say that. Can they be 10^20 digits long? Looks like
imprecise speaking to me.

and are always exact (no round-off errors ever occur).
>Decimals can be entered either as a string (e.g. "3.14159265358979") or
>as a ratio of integers (e.g. '314159265358979/10^14'), thus allowing
>for inputs of any desired length. The outputs are always ratios of
>integers.
>
>

Snip...

>(4) PI64, "Pi to 64 decimal places"
>
>To experiment with pi as an input to PDQ2, you can use PI64 to put pi
>on the stack in string form, accurate to 64 decimal places. Therefore,
>when using PI64 with PDQ2, do not specify a tolerance greater than 64,
>or trailing zeros will be assumed as significant digits, diverging away
>from pi.

You don't really mean to say that adding trailing zeroes to a 64 decimal approximation
to Pi will cause its value to diverge (further than it already is) from Pi? Won't the
value of an N-digit approximation to Pi have exactly the same value if any number of
trailing zeroes are added, rather than "diverging away from Pi" as zeroes are added?

If more than 64 digits of pi are needed (for what?!?), you'll
>have to type them in.
>

Snip...

>
>(6) \$TOF, "To Fraction" for long inputs in string form.
>
>Input: any decimal number as a string.
>Output: 'P/Q' (where P/Q is exactly equal to the input).
>
>See TOF for details.
>
>Examples:
>".123456789012345" \$TOF --> '24691357802469/200000000000000'
>which is *exactly* equal to .123456789012345, even though HP
>calculators cannot directly verify that fact, since they are limited to
>12-digit precision.

I think you mean to say something else here, rather than an unqualified "...HP
calculators cannot directly verify that fact...". It seems to me that the HP49 can, and
it's an HP calculator, isn't it?

### Rodger Rosenbaum

Oct 8, 2004, 12:15:38 AM10/8/04
to

I think it's limited by the amount of memory which the HP49 has, which means that his
algorithm, like everything else on the HP49, is limited to "calculator precision". It's
just that "calculator precision" is much bigger than it was on the HP48. It's the same
result you would get if you had a machine with arithmetic registers of 10,000 digits, or
50,000 digits, or 1,000,000 digits. That's what all arithmetic systems of so-called
"infinite precision" (metaphorically speaking, that is) do. The only way you can get true
infinite precision arithmetic is if all your numbers can be represented as rationals.
Speaking of which, I'm going to post a mini-challenge for the HP49.

>assymptotically at least. I don't know whether he has
>expanded on your "relatively small number of cases" though.

My original argument still applies. Whatever the largest fractions you can represent
in the HP49's memory, only a relatively small number of them have denominators with only 2
and 5 as prime factors.

### Tony Hutchins

Oct 8, 2004, 2:41:28 AM10/8/04
to
-=[ Fri, 8.10.04 7:26 p.m. +1300 (NZDT) ]=-

Hi The Phantom ! <pha...@aol.com>

03h03m ago, on Thu, 7.10.04 at 8:23 p.m. -0700, you wrote
in message ID <p02hm0915vjoepg8l...@4ax.com> :

> >(c) What is the simplest fraction equal to pi +/- 261 parts per
> >septillion (US) or per quadrillion (UK)? The answer is not only
> >non-obvious, it's impossible to calculate in one's head.
>
> Not if one is an Idiot Savant;

;-)

> >*ALL* questions like these, regardless of complexity, can
> >be answered by PDQ2 swiftly and exactly. The inputs and
> >outputs can be of any desired length,
>
> I don't think you really meant to say that. Can they be
> 10^20 digits long? Looks like imprecise speaking to me.

OK desired length is limited by memory.

> >To experiment with pi as an input to PDQ2, you can use
> >PI64 to put pi on the stack in string form, accurate to
> >64 decimal places. Therefore, when using PI64 with PDQ2,
> >do not specify a tolerance greater than 64, or trailing
> >zeros will be assumed as significant digits, diverging
> >away from pi.
>
> You don't really mean to say that adding trailing zeroes to
> a 64 decimal approximation to Pi will cause its value to
> diverge (further than it already is) from Pi? Won't the
> value of an N-digit approximation to Pi have exactly the
> same value if any number of trailing zeroes are added,
> rather than "diverging away from Pi" as zeroes are added?

It diverges in the sense that it becomes more and more laden
with spurious significant digits. It therefore diverges in
spuriousness, not in value.

> >Examples:
> >".123456789012345" \$TOF --> '24691357802469/200000000000000'
> >which is *exactly* equal to .123456789012345, even though HP
> >calculators cannot directly verify that fact, since they are limited to
> >12-digit precision.
>
> I think you mean to say something else here, rather than
> an unqualified "...HP calculators cannot directly verify
> that fact...". It seems to me that the HP49 can, and it's
> an HP calculator, isn't it?

Yes but
24691357802469 ENTER 200000000000000 / ->
.123456789012 not .123456789012345

And, the "directly" qualifies the statement just fine.

--
Tony Hutchins

### Tony Hutchins

Oct 8, 2004, 3:14:35 AM10/8/04
to
-=[ Fri, 8.10.04 7:56 p.m. +1300 (NZDT) ]=-

Hi Rodger Rosenbaum ! <nos...@aol.com>

02h40m ago, on Thu, 7.10.04 at 9:15 p.m. -0700, you wrote
in message ID <n44cm0tm86mq7ptqj...@4ax.com> :
[...]

> >You meant "equal" in the strictest mathematical sense. Joe's
> >PDQ1 was limited to calculator precision. His PDQ2 handles
> >unlimited precision,
>
> I think it's limited by the amount of memory which the
> HP49 has

Oops, yes, I forgot about memory.

--
Tony Hutchins

### Rodger Rosenbaum

Oct 8, 2004, 5:00:22 AM10/8/04
to

I'm trying to get people to use proper mathematical terminology, since mathematics is
what we're talking about here. "Diverges" has a standard meaning, and if it isn't
qualified then it doesn't mean "diverges in spuriousness"; I've never heard of such a
thing. How would anybody have known that he didn't mean "diverges in value"? That's what
it means if you don't say otherwise. If he wants to convey to the reader that the extra
zeroes COULD be treated as (spurious) significant digits, he should say just that.

Furthermore, extra zeroes are not necessarily spurious. Remember that recent MC, "How
many got polled?" The first percentage was 47.61% and then he asked about 47.610%; the
zero in the second number was not spurious. Joe says above, "...trailing zeros will be
assumed as significant digits..."; if that's what the user wants, then those zeroes are
not spurious. Joe was not clear.

>
>> >Examples:
>> >".123456789012345" \$TOF --> '24691357802469/200000000000000'
>> >which is *exactly* equal to .123456789012345, even though HP
>> >calculators cannot directly verify that fact, since they are limited to
>> >12-digit precision.
>>
>> I think you mean to say something else here, rather than
>> an unqualified "...HP calculators cannot directly verify
>> that fact...". It seems to me that the HP49 can, and it's
>> an HP calculator, isn't it?
>
>Yes but
>24691357802469 ENTER 200000000000000 / ->
>.123456789012 not .123456789012345

This is in approximate mode, right?

>
>And, the "directly" qualifies the statement just fine.

I disagree that it qualifies the statement "just fine". What did HE mean by the word
"directly"? Too bad he didn't tell us. If the presence or absence of the word "directly"
is that important in this sentence, then he should tell us what it means. Apparently,
with regard to the HP49, you are suggesting that it means "in approximate mode"; this is
the first time I've heard that.

And it isn't true for this example to say that "HP calculators cannot directly verify that
fact" if that's what is meant by "directly". Any Saturn based HP calculator can use my
DIGI24 package from Goodies Disk 9 (or its algorithms) to verify this example and any
others up to 24 digits, and is therefore not limited to 12-digit precision as Joe asserts
above. And in fact, Dekker's method can be extended to any multiple of the calculator's
basic floating point precision. So if "directly" means "in approximate mode", then that's
what the pre-HP49 Saturn based calcs do, and the DIGI24 package can "directly" verify the
example above. If "directly" means something else, how would I know if Joe doesn't tell
me?

Maybe he meant that if you type 24691357802469 200000000000000 and press the divide key
(or execute some other single built-in function) on the HP49 in approximate mode (or a
pre-49 Saturn based machine), then you won't get the 15 digit result. Of course, you
can't even put the 14-digit numerator on the stack of an HP49 in approximate mode (or a
pre-49 Saturn based machine). With this I would have to agree, and he should say that
this is what he means if he does. (You could use extended precision system RPL functions
I suppose, and verify with a single extended precision divide.)

But I would point out that if (in exact mode) you type the numerator on an HP49 with 15
more zeroes after it, hit enter and then type the denominator followed by IDIV2, you get
123456789012345--that's pretty "direct" verification, so it could be said that the HP49
can "directly" verify this example. This is possible because HP49 calculators are NOT
(unqualifiedly) "...limited to 12-digit precision." If they were, his PDQ2 program
wouldn't work the way it does. Had Joe said "...HP calculators...are limited to 12-digit
precision (including the HP49 in approximate mode) in their four basic arithmetic
functions.", that would have been a true, qualified, statement with which I could not
disagree.

### Tony Hutchins

Oct 8, 2004, 6:36:03 AM10/8/04
to
-=[ Fri, 8.10.04 10:31 p.m. +1300 (NZDT) ]=-

Hi Rodger Rosenbaum ! <nos...@aol.com>

31m ago, on Fri, 8.10.04 at 02:00 a.m. -0700, you wrote
in message ID <volcm0ltfn6gnvr9c...@4ax.com> :

> >> >To experiment with pi as an input to PDQ2, you can use
> >> >PI64 to put pi on the stack in string form, accurate to
> >> >64 decimal places. Therefore, when using PI64 with PDQ2,
> >> >do not specify a tolerance greater than 64, or trailing
> >> >zeros will be assumed as significant digits, diverging
> >> >away from pi.

Let me interject and ask for the last sentence to be re-read.
It is quite clear what is meant if the last phrase, which
merely qualifies "trailing zeros", is removed.

> >> You don't really mean to say that adding trailing zeroes to
> >> a 64 decimal approximation to Pi will cause its value to
> >> diverge (further than it already is) from Pi? Won't the
> >> value of an N-digit approximation to Pi have exactly the
> >> same value if any number of trailing zeroes are added,
> >> rather than "diverging away from Pi" as zeroes are added?
> >
> >It diverges in the sense that it becomes more and more laden
> >with spurious significant digits. It therefore diverges in
> >spuriousness, not in value.
>
> I'm trying to get people to use proper mathematical
> terminology, since mathematics is what we're talking about
> here. "Diverges" has a standard meaning, and if it isn't
> qualified then it doesn't mean "diverges in spuriousness";
> I've never heard of such a thing.

I made it up. How can those "trailing zeros" diverge in any
other way?

> How would anybody have known that he didn't mean "diverges
> in value"?

By reading the sentence which is quite specific.

> That's what it means if you don't say otherwise. If he
> wants to convey to the reader that the extra zeroes COULD
> be treated as (spurious) significant digits, he should say
> just that.

But that is exactly what the sentence addresses. The last
phrase is almost redundant.

> Furthermore, extra zeroes are not necessarily spurious.

Sure, but Joe was not adressing the whole universe of
possibilities.

> Remember that recent MC, "How many got polled?" The first
> the zero in the second number was not spurious. Joe says
> above, "...trailing zeros will be assumed as significant
> digits..."; if that's what the user wants, then those zeroes
> are not spurious. Joe was not clear.

Consider it without the last phrase:

> >> >Therefore, when using PI64 with PDQ2,
> >> >do not specify a tolerance greater than 64, or trailing

> >> >zeros will be assumed as significant digits.

Seems very clear to me. Maybe the "diverging away from pi"
could be left off, unless we can allow ourselves to endow
"diverge" with a merely English meaning. We are not even
talking about a series so "diverge" doesn't have a
mathematical meaning anyway.

> >> >Examples:
> >> >".123456789012345" \$TOF --> '24691357802469/200000000000000'
> >> >which is *exactly* equal to .123456789012345, even though HP
> >> >calculators cannot directly verify that fact, since they are limited to
> >> >12-digit precision.
> >>
> >> I think you mean to say something else here, rather than
> >> an unqualified "...HP calculators cannot directly verify
> >> that fact...". It seems to me that the HP49 can, and it's
> >> an HP calculator, isn't it?
> >
> >Yes but
> >24691357802469 ENTER 200000000000000 / ->
> >.123456789012 not .123456789012345
>
> This is in approximate mode, right?

Exact actually - I pressed ->NUM as well.
It seems a direct way - i.e. copy the example given directly??

> >> > '24691357802469/200000000000000' which is *exactly*
> >> > equal to .123456789012345, even though HP calculators
> >> > cannot directly verify that fact

> >

> >And, the "directly" qualifies the statement just fine.
>
> I disagree that it qualifies the statement "just fine".

I consider it meaningful in the context in which it is used.

> What did HE mean by the word

The connotation is surely just "execute the expression as given
by pressing a few keys and get the answer as quickly as
possible" - in particular use the divide key as shown in the
example.

[some snippage - ]

> But I would point out that if (in exact mode) you type
> the numerator on an HP49 with 15 more zeroes after it, hit
> enter and then type the denominator followed by IDIV2, you
> get 123456789012345--that's pretty "direct" verification,
> so it could be said that the HP49 can "directly" verify
> this example.

Pretty direct indeed ;-) Exactly what PDQ2 does.

> This is possible because HP49 calculators
> are NOT (unqualifiedly) "...limited to 12-digit precision."
> If they were, his PDQ2 program wouldn't work the way it does.
> Had Joe said "...HP calculators...are limited to 12-digit
> precision (including the HP49 in approximate mode) in their
> four basic arithmetic functions.", that would have been a
> true, qualified, statement with which I could not disagree.

Thanks Roger! Nor would I. Some HP calculators also have less
precision.

--
Tony Hutchins

### Daveh

Oct 8, 2004, 1:05:30 PM10/8/04
to
>>or trailing zeros will be assumed as significant digits, diverging away
>>from pi.

the only possible 'divergence', would be caused by the last digit due to
rounding, truncation etc. Though this could cause more errors in later
operations..

>We are not even talking about a series so "diverge" doesn't have a
mathematical meaning anyway.

An integral isn't a series, well maybe only the indefinite in the limit.
Neither is a function necessarily.

>>>'''''''''parts per septillion (US) or per quadrillion (UK)?

Daveh

"Tony Hutchins" <jus...@csi.com> wrote in message
news:1154688...@news.individual.net...

### Veli-Pekka Nousiainen

Oct 8, 2004, 1:55:34 PM10/8/04
to
"Rodger Rosenbaum" <nos...@aol.com> wrote in message
news:volcm0ltfn6gnvr9c...@4ax.com...
X

> And it isn't true for this example to say that "HP calculators cannot
> directly verify that
> fact" if that's what is meant by "directly". Any Saturn based HP
> calculator can use my
> DIGI24 package from Goodies Disk 9 (or its algorithms) to verify this
> example and any
> others up to 24 digits, and is therefore not limited to 12-digit precision
> as Joe asserts
> above. And in fact, Dekker's method can be extended to any multiple of
> the calculator's
> basic floating point precision. So if "directly" means "in approximate
> mode", then that's
> what the pre-HP49 Saturn based calcs do, and the DIGI24 package can
> "directly" verify the
> example above. If "directly" means something else, how would I know if
> Joe doesn't tell
> me?
X
I have use vectors with 6 significant numbers on one element
(at the most) for extended precision multiplication and addition
The routines were an adaption from the routines from
Jean-Michel Frerrard's book "Mastering Your HP-48 Volume 2:
Programming and Applications (Large Integers, pages 205)
[VPN]

### Tony Hutchins

Oct 8, 2004, 2:21:08 PM10/8/04
to
-=[ Sat, 9.10.04 06:45 a.m. +1300 (NZDT) ]=-

Hi Daveh ! <da...@worldonline.dk>

40m ago, on Fri, 8.10.04 at 7:05 p.m. +0200, you wrote
in message ID <GMz9d.56596\$Vf.26...@news000.worldonline.dk> :

> >>or trailing zeros will be assumed as significant digits, diverging away
> >>from pi.
>
> the only possible 'divergence', would be caused by the last digit due to
> rounding, truncation etc. Though this could cause more errors in later
> operations..
>
> >We are not even talking about a series so "diverge" doesn't have a
> mathematical meaning anyway.
>
> An integral isn't a series, well maybe only the indefinite in the limit.
> Neither is a function necessarily.

Point taken. I have almost turned into a troll<G>. The series
is the number of 'trailing zeros assumed as significant
digits' in: [PI64]0, [PI64]00, [PI64]000 etc., i.e. 1,2,3
etc. But, I diverge again...

- Tony

### Joe Horn

Oct 9, 2004, 3:11:11 AM10/9/04
to
Veky wrote:

> I want my 31bucks. :-]

You certainly earned it, since (a) you answered the challenge
precisely, and (b) in doing so you taught me something that I was dying
to know! THANK YOU!

Do you (or any other readers here) know how Mathematica's "Rationalize"
function zeros in on intermediate convergents? It astounds me that it
does it at all, let alone so quickly.

Rodger: To avoid the confusion caused by my documentation's manifest
vagaries, I'm seriously thinking of replacing the ENTIRE documentation
for PDQ2 with this single line from Mathematica's help screen:
"Rationalize[x, dx] yields the rational number with the smallest
denominator that lies within dx of x." That's exactly what PDQ2 does.
Any further "explanation" might be counterproductive.

-Joe-
"A picture is worth a thousand words, and has a higher baud rate."

### John H Meyers

Oct 9, 2004, 9:13:34 AM10/9/04
to
Joe Horn wrote:

> (a) What is the simplest fraction
> (smallest possible numerator and denominator)

Do you mean just "smallest possible denominator"?

Did you get your patent yet?

Interesting historical debate:
NEW2Q.TXT vs. DEC2FRAC on HP48 Goodies Disk #3
http://www.hpcalc.org/search.php?query=new2q
http://www.hpcalc.org/details.php?id=234
http://www.hpcalc.org/hp48/compilations/horn/horn3.zip

Want to see some remarkable "Mathematical Snapshots"?

Go here and keep clicking "next page":

About the author (quoted from source below):

[Hugo] Steinhaus is best known for his book Mathematical Snapshots written in 1937.
Kac, writing in [7] says:-

... to understand and appreciate Steinhaus's mathematical style,
one must read (or rather look at) snapshots. ... designed to appeal
to "the scientist in the child and the child in the scientist" ...
it expresses, not always explicitly and at times even unconsciously,
what Steinhaus thought mathematics is and should be. To Steinhaus
mathematics was a mirror of reality and life much in the same way
as poetry is a mirror, and he liked to "play" with numbers, sets,
and curves, the way a poet plays with words, phrases, and sounds.

http://www-gap.dcs.st-and.ac.uk/~history/Mathematicians/Steinhaus.html

"A mathematician, like a painter or a poet, is a maker of patterns.
If his patterns are more permanent than theirs,
it is because they are made with ideas.
The mathematician's patterns, like the painter's or the poet's,
must be beautiful; the ideas, like the colours or the words,
must fit together in a harmonious way. Beauty is the first test."

- G. H. Hardy, A Mathematician's Apology

With best wishes from http://www.mum.edu

### Tony Hutchins

Oct 10, 2004, 12:33:50 AM10/10/04
to
-=[ Sun, 10.10.04 5:20 p.m. +1300 (NZDT) ]=-

Hi Joe Horn ! <joe...@holyjoe.net>

21h09m ago, on Sat, 9.10.04 at 00:11 a.m. -0700, you wrote

> Do you (or any other readers here) know how Mathematica's
> "Rationalize" function zeros in on intermediate
> convergents? It astounds me that it does it at all, let
> alone so quickly.

This code by David Eppstein from Aug 1993:

http://www.ics.uci.edu/~eppstein/numth/frap.c

was posted on sci.math by him in March 1995 and he states
there that "similar routines are built into Mathmatica ..."

He's moving in a similar direction to you, and looks at some
intermediate convergents.

Enjoy :)

--
Tony Hutchins

### Joe Horn

Oct 10, 2004, 5:35:12 AM10/10/04
to
> in fact an Idiot Savant friend of mine wants to know
> what is the "simplest fraction" equal to Pi +/- 13131
> parts in 10^440 and he wants to know how long it takes
> PDQ2 to find it. It took him less than a second! Way less!

Hey, that's an excellent example of an intermediate convergent (using
10389, close to half of the partial quotient 20776). Is your friend's
name "Mathematica"? If so, you're right; he finds the answer to be the
following (broken into groups to keep Google from choking to death):

1975891706368 5406240845438 0413327813798 8557337219023 6919811855516
7226743047306 6290670359362 0215835931889 2308274160360 1397933071609
0096564056017 1119521293561 5317285063233 0284830147063 7551101789451
7380003505989 8203820427519

/

62894586416 56661036380 99443296751 77193853240 54623832389 70103364398
48926113959 96646403296 10522734293 18178568324 25836690143 45452239256
64431830927 38965162183 77776034085 40500659700 57153850182 98465893270
92348744070 13579584688

as fast as I could press the Enter key.

PDQ2, written in UBASIC, also gets the same result in no measurable
time. And it says that the diff/tolerance is
0.999987905955071639255847132....

I do not know how long the 49g+ would take, running the posted User RPL
program. Probably several minutes!

>> The inputs and outputs can be of any desired length...

>
> I don't think you really meant to say that.
> Can they be 10^20 digits long? Looks like imprecise speaking to me.

"I do not, my dear Sir, conceive you to be of that sophistical,
captious spirit, or of that uncandid dullness, as to require, for every
general observation or sentiment, an explicit detail of the correctives
and exceptions which reason will presume to be included in all the
general propositions which come from reasonable men." -- Edmund Burke,
1790.

Please also note that Mathematica's documentation even uses the same
"imprecise speaking" as I did. It says, "Mathematica can handle
approximate real numbers with any number of digits" and
"Arbitrary-precision numbers can contain any number of digits". See
both of these general propositions here:
-Joe-

### Joseph K. Horn

Oct 10, 2004, 5:49:28 AM10/10/04
to
Tony Hutchins wrote:

> http://www.ics.uci.edu/~eppstein/numth/frap.c
>
> was posted on sci.math by him in March 1995 and he states
> there that "similar routines are built into Mathmatica ..."
>

> He's moving in a similar direction to you...

Thanks, but please note that his program requires a *maximum
denominator* (not a tolerance), and uses exactly the same logic as
does my 'DEC2FRAC' program on Goodies Disk #3 (21 March 1991). That's
a simpler task than working towards a given tolerance, which is what
both PDQ2 and Mathematica's "Rationalize" do.

-Joe-

### Joseph K. Horn

Oct 10, 2004, 6:07:31 AM10/10/04
to
John H Meyers wrote:

> Joe Horn wrote:
>
> > (a) What is the simplest fraction
> > (smallest possible numerator and denominator)
>
> Do you mean just "smallest possible denominator"?

Yes. Gosper uses the term "simplest" in that sense, so I've adopted
it too. "Smallest possible numerator and denominator" is how I think
of it, but "smallest possible denominator" (or numerator) is
equivalent.

> Did you get your patent yet?

I suspect that my name will be somewhere on the patent (if one is
obtained), but I'm no longer in control of that aspect of PDQ, nor am
I at liberty to discuss the details. If Mathematica's "Rationalize"
function essentially uses the PDQ algorithm, then Wolfram has
"pre-existing art" that reduces all my work to a mere reinvention of
the patently obvious. :-(

-Joe-

### Tony Hutchins

Oct 10, 2004, 6:32:58 AM10/10/04
to
-=[ Sun, 10.10.04 11:18 p.m. +1300 (NZDT) ]=-

Hi Joe Horn ! <joe...@holyjoe.net>

3 days 19m ago, on Thu, 7.10.04 at 02:58 a.m. -0700, you wrote

> FRAC?
> @ Input: obj
> @ Output: obj, 0/1 (1 if obj is 'integer/integer', else 0)
> \<<
> IF DUP TYPE 9. -
> THEN 0.
> ELSE DUP \->LST DUP 1.

I assume that \->LST should be \->LIST
Otherwise an algebraic input results in '->LST' on the stack
and DOSUBS has a problem.

But here a \->LIST doesn't work either because the stack
doesn't have a dimension in level 1.

I'm trying it on a 49g after direct transfer from your
message.

> \<< TYPE
> \>> DOSUBS R\->I { 28 28 18 } SAME SWAP DUP SIZE GET
> { / } 1. GET SAME AND
> END
> \>>
> END
>

--
Tony Hutchins
Wellington
New Zealand

#16 Beginnings and endings are truly artificial constructs.

### Rodger Rosenbaum

Oct 10, 2004, 10:48:40 AM10/10/04
to

I doubt he was criticizing mathematical exposition.

>
>Please also note that Mathematica's documentation even uses the same
>"imprecise speaking" as I did.

"An offense may be common yet remain an offense"; Cromwell to Thomas More in "A Man For
All Seasons"

You should strive to do better than Wolfram, and not try to justify your own
shortcomings by pointing to the shortcomings of others.

Let's see how HP did it in the not too distant past. On page 82 of the HP71 Math Pac
owners manual, bottom of the page: "...solve any number of different systems, limited
only by memory...". Or page 134, bottom of the page: "The number of data points you can
use is limited only by the amount of available memory...". Maybe you've been infected by
the increasing use of hyperbole that seems common in the 21st century.

It says, "Mathematica can handle
>approximate real numbers with any number of digits" and
>"Arbitrary-precision numbers can contain any number of digits".

As you acknowledge above, this is "imprecise speaking". At least one person at
Wolfram knows the truth of it. In section 1.4.9 of the Book, you will find: "You should
have no trouble working out Expand[(1+x)^100] on any computer that can run Mathematica.
But as you increase the exponent of (1+x), the results you get will eventually become to
big for your computer's memory to hold."

### Tony Hutchins

Oct 10, 2004, 3:22:08 PM10/10/04
to
-=[ Mon, 11.10.04 08:03 a.m. +1300 (NZDT) ]=-

Hi Tony Hutchins ! <jus...@csi.com>

08h30m ago, on Sun, 10.10.04 at 11:32 p.m. +1300 (NZDT), you wrote
in message ID <1154352...@news.individual.net> :

> > ELSE DUP \->LST DUP 1.
>
> I assume that \->LST should be \->LIST

This seems to work instead of the \->LST:

OBJ\-> SWAP 1 + \->LIST

--
Tony Hutchins

### Tony Hutchins

Oct 10, 2004, 7:23:08 PM10/10/04
to
-=[ Mon, 11.10.04 11:17 a.m. +1300 (NZDT) ]=-

Hi Joe Horn ! <joe...@holyjoe.net>

12h42m ago, on Sun, 10.10.04 at 02:35 a.m. -0700, you wrote

> PDQ2, written in UBASIC, also gets the same result in no
> measurable time. And it says that the diff/tolerance is
> 0.999987905955071639255847132....

Ahha we see a memory limitation at last ;-) No no no, I see
you mean the error(diff) is 0.9999879.. times the precision
(tolernace) of 13131 E-440.

> I do not know how long the 49g+ would take, running the posted User RPL
> program. Probably several minutes!

Yup about 10 or so. It gives the same result you quoted for
the rational approximation, but unfortunately the numerator
and denominator in the error both exceed E500 so direct
floating point division gives 1. But, just using the first few
significant digits, and roughly measuring the lengths as 796
(numerator) and 1232 (denominator) - based on PC file-lengths
of the integers as strings - the difference is 436, I can see
the error is just under 13131 E-440 ;-) Probably something
like 1.313084 E-436. Ahha, yes it agrees with your UBASIC
result.

Pretty impressive! It took about 30-40 mins on the 49g. Then I
finally got the new connx working again and tried on the g+.
On win2K I had all sorts of problems installing the 30 March
2004 connx from HP's site - it kept asking for the CD to get
the driver - refused to use the new drivers - I read
all the readme's. In the end I just fed it the CD and used the
old driver, which it thought was "more recent".

I used [PI1000] instead of [PI64] ;-)

--
Tony Hutchins

### Tony Hutchins

Oct 11, 2004, 12:14:07 AM10/11/04
to
-=[ Mon, 11.10.04 5:03 p.m. +1300 (NZDT) ]=-

Hi Tony Hutchins ! <jus...@csi.com>

04h40m ago, on Mon, 11.10.04 at 12:23 p.m. +1300 (NZDT), you wrote
in message ID <1081179...@news.individual.net> :

> I can see the error is just under 13131 E-440 ;-) Probably
> something like 1.313084 E-436. Ahha, yes it agrees with

Ahha, I see the 49g can do this. BIGDIV below gives: -436 and
1.31308411931. But I'm sure it can be done better than
this<G>.

%%HP: T(3)F(.);
@ BIGDIV: Handy for evaluating error from PDQ2
@ where N or D exceeds E500.
@ Input: 'N/D'
@ Output: E,M (Exponent-like,Mantissa-like)
@ Where Evaluated error=M*10^E
\<<
OBJ\-> DROP2 DUP2
SWAP SIZE SWAP SIZE -
UNROT
SWAP \->STR 1 499 SUB STR\-> I\->R
SWAP \->STR 1 499 SUB STR\-> I\->R
/
\>>

--
Tony Hutchins
Wellington
New Zealand

#347 "Facts are a special form of fantasy." - Michael Lukas Moeller

### Joseph K. Horn

Oct 11, 2004, 12:33:43 PM10/11/04
to
Tony Hutchins wrote:

> > FRAC?
> > ...

> > ELSE DUP \->LST DUP 1.
>
> I assume that \->LST should be \->LIST

->LST is a command in library 256. It can convert an algebraic object
directly into list form. If you don't have 256 ATTACH in your
'STARTUP' program, you'll have to manually type 256 ATTACH before
keying in that program. Sorry; I am so used to having library 256
available that I forgot to mention the need for it here.

-Joe-

### Tony Hutchins

Oct 11, 2004, 4:08:51 PM10/11/04
to
-=[ Tue, 12.10.04 08:35 a.m. +1300 (NZDT) ]=-

Hi Joseph K. Horn ! <joe...@holyjoe.net>

03h01m ago, on Mon, 11.10.04 at 09:33 a.m. -0700, you wrote

> Tony Hutchins wrote:
[...]

> > I assume that \->LST should be \->LIST
>
> ->LST is a command in library 256.

Thanks Joe! BTW your CF~I2 gets the 432 items in the CF for
that PI rational approximation with tolerance 13131 E-440 in
20 seconds. I wish I had Mathematica here - I'd try to find
an example where it produces a different result to PDQ2.

--
Tony Hutchins

### Rodger Rosenbaum

Oct 11, 2004, 6:22:27 PM10/11/04
to

Why would you expect to find such an example? As long as the 49 doesn't run out of
memory, I would expect PDQ2 to work properly. Joe is a good programmer, even if his
documentation is some what enthusiastic.

FYI, these interesting tests of PDQ2 occur where the partial quotients of the CF of PI
are big. After the one at 432 terms (the PQ is 20776 there), the next four in the
expansion of Pi are:

Partial Quotient Location

78629 28422
179136 156382
528210 267314
12996958 453294

It's beginning to take even Mathematica a while (minutes) to find these.

### Tony Hutchins

Oct 11, 2004, 8:45:28 PM10/11/04
to
-=[ Tue, 12.10.04 1:06 p.m. +1300 (NZDT) ]=-

Hi Rodger Rosenbaum ! <nos...@aol.com>

01h44m ago, on Mon, 11.10.04 at 3:22 p.m. -0700, you wrote
in message ID <l81mm0tsbjatc6smd...@4ax.com> :

> >> ->LST is a command in library 256.
> >
> >Thanks Joe! BTW your CF~I2 gets the 432 items in the CF for
> >that PI rational approximation with tolerance 13131 E-440 in
> >20 seconds. I wish I had Mathematica here - I'd try to find
> >an example where it produces a different result to PDQ2.
>
> Why would you expect to find such an example?

Did I say that?<G>. It would indeed be surprising to find one.

> As long as the 49 doesn't run out of memory, I would expect
> PDQ2 to work properly. Joe is a good programmer, even if
> his documentation is some what enthusiastic.

Yes, I'm sure PDQ2 works as I can see the algorithm in source
code. And we have documentation ;-). The Mathematica algorithm
for Rationalize[x,dx] does not seem to be published on the net.
AFAIK the only reason we have for thinking it is the same as
Joe's is its surprising speed.

> FYI, these interesting tests of PDQ2 occur where the partial quotients of the CF of PI
> are big. After the one at 432 terms (the PQ is 20776 there), the next four in the
> expansion of Pi are:
>
> Partial Quotient Location
>
> 78629 28422
> 179136 156382
> 528210 267314
> 12996958 453294
>
> It's beginning to take even Mathematica a while (minutes)
> to find these.

Thanks Roger - very interesting. That last one is 1000*
further out than the 432nd and may take 6 hours or so on the
49g+. Hmm I might see how far CF~I2 can get before we hit
the space-time limitations ;-) Oops I mean memory-battery_life
limitations.

Joe - all your documentation needs at the beginning is a
little poem about memory and battery-life limitations ;-)

Roger - I just twigged that you MUST be "The Phantom"!! It
could only be you ;-) LOL! Excellent example BTW - thanks.
Nice to see your "idiot savant" friend finally taking
measurable minutes to compute though - was he computing
Rationalize[x,dx] out near these peaks, or just doing the
prolonged CF?

--
Tony Hutchins

### Rodger Rosenbaum

Oct 11, 2004, 10:04:06 PM10/11/04
to
On Tue, 12 Oct 2004 13:45:28 +1300 (NZDT), Tony Hutchins <jus...@csi.com> wrote:

> -=[ Tue, 12.10.04 1:06 p.m. +1300 (NZDT) ]=-
>
>Hi Rodger Rosenbaum ! <nos...@aol.com>
>
>01h44m ago, on Mon, 11.10.04 at 3:22 p.m. -0700, you wrote
>in message ID <l81mm0tsbjatc6smd...@4ax.com> :
>
>> >> ->LST is a command in library 256.
>> >
>> >Thanks Joe! BTW your CF~I2 gets the 432 items in the CF for
>> >that PI rational approximation with tolerance 13131 E-440 in
>> >20 seconds. I wish I had Mathematica here - I'd try to find
>> >an example where it produces a different result to PDQ2.
>>
>> Why would you expect to find such an example?
>
>Did I say that?

A few lines above you said: "I wish I had Mathematica here - I'd try to find

an example where it produces a different result to PDQ2."

Why would you bother if you didn't expect (at least a little bit) to find such an
example?

Maybe you're a little like me. When I got my first HP calculator (an HP 45) I used to go
through log tables and trig tables and pick a value at random to see if the HP got the
right result. I didn't really expect it to fail, but it just amazed me (at that time; I'm
spoiled now) to see all those digits match the numbers in the table. And this was in a
book of 7 place (!!!) tables; the HP got 10 digits!!! WOW!!

So I understand how much fun it is to see your HP do all this stuff almost as well as the
big boys. Joe would probably tell us that the Hurwitz accuracy of the HP calcs is
superior since you are getting so much more "bang" for the buck. The HP runs on penlight
batteries and can still do this kind of thing in a reasonable time.

I was looking at a paper in the IEEE Transactions on Circuit Theory from the 1950's where
the author (at Bell Labs) discussed the need to find the roots of a 12th degree polynomial
and separate the roots with positive real parts from those with negative real parts. She
said that one way would be to just find the roots and manually separate them, but she went
on to describe this really convoluted algorithm that she said gave better performance. I
typed in the coefficients on my HP48G and found the roots in 21 seconds, saved them to a
variable, OBJ-> on the stack, deleted the ones with positive real parts, saved them, did
it again to get the ones with negative real parts. This took about 30 seconds each. She
said her fancy algorithm took 80 (!) minutes on an IBM 650, and her results weren't as
accurate as mine! It's fun and good practice on your calc to see how fast and easy it is
to redo things like that from old journals; and I usually get better results.

Have you looked at my nearby "big numbers" MC? Another example where the HP exhibits
its astounding performance. But, as I said, I'm spoiled now and wouldn't expect anything
less.

Just a prolonged CF. I've been trying to find some irrational number that has a really
huge PQ after only a few terms, but I think it will be necessary to use come clever
artifice to get one.

The "best" one I found so far is 7^(1/5); it has a PQ of 55318 at position 125 in the
expansion.

### Rodger Rosenbaum

Oct 11, 2004, 10:09:31 PM10/11/04
to
On Tue, 12 Oct 2004 13:45:28 +1300 (NZDT), Tony Hutchins <jus...@csi.com> wrote:

Snip...

>> >Thanks Joe! BTW your CF~I2 gets the 432 items in the CF for
>> >that PI rational approximation with tolerance 13131 E-440 in
>> >20 seconds. I wish I had Mathematica here - I'd try to find
>> >an example where it produces a different result to PDQ2.

I know why you want to do that!! You want to find a bug so you can fix it!! :=)

### Tony Hutchins

Oct 12, 2004, 5:21:43 AM10/12/04
to
-=[ Tue, 12.10.04 9:43 p.m. +1300 (NZDT) ]=-

Hi Rodger Rosenbaum ! <nos...@aol.com>

06h39m ago, on Mon, 11.10.04 at 7:04 p.m. -0700, you wrote
in message ID <mtcmm0p30iejlo6i1...@4ax.com> :

> A few lines above you said: "I wish I had Mathematica here
> - I'd try to find an example where it produces a different
> result to PDQ2."
>
> Why would you bother if you didn't expect (at least a little
> bit) to find such an example?

I think I was hoping that a difference in algorithms would
somehow show - maybe a case that took Mathematica 2 seconds
where it should take 1 :-)

> Maybe you're a little like me. When I got my first HP
> calculator (an HP 45) I used to go through log tables and
> trig tables and pick a value at random to see if the HP
> got the right result. I didn't really expect it to fail,
> but it just amazed me (at that time; I'm spoiled now)
> to see all those digits match the numbers in the table.
> And this was in a book of 7 place (!!!) tables; the HP got
> 10 digits!!! WOW!!

Indeed. My first HP was an HP 80. The financial tables were
then useful only to test the HP-80 ;-)

> So I understand how much fun it is to see your HP do all
> this stuff almost as well as the big boys.

And we didn't know till recently that Mathematica did the
Rationalize[x,dx].

> Joe would probably tell us that the Hurwitz accuracy of the
> HP calcs is superior since you are getting so much more
> "bang" for the buck. The HP runs on penlight batteries
> and can still do this kind of thing in a reasonable time.

That means the hp49g+ probably has the highest "HP" Hurwitz
Factor. It's quite interesting to find its limitations. I
tried a PI to 5000 digits this afternoon and naturally got a
partial quotients in the stack. Also, using the 'PI5000'
really slows down that 13131 E-440 example - it took over an
hour to get the result that took just a few minutes using
'PI1000'.
[...]

> Have you looked at my nearby "big numbers" MC? Another
> example where the HP exhibits its astounding performance.
> But, as I said, I'm spoiled now and wouldn't expect
> anything less.

I will have a closer look.

> >measurable minutes to compute though - was he computing
> >Rationalize[x,dx] out near these peaks, or just doing the
> >prolonged CF?
>
> Just a prolonged CF. I've been trying to find some
> irrational number that has a really huge PQ after only a
> few terms, but I think it will be necessary to use come
> clever artifice to get one.

Very interesting challenge.

> The "best" one I found so far is 7^(1/5); it has a PQ of
> 55318 at position 125 in the expansion.

On the native 49g+ I only get to position 28 in the CF
expansion of that one.

--
Tony Hutchins

### Rodger Rosenbaum

Oct 12, 2004, 6:42:06 AM10/12/04
to

Joe only had a PI64 function in PDQ2. What have you done to get a
PI1000 and PI5000? Did you compute it on the 49 itself, or get the
digits from another, bigger, computer and import into the 49?

>
>> Have you looked at my nearby "big numbers" MC? Another
>> example where the HP exhibits its astounding performance.
>> But, as I said, I'm spoiled now and wouldn't expect
>> anything less.
>
>I will have a closer look.
>
>> >measurable minutes to compute though - was he computing
>> >Rationalize[x,dx] out near these peaks, or just doing the
>> >prolonged CF?
>>
>> Just a prolonged CF. I've been trying to find some
>> irrational number that has a really huge PQ after only a
>> few terms, but I think it will be necessary to use come
>> clever artifice to get one.
>
>Very interesting challenge.
>
>> The "best" one I found so far is 7^(1/5); it has a PQ of
>> 55318 at position 125 in the expansion.
>
>On the native 49g+ I only get to position 28 in the CF
>expansion of that one.

What was the limiting factor? Why could you only get to position
28?

>

### Tony Hutchins

Oct 12, 2004, 12:44:50 PM10/12/04
to
-=[ Wed, 13.10.04 04:23 a.m. +1300 (NZDT) ]=-

Hi Rodger Rosenbaum ! <nos...@aol.com>

04h41m ago, on Tue, 12.10.04 at 05:42 a.m. -0500, you wrote
in message ID <15dnm0purih79nvf3...@4ax.com> :

[...]

> >Also, using the 'PI5000' really slows down that 13131
> >E-440 example - it took over an hour to get the result
> >that took just a few minutes using 'PI1000'. [...]
>
> Joe only had a PI64 function in PDQ2. What have you done
> to get a PI1000 and PI5000? Did you compute it on the 49
> itself, or get the digits from another, bigger, computer
> and import into the 49?

The latter. Thanks to Google ;-)
[...]

> >> The "best" one I found so far is 7^(1/5); it has a PQ of
> >> 55318 at position 125 in the expansion.
> >
> >On the native 49g+ I only get to position 28 in the CF
> >expansion of that one.
>
> What was the limiting factor? Why could you only get to
> position 28?

I started with the 12 digit 1.47577316159. Adding another 4
digits from the 200LX palmtop didn't help. I'd probably need at
least 100 digits to find the 55318.

--
Tony Hutchins

### Gjermund Skailand

Oct 12, 2004, 1:27:10 PM10/12/04
to
Hi,
Rodger Rosenbaum <nos...@aol.com> wrote in message news:<15dnm0purih79nvf3...@4ax.com>...

> On Tue, 12 Oct 2004 22:21:43 +1300 (NZDT), Tony Hutchins
> <jus...@csi.com> wrote:
>
> > -=[ Tue, 12.10.04 9:43 p.m. +1300 (NZDT) ]=-
> >
> >Hi Rodger Rosenbaum ! <nos...@aol.com>
> >
> >06h39m ago, on Mon, 11.10.04 at 7:04 p.m. -0700, you wrote
> >in message ID <mtcmm0p30iejlo6i1...@4ax.com> :
> >
> >> A few lines above you said: "I wish I had Mathematica here
> >> - I'd try to find an example where it produces a different
> >> result to PDQ2."
> >>
> Joe only had a PI64 function in PDQ2. What have you done to get a
> PI1000 and PI5000? Did you compute it on the 49 itself, or get the
> digits from another, bigger, computer and import into the 49?

Sorry to pop in, hope this is relevant:
Just for information: my longfloat library at hpcalc.org has functions
for PI to arbitrary length and conversion between floating / simple
fraction / continuous fraction representation. Could be handy when
doing this kind of calculations.
Best regards
Gjermund Skailand
PS Try continuous fraction expansion of this fraction ( 1.3 sec on
hp49g+ with longfloat lib):

1829151464613551845229766684240085457977114066958132611028137640130366249249264703822076231310349502716995004329909838833967799149726538729584636544830491123626011712281251039707242014594449383181594119455905756335372184341285442894833577777675820454061486809900044466151484538203903862528738851534878481429062476876438743657726379178432129015903949079323534884585879968894630898002489789447203843435
000
/
1000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000

Continuous fraction expansion:
{ 1829 6 1 1 1 1 17 1 1 3 2 2 3 7 1 7 1 428 1 1 3 52 1 5 1 4 1 2 1 1 2
4 3 3 8 1 1 2 1 2 7 5 1 17 2 5 2 3 4 8 2 5 1 13 35 5 1 2 2 1 102 1 1 3
2 1 6 1174432375356999418156906759746231648661456610679530520916425437057426881793642718189281666359914553831730104240207954136420939696280455479951978198341209475388354716281725861581739201751862049650448139929231944021896869379111614535260523870011376825005238855165204139267124381959584353592123424607279088687963807804342933494
1 1 2 1 1 3 41 1 5 1 1 3 1 1 1 13 1 1 1 1 1 1 2 10 2 2 1 2 2 1 2 1 9 4
1 1 1 8 1 6 2 1 13 315 2 2 2 2 2 4 323 28 2 12 5 1 1 1 3 2 2 2 2 2 3
15 1 5 207 2 }

### Tony Hutchins

Oct 12, 2004, 6:43:30 PM10/12/04
to
-=[ Wed, 13.10.04 11:11 a.m. +1300 (NZDT) ]=-

Hi Gjermund Skailand ! <gjermund...@yahoo.no>

04h43m ago, on Tue, 12.10.04 at 10:27 a.m. -0700, you wrote

> > Joe only had a PI64 function in PDQ2. What have you done to get a
> > PI1000 and PI5000? Did you compute it on the 49 itself, or get the
> > digits from another, bigger, computer and import into the 49?
>
> Sorry to pop in, hope this is relevant:

It is, it is! :)

> Just for information: my longfloat library at hpcalc.org
> has functions for PI to arbitrary length and conversion
> between floating / simple fraction / continuous fraction
> representation. Could be handy when doing this kind of
> calculations.

Thanks Gjermund! I just tried Roger's 7^(1/5) with DIGITS=200
(took under a minute) in your library and then [Fl->SF]
[SF->CF] and can see Roger's 55318 at position 125. What fun!
Mirabile dictu!

Next I will experiment with the [FPI] to a few thousand digits
to test the limits of what is possible on the 49g+. I will
work back from DIGITS=5000, where I got "positive underflow"
after a few minutes :) No chance to find Roger's PQ of 78629
at position 28422.

> Best regards
> Gjermund Skailand
> PS Try continuous fraction expansion of this fraction ( 1.3 sec on
> hp49g+ with longfloat lib):

Amazing. I see a big partial quotient early in the list. What
is the origin of the number input? Is it irrational?

BTW minor thing - CF= "continued", rather than
"continuous", fraction. It's in the longfloatv31b.doc a couple
of times near the end. I tried to go to www.praxelius.de as
mentioned on page 1 of the .doc but alas it was not
contactable.

Cheers - I look forward to your matrix enhancement in
longfloat.

--
Tony Hutchins

### Veli-Pekka Nousiainen

Oct 13, 2004, 6:49:30 AM10/13/04
to
"Tony Hutchins" <jus...@csi.com> wrote in message
news:1076850...@news.individual.net...
X

> Cheers - I look forward to your matrix enhancement in
> longfloat.
This would be one of the most important additions to the 49g+
When a matrix is ill-conditioned in such a way
that the elements have exponents range exeeding
even the internal 15 digit accuracy, it's nice to have
1000 digits of floating point range to (possibly) overcome
this lost accuracy during matrix inversions, etc...

Are you going to use
A) the "old" symbolic matrices (internally as a list)
B) a new type of SymMat with slighly different prolog
C) just use lists and let people set the proper flag (-91)
in order to use the Matrix Writer for input
[VPN]
hpgcc 4 ever

### Rodger Rosenbaum

Oct 13, 2004, 9:32:42 PM10/13/04
to

I got one more after an all night run:
878783625 11504931

### Joseph K. Horn

Oct 13, 2004, 9:53:27 PM10/13/04
to
Rodger Rosenbaum writes:

> I've been trying to find some irrational number
> that has a really huge PQ after only a few terms,
> but I think it will be necessary to use come clever
> artifice to get one.

A decidedly UNclever technique is a web search, but it reveals some
goodies, such as these (all on the first page of an AltaVista search
on "partial quotients" and "pi"):

(A) (1-euler)^2 has 372792 as term 8 (if the leading term, the integer
part, is called "term zero"). "euler" is Euler's constant, approx.
0.5772156649.

(B) pi^4 which has 16539 as term 5 (if the leading term is considered
"term zero").

(C) "The number known to recreational mathematicians as Champernowne
Constant has a decimal expansion consisting of the successive digits
of all the integers: 0.123456789101112131415161718192021... It has a
weird continued fraction expansion, which is fairly delicate to
obtain:
[0; 8, 9, 1, 149083, 1, 1, 1, 4, 1, 1, 1, 3, 4, 1, 1, 1, 15, A18, 6,
1, 1, 21, 1 ...]
A18 is a 166 digit integer!"

(D) "The real root of x3 -3600 x2 +120 x -2 :
x = [3599; 1, 28, 1, 7198, 1, 29, 388787400, 23, 1, 8998, 1, 13,
1, 10284, 1, 2, 25400776804, 1, 1, 3, 4, 93, 3, 1, 2, 11, 1, 9, 1, 99,
1, 3, 1, 3, 9, 1, 603118914 ...]"

(E) "the real zero of y3 - 8y - 10 (y = 3.31862821775...) :
y = [3; 3, 7, 4, 2, 30, 1, 8, 3, 1, 1, 1, 9, 2, 2, 1, 3, 22986, 2,
1, 32, 8, 2, 1, 8, 55, 1, 5, 2, 28, 1, 5, 1, 1501790, 1, 2, 1, 7, 6,
1, 1, 5, 2, 1, 6, 2, 2, 1, 2, 1, 1, 3, 1, 3, 1, 2, 4, 3, 1, 35657, 1,
17, 2, 15, 1, 1, 2, 1, 1, 5, 3, 2, 1, 1, 7, 2, 1, 7, 1, 3, 25, 49405,
1, 1, 3, 1, 1, 4, 1, 2, ...]"

-Joe-

### Gjermund Skailand

Oct 13, 2004, 5:34:31 PM10/13/04
to
"Veli-Pekka Nousiainen" <DROP...@THIS.welho.com> wrote in message news:<ckj17q\$3ck\$1...@nyytiset.pp.htv.fi>...

> "Tony Hutchins" <jus...@csi.com> wrote in message
> news:1076850...@news.individual.net...
> X
> Are you going to use
> A) the "old" symbolic matrices (internally as a list)
> B) a new type of SymMat with slighly different prolog
I have uploaded a version 3.5b to hpcalc.org today. The library uses
the old arry type, so there are some restrictions. You have to see the
manual for that. This is simpler and a bit faster, but not really
supported by the OS. Thus it will always be a beta-version.

Tony: A nice manual is partly included in the new package. If you want
to calculate eg. 7^(1/5) use FXROOT it is much faster that FY^X which
uses logarithms, for 200 Digits the timing of FXROOT is about 3 sec.
Some constants like LN2, PI are kept in library data, to avoid
unnecesary recalculation. (Try FXROOT a second time).
Regarding "continued fractions", the prevously mentioned fraction is a
solution to a recent post, using both exact and longreal matrix
calculation, calculated with 400 Digits. Send me a mail if you're
interested in more details.

Best regards
Gjermund Skailand
PS Using the emu49 I was able to calculate PI at upper limit to 9999
digits with version 3.5.

### Tony Hutchins

Oct 13, 2004, 9:10:46 PM10/13/04
to
-=[ Thu, 14.10.04 1:26 p.m. +1300 (NZDT) ]=-

Hi Gjermund Skailand ! <gjermund...@yahoo.no>

02h52m ago, on Wed, 13.10.04 at 2:34 p.m. -0700, you wrote

> Tony: A nice manual is partly included in the new package.

I'll look out for 3.5b on hpcalc - just checked - it'snot
there yet.

> If you want to calculate eg. 7^(1/5) use FXROOT it is much
> faster that FY^X which uses logarithms, for 200 Digits the
> timing of FXROOT is about 3 sec.

Great. Thanks!

> Some constants like LN2, PI are kept in library data, to
> avoid unnecesary recalculation. (Try FXROOT a second time).

I noticed :). Very nice feature.

> Regarding "continued fractions", the prevously mentioned
> fraction is a solution to a recent post, using both exact and
> longreal matrix calculation, calculated with 400 Digits. Send
> me a mail if you're interested in more details.

Ahha. Thanks

> PS Using the emu49 I was able to calculate PI at upper
> limit to 9999 digits with version 3.5.

Ahha that 9999 I keep coming across as an upper limit for some
operations. I look forward to trying that with v3.5 on my 49g+.

--
Tony Hutchins

### Joseph K. Horn

Oct 19, 2004, 7:55:24 PM10/19/04
to
Hey, Rodger! I finally sat down and started to learn Mathematica!
(Ok, ok, so maybe I'm a little slow on the uptake. Better late than
never!)

There are an infinite number of ways of generating increasingly large
partial quotients, of course, but here's one that tickles me: the 3rd
partial quotient of the cosines of the increasingly "better" fractions
for pi. So easy to do in Mathematica! Out[2] below is the list of
numerators (of 3/1, 22/7, etc) and Out[3] is the list of the 3rd
partial quotients of their cosines.

In[1]:=
<< NumberTheory`ContinuedFractions`

In[2]:=
Numerator[Convergents[Pi, 25]]

Out[2]=
{3, 22, 333, 355, 103993, 104348, 208341, 312689, 833719, 1146408,
4272943, \
5419351, 80143857, 165707065, 245850922, 411557987, 1068966896,
2549491779, \
6167950454, 14885392687, 21053343141, 1783366216531, 3587785776203, \
5371151992734, 8958937768937}

In[3]:=
ContinuedFraction[Abs[Cos[%]], 3][[All, 3]]

Out[3]=
{98, 25526, 25701, 2200989908, 5465503977, 16483886140, 30375674292, \
237697464198, 373859991506, 5788958372881, 6621691598463,
1370542374992698, \
9164348719628591, 26700401132899689, 53432041505819839,
310803772468965157, \
1832746004640389165, 9989458647919860071, 89204641469760295779, \
91331082019550175280, 650545454929017880194800,
4117455774884563343654118, \
15476430462882332335591161, 17562023095942885732083182, \
4124798882403720305240845429}

The sines' 2nd PQ likewise increase without bound as the input
approaches a multiple of pi. Fun! Mathematica is so fast!!!

However, Mathematica's "Rationalize" function takes about 2.5 times
LONGER than my UBASIC program to perform PDQ2 on 2000-digit inputs.
Strange. Methinks they missed a fundamental shortcut somewhere.

-Joe-

### Rodger Rosenbaum

Oct 20, 2004, 6:58:23 AM10/20/04
to
On 19 Oct 2004 16:55:24 -0700, joe...@holyjoe.net (Joseph K. Horn) wrote:

>Hey, Rodger! I finally sat down and started to learn Mathematica!
>(Ok, ok, so maybe I'm a little slow on the uptake. Better late than
>never!)
>
>There are an infinite number of ways of generating increasingly large
>partial quotients, of course, but here's one that tickles me: the 3rd
>partial quotient of the cosines of the increasingly "better" fractions
>for pi.

I'm puzzled; here you seem to be speaking of the convergents themselves--

So easy to do in Mathematica! Out[2] below is the list of
>numerators (of 3/1, 22/7, etc) and Out[3] is the list of the 3rd
>partial quotients of their cosines.

But here you are using the *numerators* of the convergents. If you do this same thing
with the convergents themselves rather than their numerators, the 3rd PQ gets bigger a lot
faster!

When I run Timing[Rationalize[Pi,10^-2000];]

It takes 3.46 seconds under Mathematica 4.1
but only .88 seconds under Mathematica 5.01

Later versions have significantly faster arbitrary precision arithmetic.

This is on my 600 MHz laptop Pentium III
>
>-Joe-

### Bhuvanesh

Oct 20, 2004, 9:19:52 AM10/20/04
to
joe...@holyjoe.net (Joseph K. Horn) wrote:
> However, Mathematica's "Rationalize" function takes about 2.5 times
> LONGER than my UBASIC program to perform PDQ2 on 2000-digit inputs.
> Strange. Methinks they missed a fundamental shortcut somewhere.

Thanks for pointing this out. I haven't tried this yet, but I'll look
into it. Is the source for your UBASIC program available?

Cheers,
Bhuvanesh.

### Joseph K. Horn

Oct 21, 2004, 2:47:22 PM10/21/04
to
Bhuvanesh wrote!

> Is the source for your UBASIC program available?

http://holyjoe.net/hp/pdq2.txt
Just rename it from .TXT to .UB and it'll load directly into UBASIC.

If you want to time huge inputs, be sure to comment out all but the final output.

-Joe-