Post Reply 
HP 50g Double factorial
05-01-2019, 11:19 AM (This post was last modified: 05-01-2019 11:29 AM by joeres.)
Post: #1
HP 50g Double factorial
Hi,

currently I extend my 50g user's guide with additional information. Among others, I would like to add a User RPL program to calculate the double factorial (https://en.wikipedia.org/wiki/Double_factorial).
My version:
Code:

@program DF
@in: number > 0
@out: double factorial of number
\<<
  \-> N
  \<<
     IF 1 > THEN @valid number?
       N 2 MOD @even or odd?
       IF 0 == THEN @even number
          2 N 2 / ^ N 2 / ! *
       ELSE @odd number
          N 1 + ! 2 N 1 + 2 / ^ N 1 + 2 / ! * /          
       END
     ELSE @not valid
       1
     END 
  \>>
\>>
It based on the fourmulars n!!=2^(n/2)*(n/2)! for even numbers and n!!=(n+1)!/2^(n+1/2)*(n+1/2)! for odd numbers and it works, but can you tell me where you could optimize it or how you would solve this task?
But no additional libraries should be used and it should be a User RPL program.

Many thanks and kind regards
Joerg
Find all posts by this user
Quote this message in a reply
05-01-2019, 01:40 PM (This post was last modified: 05-01-2019 02:00 PM by Gilles.)
Post: #2
RE: HP 50g Double factorial
Hi Joerg , another way is :

Code:
 « -> n 'IFTE(n MOD 2,(n+1)!/(2^((n+1)/2)*((n+1)/2)!),2^(n/2)*(n/2)!)'»

returns 1 for zero, and an error for n<0
And works in symbolic mode with a symbolic variable
Find all posts by this user
Quote this message in a reply
05-01-2019, 03:57 PM
Post: #3
RE: HP 50g Double factorial
Hi Gilles,

(05-01-2019 01:40 PM)Gilles Wrote:  
Code:
 « -> n 'IFTE(n MOD 2,(n+1)!/(2^((n+1)/2)*((n+1)/2)!),2^(n/2)*(n/2)!)'»
returns 1 for zero, and an error for n<0
And works in symbolic mode with a symbolic variable

thanks for the answer with the short and sweet version.
I like that a question mark is displayed for negative numbers. However, at -1 and -2 not. This is certainly related to the formula?

Kind regards
Joerg
Find all posts by this user
Quote this message in a reply
05-01-2019, 04:14 PM (This post was last modified: 05-01-2019 05:52 PM by Gilles.)
Post: #4
RE: HP 50g Double factorial
(05-01-2019 03:57 PM)joeres Wrote:  thanks for the answer with the short and sweet version.
I like that a question mark is displayed for negative numbers. However, at -1 and -2 not. This is certainly related to the formula?

Kind regards
Joerg

Yes it is. We can do :

Code:
«  -> n 
 «
  IF 'n<0 OR FP(n)' THEN  #203h DOERR ELSE
   'IFTE(n MOD 2,(n+1)!/(2^((n+1)/2)*((n+1)/2)!),2^(n/2)*(n/2)!)' ->NUM
  END 
 »
»

A negative or not integer input will result an error "Bad argument value".
Find all posts by this user
Quote this message in a reply
05-01-2019, 04:59 PM (This post was last modified: 05-01-2019 04:59 PM by grsbanks.)
Post: #5
RE: HP 50g Double factorial
Gilles, I'm curious. What are these rectangles in your post supposed to be?

Ces petits rectangles sont censés représenter quoi au juste (voir image attachée) ?

   
Find all posts by this user
Quote this message in a reply
05-01-2019, 05:05 PM (This post was last modified: 05-01-2019 09:26 PM by Giuseppe Donnini.)
Post: #6
RE: HP 50g Double factorial
Here's another solution:

Considering that  n!!  is equivalent to  n·(n−2)·(n−4)· ... ·1  for every positive integer  n ≠ 0  (odd or even), all it needs is a simple loop which calculates the next factor and multiplies it with the running product.

1             @ Initialize result.
SWAP 2 FOR f  @ Loop down from n to 2.
  f *         @ Multiply with current factor.
-2 STEP       @ Decrement factor by 2 ; repeat until factor < 2.

Taking into account the special case where  n = 0 , we finally have:

Code:

