factorial

[fuzzy wozzy]

(define factorial

     N -> (factorial-h (* N (- N 1) (/ N 2)) (* (/ N 2) (/ N 2)) (- (/ N 2) 2) 1)

                                                               where (integer? (/ N 2))

     N -> (factorial-h (* N (/ (+ N 1) 2)) (* (/ (+ N 1) 2) (/ (+ N 1) 2))

                                                    (- (/ (+ N 1) 2) 2) 1))

(define factorial-h

     X S 0 A -> (* X A)

    X S C A -> (factorial-h X S (- C 1) (* (- S (* C C)) A) ))

[Robert Koeninger]

Interesting - is this your formulation? Does it have any advantages over just multiplying the numbers 1..N?

Also, it recurses infinitely for the inputs 0, 1 and 2 for me on Shen/SBCL 2.2. I tried switching the where condition to (= 0 (shen.mod N 2)), but that didn't help.

I tried refactoring it a bit to see if I could get some insight into how it works, pen and paper would probably be more useful in doing that.

(define factorial-h

  X S 0 A ->

    (* X A)

  X S C A ->

    (factorial-h

      X

      S

      (- C 1)

      (* (- S (* C C)) A)))

(define factorial

  0 -> 1

  1 -> 1

  2 -> 2

  N ->

    (let R (shen.mod N 2)

         L (if (= 0 R) (- N 1) 1)

         M (/ (+ N R) 2)

      (factorial-h

        (* N L M)

        (* M M)

        (- M 2)

        1)))

[fuzzy wozzy]

you're right, 0 -> 1, 1 -> 1, 2 -> 2 would have to be there to cover all cases,
I was thinking of a different way to calculate factorials of large numbers,
this method is not faster, but it gives a way of representing factorials 
in an equation form.

 

to find 10!

z = (/ 10 2) = 5

 

10! = 10 x 9 x 8 x 7 x 6 x 5 x 4 x 3 x 2 x 1
    = (10 x 1) x (9 x 2) x (8 x 3) x (7 x 4) x (6 x 5)
    = 10 x (z + 4) x (z - 3) x (z + 3) x (z - 2) x (z + 2) x (z - 1) x (z + 1) x z
    = 10 x 9 x 5 x (z**2 - 3**2) x (z**2 - 2**2) x (z**2 - 1**2)

 

except for 10, 9, and 5, the rest of calculation is a recursion of (z**2 - c**2),

n! = n x (n - 1) x (/ n 2) x II(z**2 - c*2),

if II is a symbol for a product of sequence.

[Antti Ylikoski]

How much does your method help in the practice?

It is an old wisdom in the Computer Science, that sophisticated coding to save some milliseconds is not worthwhile.

For example, with the SBCL bignums, if you calculate 2000!, or so, how much will this one help?

AJY

The EU

[Mark Tarver]

The time to work out 2000! using the simple factorial function in the maths-lib is 0.0 secs. The number has 5736 digits.

On Tuesday, 23 January 2018 03:40:53 UTC, fuzzy wozzy wrote:

