PI (Product) function for 48/49/50
John H Meyers
we can turn the summation function (Sigma)
into a product function,
at least for numeric arguments and limits:
To calculate the factorial of 253, for example,
(using the built-in factorial function
as a "control" for comparison), we can write:
'\GS(k=1,253,PI(k))' ==> 5.17346099269E499
Where we have applied the PI() function
to the expression whose product we actually desire,
which transforms the "summation" into a "product."
This is only slightly less accurate than the built-in
factorial function, which gives 5.17346099264E499 [correct
to 12 digits, so our error is in the last digit only],
while a method using logarithms would be much less accurate,
yielding 5.17346106391E499, in addition to having
some difficulty with any negative factors.
The use of a function which the built-in compiler understands,
rather than creating a completely new user-defined function,
allows us to use the same syntax for summation and products;
the compiler also knows that it must quote both the first
and last arguments, which it won't do for any user function's
arguments, and this preserves the original symbolic expression
even if EVAL is applied with symbolic limits.
The helper function PI(x) is provided below for 48/49/50 series;
a flaw in this method is that it will not handle
a result which becomes zero before evaluation is complete,
due either to a zero factor or an underflow,
because a previous "running product" value of zero
is what we use to detect that the very first factor
is being presented, at which point we start over again;
invalid limits (e.g. n=1 to -8) are also not acceptable.
HP48[S/G][X][+] series product helper function:
\<< \-> x. \<< x. OVER { * 0 } IFT \>> \>> 'PI' STO
HP49G/50G/48Gii series product helper function (SysRPL):
:: CK1NOLASTWD ( with two up-front safety checks )
GetVStackProtectWord #0=?SEMI
TWO GetElemBotVStack
DUP { x* %0 } xIFT
TWO PutElemBotVStack ;
An alternate way to create this program
via UserRPL (no assembly required :)
"D9D207926233262E7A43B1133422628813047A20"
"8ED93739F2B21304A693B1133A6262B2130" +
IF 64 STWS DUP BYTES DROP #3C00h ==
THEN #100001h LIBEVAL END 'PI' STO
You may notice that a product of Integers (zints)
will produce a "real" result; this would be the case
even if PI were replaced with the *empty* program :: ;
(in which case it would of course revert back to summation)
so it just doesn't pay to try to use ZINT 0 in place of %0
(in most cases it would even slow down the calculation).
Why does using even an *empty* PI program
produce a "real" sum of zints, with decimal point?
This program was greatly facilitated by
JYA's complete explanation of the "Virtual Stack"
(which must be used in the 49/50/48Gii series version):
http://groups.google.com/group/comp.sys.hp48/msg/871ede0942adffe4?dmode=source
http://www.hpcalc.org/details.php?id=2992
http://www.hpcalc.org/hp49/docs/programming/stackapi.zip
With best wishes from http://www.mum.edu
and http://www.maharishischooliowa.org
Fabian
one could expect since there is a superfluous addition involved in each
step.
Cheers,
Fabian
John H Meyers
> Really unconventional approach.
It has many fewer problems, however,
than any "conventional" approach involving a "user" function.
> I fear that the timing is not what one could expect
> since there is a superfluous addition involved in each step.
The extra addition is really almost completely insignificant;
a clue to what is more significant might come from the fact
that on a 48GX, even with a UserRPL "helper" function,
the result computes much faster than on a 49G
having the same CPU speed, even using a SysRPL helper function
(a 49G+ has a much faster effective emulated CPU speed,
but still may have a hard time competing with an old 48GX
in this particular race).
I'm afraid that it's the CAS overhead; it might be fair
to say that no CAS could be a speed demon when not only
forced to use this processor, but also the entire
limited environment it has to run in here;
that may be why an ISOL command (isolate a variable
in an algebraic expression) takes 15 times as long
to execute on a 49G as on my 48GX, which at one time
put a small crimp into my once speedy triangle solver :)
Upon every single invocation of the helper function,
a command known as "xFCNAPPLY" within the overall algebraic
invokes new CAS activities, which include the stashing away
of the entire stack before calling the helper function,
and the restoration of the entire stack afterwards;
you can get an idea of how little the helper function itself
contributes to the overall time by simply replacing it
with an empty program -- this will give you the same
old sum again, rather than a product, but it will also show you
that nearly all the execution time and overhead is spent in the CAS
(the "core" execution time for a mere 252 floating point
multiplications *and* additions amounts to very little).
I will post for you later an independent UserRPL product function,
even though you will see that it poses various other problems
to even attempt to rely on "user" functions for this sort of thing,
even if they were housed in libraries, which is why I devised
this trick of "piggybacking" onto the Summation function instead.
C'ya later.
[r->] [OFF]
John H Meyers
UserRPL programs for the \PI (Product) function,
using two more "really unconventional approaches."
As we've noted earlier, however, when you invoke
a UserRPL program as a "user defined function,"
called from within an algebraic expression,
you must use standard algebraic syntax (all arguments
separated by commas or by semicolons), and you must
*quote* any symbolic arguments yourself.
For example, instead of '\GS(n=1,253,PI(n))'
[where we had transformed Summation into a Product]
you must here write '\PI(QUOTE(n),1,253,QUOTE(n))'
unless you are sure that 'n' EVAL will remain 'n'
and won't find anything previously defined.
Otherwise, the calculator will, by default, attempt
to pre-evaluate all the arguments if at all possible,
before even passing their values to the user function;
if the summation/product variable 'n' already exists,
this means that the symbolic arguments
would not reach our program intact,
resulting in either an error or a wrong answer
(with the built-in Summation function, on the other hand,
the compiler *automatically* quotes the symbolic arguments,
which is one reason why we suggested the unconventional idea
of tricking the Summation function into doing products for us).
To avoid any mixups, note that our earlier "helper" program
to transform the Summation function into a Product function
was named 'PI', while our directly executable stand-alone
Product functions below are named '\PI' (the Greek Capital
letter \PI, Alpha right-shift P on 49G[+]/50G/48Gii
or Alpha right-shift Z on 48G[X][+])
Here is a UserRPL "Product" function based on
the behavior of the internal graph-plotting engine,
which needs to evaluate a formula to be plotted
using changing values of a stored independent variable,
just as we need to evaluate a formula repeatedly using
different values, to get the individual terms in our product.
To avoid potential conflicts with pre-existing variables,
it creates a brand *new* variable, even if an old variable
having the same name already exists! -- of course, it also
traps errors, so that it will be sure to delete the new variable
before quitting, whether or not the plot is prematurely canceled,
in all cases leaving any old variable completely undisturbed:
@ HP48/49/50 [all models]
\<< \-> ?n ?a ?b ?f \<< 1 ?a DUP
?n DUP TYPE 6 \=/ { 514 DOERR } IFT
#8696h SYSEVAL ?b @ create new variable
IFERR FOR ?i ?i ?n STO ?f EVAL * NEXT
THEN 1 ELSE 0 END ?n PURGE
{ ERRN DOERR } IFT \>> \>> '\PI' STO
Try it with this test:
12345 'n' STO @ create a conflicting variable
'\PI(QUOTE(n),1.,253.,QUOTE(n))' EVAL
'n' RCL @ is the original variable still there?
The next version of a UserRPL "Product" function
is based on a radically different idea,
which is that the fastest UserRPL *program*,
rather than a formula, to compute '\PI(n,a,b,f(n))'
would be: \<< 1 'a' 'b' FOR n 'f' EVAL * NEXT \>>
So okay, what we'll do is to actually "manufacture"
that program (in *text* form) from the arguments
we receive, and then we'll compile and run it!
Well, that wasn't my own idea -- it happens to be
the actual way in which the built-in SEQ command
has always worked in the HP48G/49G/50G series,
so if HP is that weird, why can't we be too?
\<<
\-> n a b f \<< n a b f @ HP48[S/G][X][+]
RCLF \-> ?f \<< IFERR STD 64 STWS
ROT \->STR " " + ROT \->STR + " FOR " + ROT + " " +
SWAP \->STR + " EVAL * NEXT" + "1 " SWAP +
STR\-> THEN 1 ELSE 0 END ?f STOF { ERRN DOERR } IFT
\>> @ HP48[S/G][X][+]
\>> \>> '\PI' STO
Please be sure to provide the blank spaces in strings
exactly as shown; an extra space would not be harmful,
but a missing space might possibly be so.
The several ->STR commands are also important, in case any
argument happens to be a quoted variable name (as is allowed).
The odd-looking local variable names are unlikely
to conflict with any variable names in the stack arguments.
Now give this a try, using the same set-up as before.
BTW, lines with "@ HP48..." are necessary on 48 series,
but optional on 49/50 series, because the latter
do *not* check for "Invalid user functions"!
Had you noticed that the "piggyback" approach
(tricking the Summation function into doing products)
was a bit slowed down by the CAS on the 49/50 series,
while not too bad at all on the original 48 series?
The stand-alone "plotting engine" approach gives a speed
more like what you'd expect, but the last method,
as strange as it may appear, is the fastest of all
on a very long calculation, using only
a temporary variable during evaluation
(if the CAS happens to do this as well for summations,
any advantage is somewhat neutralized by other effects
when we introduce our own "helper" function).
If the result of a Product function would end up symbolic,
however, either because of symbolic start/end limits
or a "summand" which remains symbolic when evaluated,
you won't be happy with the programs introduced today
(they will either error or produce a humongously large
algebraic result), and you might have been better off
with the "piggyback" approach in such cases (also you
wouldn't have had to awkwardly quote your arguments).
So what's the "best way"?
I think it's good to have a broad choice of tools
in your kit, because even tools that ostensibly
"do the same thing" may each excel in turn
under different circumstances, just as do people.
With best wishes (and Happy Thanksgiving day in USA)