@ NAME  : DF (Double factorial)
@ STACK : ( n --> n!! )
\<<
  1                 @ Initialize result.
  IF SWAP DUP 0 ==  @ If n = 0,
  THEN DROP         @   then n!! = 1 (just drop n, 1 already on stack).
  ELSE 2 FOR f      @ Otherwise: Loop down from n to 2.
    f *             @   Multiply with current factor.
  -2 STEP END       @   Decrement factor by 2 ; repeat until factor < 2.
\>>
Find all posts by this user
Quote this message in a reply
05-01-2019, 05:50 PM
Post: #7
RE: HP 50g Double factorial
Hi joeres, the program you posted contains a small bug:

After having created and initialized the local variable N, you need to put it on the stack again for comparison, because \-> consumes its argument.

So, instead of:

\<<
  \-> N
  \<<
     IF 1 > THEN

(...)

you have to do:

\<<
  \-> N
  \<<
     IF N 1 > THEN

(...)
Find all posts by this user
Quote this message in a reply
05-01-2019, 06:01 PM
Post: #8
RE: HP 50g Double factorial
(05-01-2019 04:14 PM)Gilles Wrote:  Yes it is. We can do :

Code:
«  -> n 
 «
  IF 'n<0 OR FP(n)' THEN  #203h DOERR ELSE
   'IFTE(n MOD 2,(n+1)!/(2^((n+1)/2)*((n+1)/2)!),2^(n/2)*(n/2)!)' ->NUM
  END 
 »
»

A negative or not integer input will result an error "Bad argument value".

Thanks Gilles, it works fine.

Kind regards
Joerg
Find all posts by this user
Quote this message in a reply
05-01-2019, 06:17 PM
Post: #9
RE: HP 50g Double factorial
(05-01-2019 05:05 PM)Giuseppe Donnini Wrote:  
1          @ Initialize result.
SWAP 2     @ Loop down from n to 2.
FOR f f *  @ Multiply with current factor.
-2 STEP    @ Decrement factor by 2 ; repeat until factor < 2.

Taking into account the special case where  n = 0 , we finally have ...

I am not familiar with RPL code, but is there a loop construct where you test the condition first ?
That way, the 0 special case need not be added.

BTW, I like your code. Doing odd n!! as n! / (n-1)!! might overflow the result.
Find all posts by this user
Quote this message in a reply
05-01-2019, 06:26 PM
Post: #10
RE: HP 50g Double factorial
(05-01-2019 05:05 PM)Giuseppe Donnini Wrote:  Here's another solution:

Considering that  n!!  is equivalent to  n·(n-2)·(n-4)· ... ·1  for every positive integer  n ≠ 0  (odd or even), all it needs is a simple loop which calculates the next factor and multiplies it with the running product. Such a solution could do entirely without local variables and algebraic objects, requiring only stack manipulations, and very few in number:

1          @ Initialize result.
SWAP 2     @ Loop down from n to 2.
FOR f f *  @ Multiply with current factor.
-2 STEP    @ Decrement factor by 2 ; repeat until factor < 2.

Taking into account the special case where  n = 0 , we finally have:

Code:

@ NAME  : DF (Double factorial)
@ STACK : ( n --> n!! )
\<<
  1                 @ Initialize result.
  IF SWAP DUP 0 ==  @ If n = 0,
  THEN DROP         @   then n!! = 1 (just drop n, 1 already on stack).
  ELSE 2            @ Otherwise: Loop down from n to 2.
  FOR f f *         @   Multiply with current factor.
  -2 STEP           @   Decrement factor by 2 ; repeat until factor < 2.
\>>

Hi Giuseppe,

Thank you, very interesting and nice solution. I tested it (END added). Works fine, but if I use very large numbers, the program runs in an infinite loop. The other programs output the message "Integer too large". But I can easily query this.

Kind regards
Joerg
Find all posts by this user
Quote this message in a reply
05-01-2019, 06:35 PM
Post: #11
RE: HP 50g Double factorial
(05-01-2019 05:50 PM)Giuseppe Donnini Wrote:  Hi joeres, the program you posted contains a small bug:

After having created and initialized the local variable N, you need to put it on the stack again for comparison, because \-> consumes its argument.

So, instead of:

\<<
  \-> N
  \<<
     IF 1 > THEN

(...)

you have to do:

\<<
  \-> N
  \<<
     IF N 1 > THEN

(...)

Many thanks for the hint. I have corrected it.