(define factoria

[Antti Ylikoski]

Very very good, but can you find a (much) larger integer, such that the simple factorial calculation, and fuzzy wozzy's function, will produce significantly different computation times?   10000!  or 100000! or even more?

AJY

The EU

--
You received this message because you are subscribed to a topic in the Google Groups "Shen" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/qilang/0YVAmm2QW6Q/unsubscribe.
To unsubscribe from this group and all its topics, send an email to qilang+unsubscribe@googlegroups.com.
To post to this group, send email to qil...@googlegroups.com.
Visit this group at https://groups.google.com/group/qilang.
For more options, visit https://groups.google.com/d/optout.

[Mark Tarver]

The time to work out 2000! using the simple factorial function in the maths-lib is 0.0 secs. The number has 5736 digits.

 

The above is not my post btw - but Willi's - somehow he has wandered into my account!

Mark

[Mark Tarver]

I hope it is me (Willi) replying, and not Mark!

I have tested fuzzy's code, and it does not perform significantly better for up to 10000!

I notice though that in his code there are lots of divisions by 2. They are most likely integer divisions, such as 13/2 = 6.

Shen does not provide integer division. In the maths library,there is a function rsh which simulates a right shift, i.e. integer  division by 2.

In Shen/SBCL  rsh can be replaced by the equivalent native function. I suspect that fuzzy's version will then be faster. 

I have not enough time now to translate fuzzy's code in terms of rhs.

Willi

.

On Tuesday, 23 January 2018 03:40:53 UTC, fuzzy wozzy wrote:

[fuzzy wozzy]

in the example of 10!, since 9 = 5 + 4 and 1 = 5 - 4, these can go in the recursion, simplifying the equation,

if n is an even number greater than 2, and TT (two capital T's) is a symbol pi for product of sequence,

n! = ((n^2) / 2) x TT((n / 2)^2 - c^2), where c iterates from 1 to (n / 2) - 1

if n is an odd number greater than or equal to 3,

n! = ((n + 1) / 2) TT((n + 1) / 2)^2 - c^2), where c iterates from 1 to ((n + 1) / 2 - 1)

my computer goes into stack overflow at around 4000!, are there any significant difference after 10000!?

shen is versatile and flexible that anyone from beginner to expert could use any method to calculate 2000!

but what higher number factorial would be high enough that it wouldn't be trivial?

I think evne 1000000! might be trivial to some?

[fuzzy wozzy]

this is just a hunch, but now I believe that 10^100 can be done in shen, easily and within timely manner,

provided that the limitation of shen (because there eventually is, usually) can tolerate it...

but 10^1000000... 

[fuzzy wozzy]

the integer division is for dividing an even number only,

there is a function posted here, that can divide a number with millions of digits by two,

I forgot to add that, for an odd number, one can either do (n + 1) / 2, or (n - 1) / 2,

someone wrote a paper on how to find 20000! in lisp, several pages long, the answer has 77000+ digits,

is that a high number factorial? there is a number for 1000000!, and even 10^100!, so why only ask for 100000!?

[Antti Ylikoski]

Sorry to say, but Willi sometimes has been known to have made some problems.

AJY

The Eu

[fuzzy wozzy]

by limitation is meant stack overflow, which every language has, every computer too...

the recursion method shown here doesn't work with previously known factorial, such as 5 x 4 x 3!

but the equation or formula can show the enormity of the scope of solving something like (10^100)! rather nicely,

(10^100)! = ((10^100)^2) / 2 x TT( ((10^100) / 2)^2 - c^2), where c iterates from 1 to (10^100) / 2 - 1

right away, you can tell that this scale of calculation belongs to the computing giants, with easy access to

supercomputers and such, you wonder if you should stick with calculating 10000! or 100000! instead,

even those problems aren't all that trivial using off the shelf machines, but also there's the benefit of 100000!

being considered as trivial once you solve (10^100)! or something close to it, but an ordinary method of doing

n x (n - 1) x (n - 2)... 3 x 2 x 1 won't be enough even thought that's what it will take ultimately...

[fuzzy wozzy]

I appreciate willi's comment, always good to hear his input...

is rsh good for multiplying large numbers with millions of digits?

that is very hard to do, even if it could be done it would take days, possibly,

to multiply two numbers of millions of digits, and imagine if you had to do it

almost (5 x 10^198) times... it's amazing that technology has developed to a point

where something like that is a possibility and a reality...

no one much mentions about how fast some factorial can be calculated,

people seem more interested in how large of a number can be used...

if anyone could share some fast multiplication method of large numbers that

they know of, it would be nice...

[Antti Ylikoski]

https://en.wikipedia.org/wiki/Multiplication_algorithm#Karatsuba_multiplication

AJY

[fuzzy wozzy]

thank you... is that the method used by professional researchers who calculate (10^100)! and beyond?

it would be interesting to know how they do it...

the recursion method shows patterns that weren't obvious with n x (n - 1)... method,

such patterns might contribute 0.1% toward the solution, or even 0.001%...

something some shenturians might be already using, even with the recursion, the number is so large

it's kind of hard to imagine what complexity level of calculation would be entailed...

how many digits is the answer for (10^100)!?

if nothing else, this would be a good learning opportunity, possibly

[Antti Ylikoski]

"Professional researchers who calculate (10^100)! and beyond?"

I don't know what to say.

But, there exist methods to calculate (10^100)! type numbers.  Use logarithms.

AJY

[Antti Ylikoski]

One of my mathematics teachers said that for a mathematician, any finite number is small.

AJY

torstai 1. helmikuuta 2018 12.11.08 UTC+2 fuzzy wozzy kirjoitti:

[fuzzy wozzy]

for an even number n>5, if (n-2) is evenly divided by 4, this reduces the number of iteration 

to 1/4 of what would be required by the multiplication method, and it can use previously

known factorials...

(define fs

        N -> (fs-h (/ (* N N) 2) N (* (/ N 2) (/ N 2)) (/ (- N 2) 4) 1))

define fs-h

        X N S 0 A -> (* X A)

        X N S C A -> (fs-h X N S (- C 1)

                                                    (* C (- N C)

[Antti Ylikoski]

Textbook algorithms almost invariably are the best ones.

AJY

The EU

[fuzzy wozzy]

well, if it's true that one pattern begets another that begets another, ad infinitum...

with each pattern reducing the number of iteration, then is it any longer appropriate

to say that calculating (10^100)! is a hard problem... there is the actual task of

multiplying millions and billions of numbers, but I think there may be shenturians

who know the answer, but not willing to say it just yet...?

[fuzzy wozzy]

I don't think the many patterns that can be derived mean much, maybe they don't do much to reduce the number of iterations,

not any more than the two patterns that were posted above, which may be more than enough to recursively reduce the iteration

for ever smaller portions of calculation, the two patterns, c x (n - c) and (n/2)^2 - c^2, complement each other by iterating in the

reverse way to each other, which is why they can be used together in the same number of iteration, for example if c iterates from 1 to 20...

c x (n - c)     = 1 x (n - 1), 2 x (n - 2) , ..., 19 x (n - 19), 20 x (n - 20)

(n/2)^2 - c^2 = 20 x (n - 20), 19 x (n - 19), ... 3 x (n - 3), 2 x (n - 2), 1 x (n - 1) , in its equivalent form,

there may be other patterns that work more efficiently, but I'm not familiar with them, I don't know much about this, but I believe

that these two patterns can be used recursively to reduce the number of iteration more reliably than anything I've found so far...

  

[Antti Ylikoski]

Can you try that one in the practice?  Does it help significantly?

AJY

The EU

--

[fuzzy wozzy]

patterns need not be derived, but we can let the computer find it and use it directly within the calculation,

for example, instead of the pattern, c x (n - c), in the (fs) function, you can use the upper bound of c in the

original pattern, (n/2) - 1, and change it to (n/2) - c and use that in place of (c x (n - c)), such that

n! = (n^2)/2 x TT(((n/2)^2 - c^2) ((n/2)^2 - ((n/2) - c)^2)), c iterates from 1 to (n-2)/4

[Antti Ylikoski]

And how actually useful is this one?

Can you build an algorithm out of it, or, maybe a theorem?  

AJY

The EU

[fuzzy wozzy]

just trying out this or that...

since (10^100)! or some such high number factorials require a for-loop from 1 to (5 x 10^99) - 1, at least...

shen can do it easily, using list, of course...  

I feel kind of stuck at the moment, but I wasn't expecting to reduce the number of iteration to a quarter of

the usual method, so it's kinda fun to figure out different ways about it...

and if all I get is to reduce the for-loop, or its equivalent recursive loop, to a quarter of (5 x 10^99) - 1, I'm ok with that too...

[Antti Ylikoski]

"Reduce the for-loop, or its equivalent recursive loop, to a quarter of (5 x 10^99) - 1"?

Even a QUARTER of that number is pretty large.   A Shen FOR loop up to that number would take many nanoseconds.

Dr A. J. Y.

The EU

--

[fuzzy wozzy]

I'm not sure about significantly...

but my old machine goes into stack overflow at around 4500! with textbook method,

but with the recursion method here it can do 5500! with less halting, at least once...

all these patterns are new to me, I was surprised when one time willi mentioned that

he uses a different method to find factorials, I thought for sure that the textbook approach

was the ultimate and the most sensible one to use, but maybe this is old news to you...

then, let me ask you the same questions you are asking...

is there an algorithm or a theorem to be had from all this? or is this all an exercise in futility,

as the old computing wisdom once said, as you assert, to find fancy ways to save mere milli seconds?

somehow you might know...

[Antti Ylikoski]

See

https://en.wikipedia.org/wiki/Factorial#Computation

What is Willi's method?  He is a professional mathematician.

AJY

The EU

--

[fuzzy wozzy]

it can now be reduced to 1/8, and to 1/k where k approaches infinity 

[Antti Ylikoski]

What can be reduced, and how?

AJY

The EU

[fuzzy wozzy]

pretty amazing that ((n/2)^2 - c^2) may be all you need for n! where n gets to infinity 

[Antti Ylikoski]

Yes, but can you make a working divide and conquer algorithm out of this one?

Can you exemplaris gratis code the algorithm in Shen?

AJY

The EU

[fuzzy wozzy]

some say there's nothing new under the sun...

all this may be old news to some, perhaps to you too,

in which case all this would just be a rhetorical question...

an algorithm that wasn't workable and wasn't working would be an oxymoron or such,

now the (10^100)! doesn't seem so undividable and unsolvable...

I believe anyone who's been reading this short thread can write their own shen code easily,

since they would have 99.999% of the information needed to do so...

besides, the whole bulk of the program would just be a line or two long, mainly consisting of

the pattern, ((n/2)^2 - c^2), with shen finding millions and billions of its derivatives automatically...

 

[Antti Ylikoski]

The formula ((n/2)^2 - c^2)does not compute the factorial n!.

AJY

The EU

[fuzzy wozzy]

((n/2)^2 - c^2) is more a pattern, of which there are as many variations as there are numbers... 

if n is an even number greater than 2, and TT (two capital T's) is a symbol pi for product of sequence,

n! = ((n^2) / 2) x TT((n / 2)^2 - c^2), where c iterates from 1 to (n / 2) - 1

to reduce the number of iteration to 1/8 of textbook method, one could change the upper bound in (fs) to (/ (- N 2) 8)

add the following in (fs-h)

(/ (* (+ N (* 2 (- (* 2 C)  1))) (- N (* 2 (- (* 2 C) 1))) (+ (* 3 N) (* 2 (- (* 2 C) 1))) (- (* 3 N) (* 2 (- (* 2 C) 1))) ) 256)

finding actual patterns are ineffective because they do not show how the upper bound moves...

[Antti Ylikoski]

You write:

"n! = ((n^2) / 2) x TT((n / 2)^2 - c^2), where c iterates from 1 to (n / 2) - 1"

Can you prove this theorem?

AJY

The EU

[fuzzy wozzy]

patterns are nice because, once you find them, they're ready to use with no frills and they show that the method works as expected,

automating the process would be nice too, but they're not easy but just as complicated and complex, if not more so,

as easy as it is to find and use the patterns, they are powerful in their own right, for example it would take 320+ patterns

to reduce the iteration for (10^100)! to under 50, that seems easier to me than the alternative, but of course for factorial

of large numbers it would be less practical to rely on patterns alone, and automation of the whole process would have to be

utilized out of necessity rather than luxury...

[Antti Ylikoski]

How do you compute a factorial with "patterns"?

Can you make a theorem out of this one, and prove it?

AJY

The EU

[fuzzy wozzy]

 

I wouldn't say that automating the process is complicated or complex so much as it's convoluted,

more so than I would have expected it to be...

for a minute there I thought there might be some kind of meta pattern to be found in the patterns,

even with some sort of symmetry in the mix, which might be possible since we're working with

well ordered numbers rather than random ones?

but if I had a supercomputer I would go with the textbook method, easy and reliable...

[Mark Tarver]

Think 39 posts on the factorial function is probably enough for the moment.

Mark

--
You received this message because you are subscribed to the Google Groups "Shen" group.
To unsubscribe from this group and stop receiving emails from it, send an email to qilang+unsubscribe@googlegroups.com.

[fuzzy wozzy]

I could probably come up with something... you could too, right? 

[fuzzy wozzy]

Now the automation process for the reduction of iteration for n! works without a glitch.

This is a code segment for n! that reduces iteration to (n-2)/64.

It's obvious that (d) would be defined as, d = (n - (2 * (2 * c - 1))) / 4,

and defining (e) would be just as simple - can you define it?

(define fs-h

   X N S 0 A -> (* X A)
   X N S C A -> (fs-h X N S (- C 1)

                    (* (- S (* C C)) C (- N C)
                       (- S (* (d N C) (d N C))) (d N C) (- N (d N C))

                       (- S (* (e N C) (e N C))) (e N C) (- N (e N C))
                       (- S (* (d N (e N C)) (d N (e N C))))
                       (d N (e N C)) (- N (d N (e N C)))

                       (- S (* (f N C) (f N C))) (f N C) (- N (f N C))
                       (- S (* (d N (f N C)) (d N (f N C)))) (d N (f N C))
                       (- N (d N (f N C)))

                       (- S (* (e N (f N C)) (e N (f N C)))) (e N (f N C))
                       (- N (e N (f N C)))
 
                       (- S (* (d N (e N (f N C))) (d N (e N (f N C)))))
                       (d N (e N (f N C))) (- N (d N (e N (f N C))))

                       (- S (* (g N C) (g N C))) (g N C) (- N (g N C))
                       (- S (* (d N (g N C)) (d N (g N C)))) (d N (g N C))
                       (- N (d N (g N C)))

                       (- S (* (e N (g N C)) (e N (g N C))))
                       (e N (g N C)) (- N (e N (g N C)))
 
                       (- S (* (d N (e N (g N C))) (d N (e N (g N C)))))
                       (d N (e N (g N C))) (- N (d N (e N (g N C))))  

                       (- S (* (f N (g N C)) (f N (g N C)))) (f N (g N C))
                       (- N (f N (g N C)))
 
                      (- S (* (d N (f N (g N C))) (d N (f N (g N C)))))
                      (d N (f N (g N C))) (- N (d N (f N (g N C))))

                       (- S (* (e N (f N (g N C))) (e N (f N (g N C)))))
                       (e N (f N (g N C))) (- N (e N (f N (g N C))))

                       (- S (* (d N (e N (f N (g N C)))) (d N (e N (f N (g N C)))) ))
                       (d N (e N (f N (g N C)))) (- N (d N (e N (f N (g N C)))))          
                       A)))

[Mark Tarver]

Do you have any timings on this wrt the standard def of factorial?  

Mark 

[fuzzy wozzy]

 \* (fa) can find a portion of the factorial such that,
   (fa 11 21) = 11 x 12 x... 19 x 20 x 21, which is useful for using a
   known factorial value such as, 21! = 21 x 20 x... x 12 x 11 x 10!.
   (fa) works with the reduction method. The comment lines in (fa) are
   for 1/2 reduction of iteration, no need for Diff or (/ (...) 4) lines.
   d = (Diff - 2*(2*C - 1))/4.
   Try (fa 11 52), (fa 11 53), (fa 1 403),... *\
 
(define fa
     0 0 -> 1
     0 1 -> 1
     0 UB -> (fa 1 UB)      
   LB UB -> [upper bound > lower bound, lower bound >= 0, upper bound > 0]
                 where (or (>= LB UB) (< LB 0) (<= UB 0))

   LB UB -> (fa-h (* LB UB (midp LB UB)) (* (midp LB UB) (midp LB UB))
            \\ (- (/ (- UB LB) 2) 1) 1)
                   (/ (- UB (+ LB 2)) 4) (- UB LB) 1)
                           where (integer? (/ (- UB LB) 2))

   LB UB -> (fa-h1 (* UB (midq LB UB)) (* (midq LB UB) (midq LB UB))
             \\ (/ (- UB (+ LB 1)) 2) 1))
                    (/ (- UB (+ LB 1)) 4) (- UB LB) 1)) 

(define fa-h
   X S 0 Diff A -> (* X A)
   X S C Diff A -> (fa-h X S (- C 1) Diff

    (* (- S (* C C))

                                   (- S (* (- (/ Diff 2) C) (- (/ Diff 2) C)))
     A)))

(define fa-h1
   X S 0 Diff A -> (* X A)
   X S C Diff A -> (fa-h1 X S (- C 1) Diff

                       (* (- S (* C C))

                       (- S (* (/ (- Diff (- (* 2 C) 1)) 2) (/ (- Diff (- (* 2 C) 1)) 2))) 
          A)))

(define midp
 LB UB -> (+ LB (/ (- UB LB) 2)) )

(define midq
 LB UB -> (+ LB (/ (- UB (+ LB 1)) 2)))

[fuzzy wozzy]

Normally, (fa 1 6000) would go into stack overflow in the repl.

(fa 1000000103 1000001286) works okay, but (fa 1000000103 10000001286) doesn't go into stack overflow,

instead the cursor blinks for a long time as if busy calculating, until the repl is closed.

I'm not sure why that is...