Kind regards
Joerg
Find all posts by this user
Quote this message in a reply
05-01-2019, 07:00 PM (This post was last modified: 05-01-2019 07:08 PM by John Keith.)
Post: #12
RE: HP 50g Double factorial
(05-01-2019 04:14 PM)Gilles Wrote:  A negative or not integer input will result an error "Bad argument value".

Not necessarily. Smile See this thread.

Here is a version for the HP-48 or HP 50:

Code:

\<< \pi EVAL DUP2 * COS 1. - .25 * SWAP .5 * SWAP ^ SWAP .5 * DUP 2. SWAP ^ SWAP 1. + GAMMA * *
\>>

It will handle non-integer arguments, and complex numbers if complex mode is enabled. Negative even integers return results near negative infinity.

See the graphs here.
Find all posts by this user
Quote this message in a reply
05-01-2019, 07:01 PM
Post: #13
RE: HP 50g Double factorial
(05-01-2019 06:26 PM)joeres Wrote:  If I use very large numbers, the program runs in an infinite loop.
I cannot confirm this behavior. If n!! is too large, I get MAXR, i.e. 9.99999999999E499.
You must be using ZINTs, then. As a matter of fact, I wrote the code with an HP-48 in mind, so it is primarily intended for floats.

(05-01-2019 06:26 PM)joeres Wrote:  The other programs output the message "Integer too large".
Again, this must be a ZINT-specific message.

P.S. Thanks for reminding me of the missing "END", I just added it.
Find all posts by this user
Quote this message in a reply
05-01-2019, 07:06 PM (This post was last modified: 05-01-2019 07:15 PM by Gilles.)
Post: #14
RE: HP 50g Double factorial
(05-01-2019 04:59 PM)grsbanks Wrote:  Gilles, I'm curious. What are these rectangles in your post supposed to be?
Ces petits rectangles sont censés représenter quoi au juste (voir image attachée) ?

I think it's because I use HPUserEdit and drag and drop from it and then, I change the special chars by hand :/

[Image: fske.png]
Find all posts by this user
Quote this message in a reply
05-01-2019, 09:26 PM (This post was last modified: 05-02-2019 02:55 PM by Giuseppe Donnini.)
Post: #15
RE: HP 50g Double factorial
(05-01-2019 06:17 PM)Albert Chan Wrote:  I am not familiar with RPL code, but is there a loop construct where you test the condition first ? That way, the 0 special case need not be added.

Thanks, Albert, for your valuable comments. Sure there is! Here we go:

\<<
  1 SWAP          @ Initialize result.
  WHILE
    DUP 2 \>=     @ Repeat while factor >= 2.
  REPEAT
    DUP ROT *     @ Multiply product with current factor.
    SWAP 2 -      @ Decrement factor by 2.
  END DROP        @ Drop remaining factor (by now invalid, i.e. < 2).
\>>

Saves 9 bytes, but negative arguments will return 1, too, whereas my former implementation left them unchanged.
Find all posts by this user
Quote this message in a reply
05-02-2019, 08:33 PM
Post: #16
RE: HP 50g Double factorial
(05-01-2019 07:00 PM)John Keith Wrote:  Not necessarily. Smile See this thread.

Here is a version for the HP-48 or HP 50:

Code:

\<< \pi EVAL DUP2 * COS 1. - .25 * SWAP .5 * SWAP ^ SWAP .5 * DUP 2. SWAP ^ SWAP 1. + GAMMA * *
\>>

It will handle non-integer arguments, and complex numbers if complex mode is enabled. Negative even integers return results near negative infinity.

See the graphs here.

Hi John,
Thank you for this information. I needed something to understand everything reasonably. The connection between the gamma function and the double faculty I did not know. Mathematics and RPL, two really incredible things ....

many Greetings
Joerg
Find all posts by this user
Quote this message in a reply
05-02-2019, 09:34 PM
Post: #17
RE: HP 50g Double factorial
(05-01-2019 07:01 PM)Giuseppe Donnini Wrote:  I cannot confirm this behavior. If n!! is too large, I get MAXR, i.e. 9.99999999999E499.
You must be using ZINTs, then. As a matter of fact, I wrote the code with an HP-48 in mind, so it is primarily intended for floats.

Again, this must be a ZINT-specific message.

Sorry Giuseppe, you are right. The answer will be given after about 21 seconds for number 9999. I was just confused, the other program versions end immediately.

Kind regards,
Joerg
Find all posts by this user
Quote this message in a reply
Post Reply 




User(s) browsing this thread: 1 Guest(s)