Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss

Number of digits in factorial

73 views
Skip to first unread message

mukesh tiwari

unread,
Jan 22, 2010, 4:47:34 PM1/22/10
to
Hello everyone. I am trying to solve a programming problem (https://
www.spoj.pl/problems/LENGFACT/) which ask for counting the number of
digits in factorial N, N<=5*10^9. I used Stirling approximation and
gamma function approximation. Here is my python code but i am
interested in mathematical part. How to solve this problem.

from cmath import *
import math
g = 7
p = [0.99999999999980993, 676.5203681218851, -1259.1392167224028,
771.32342877765313, -176.61502916214059, 12.507343278686905,
-0.13857109526572012, 9.9843695780195716e-6,
1.5056327351493116e-7]

def gamma(z):
z -= 1
x = p[0]
for i in range(1, g+2):
x += p[i]/(z+i)
t = z + g + 0.5
return log10(2*pi)+(z+0.5)*log10(t)-t*log10(e)+log10(x);

if __name__ == '__main__':
num = int(raw_input())
for i in range(num):
inp = int (raw_input())
if(inp==0 or inp==1 or inp==2):
print 1;
continue;
print int (math.floor(gamma(inp+1).real+0.5));


~

Dan Cass

unread,
Jan 25, 2010, 10:39:27 AM1/25/10
to
> return

> n log10(2*pi)+(z+0.5)*log10(t)-t*log10(e)+log10(x);
>
>
>
> if __name__ == '__main__':
> num = int(raw_input())
> for i in range(num):
> inp = int (raw_input())
> if(inp==0 or inp==1 or inp==2):
> print 1;
> continue;
> print int
> print int (math.floor(gamma(inp+1).real+0.5));
>
>
> ~

In base 1, the factorial n! has n! digits.
[OK I realize there's no "base 1"...]

Chip Eastham

unread,
Jan 25, 2010, 11:05:02 AM1/25/10
to
On Jan 22, 4:47 pm, mukesh tiwari <mukeshtiwari.ii...@gmail.com>
wrote:

> Hello everyone. I am trying to solve a programming problem
> (https://www.spoj.pl/problems/LENGFACT/) which ask for

> counting the number of digits in factorial N, N<=5*10^9.
> I used Stirling approximation and gamma function
> approximation. Here is my python code but i am interested
> in mathematical part. How to solve this problem.

I think you are on the right track with Stirling's
approximation. There are some tweaks to get more
correct digits, but this is unlikely to affect the
number of digits.

Note that the number of digits in (a base ten representation
of) a positive integer x is floor(log10(x)) + 1.

Some remarks you may find useful (or overkill):

[Q: Factorial for large numbers -- Google Answers]
http://answers.google.com/answers/threadview/id/33709.html

[Q: Large Factorial -- Google Answers]
http://answers.google.com/answers/threadview/id/509662.html

The latter in particular will point you to more than
you ever wanted to know about Stirling's approximation
and its refinements.

regards, chip

David W. Cantrell

unread,
Jan 25, 2010, 11:25:47 AM1/25/10
to
Chip Eastham <hard...@gmail.com> wrote:
> On Jan 22, 4:47=A0pm, mukesh tiwari <mukeshtiwari.ii...@gmail.com>

> wrote:
> > Hello everyone. I am trying to solve a programming problem
> > (https://www.spoj.pl/problems/LENGFACT/) which ask for
> > counting the number of digits in factorial N, N<=3D5*10^9.

> > I used Stirling approximation and gamma function
> > approximation. Here is my python code but i am interested
> > in mathematical part. How to solve this problem.
>
> I think you are on the right track with Stirling's
> approximation. There are some tweaks to get more
> correct digits, but this is unlikely to affect the
> number of digits.
>
> Note that the number of digits in (a base ten representation
> of) a positive integer x is floor(log10(x)) + 1.
>
> Some remarks you may find useful (or overkill):
>
> [Q: Factorial for large numbers -- Google Answers]
> http://answers.google.com/answers/threadview/id/33709.html
>
> [Q: Large Factorial -- Google Answers]
> http://answers.google.com/answers/threadview/id/509662.html
>
> The latter in particular will point you to more than
> you ever wanted to know about Stirling's approximation
> and its refinements.

Also see the OEIS entry for the number of decimal digits in n! :
<http://www.research.att.com/~njas/sequences/A034886>,
in particular, the last item in the formula section. I wonder if
the formula given in that item is correct.

David

David W. Cantrell

unread,
Jan 25, 2010, 11:11:51 PM1/25/10
to
David W. Cantrell <DWCan...@sigmaxi.net> wrote:
> Also see the OEIS entry for the number of decimal digits in n! :
> <http://www.research.att.com/~njas/sequences/A034886>,
> in particular, the last item in the formula section. I wonder if
> the formula given in that item is correct.

Now that I've had some time to think about it, I must say that I doubt the
correctness of that formula:

floor((log(2 pi n)/2 + n (log(n) - 1))/log(10)) + 1 [1]

Indeed, I now have a heuristic argument that [1] should fail to give the
number of decimal digits in n! correctly for infinitely many n, with the
first failure occurring before n = 600000 with 50% probability and before
n = 6 * 10^11 with 100% probability. (Of course, the first n for which [1]
fails may be difficult to find.)

But there is a very simple modification to the formula which should,
according to a heuristic argument, make it valid, with high probability,
for all n >= 1:

floor((log(2 pi n)/2 + n (log(n) - 1) + 1/(12 n))/log(10)) + 1 [2]

David W. Cantrell

David W. Cantrell

unread,
Feb 6, 2010, 1:49:16 AM2/6/10
to
On Jan 26, 4:11 am, David W. Cantrell <DWCantr...@sigmaxi.net> wrote:
> David W.Cantrell<DWCantr...@sigmaxi.net> wrote:
>
> > Also see the OEIS entry for the number of decimaldigitsin n! :

> > <http://www.research.att.com/~njas/sequences/A034886>,
> > in particular, the last item in the formula section. I wonder if
> > the formula given in that item is correct.

The author of [1] recently informed me that, as indicated in the OEIS
entry, it is merely based on Stirling's approximation, without proof
that it is correct for all n >= 2. But see below...

> Now that I've had some time to think about it, I must say that I
> doubt the correctness of that formula:
>
> floor((log(2 pi n)/2 + n (log(n) - 1))/log(10)) + 1 [1]
>
> Indeed, I now have a heuristic argument that [1] should fail to give

> the number of decimaldigitsin n! correctly for infinitely many n,


> with the first failure occurring before n = 600000 with 50%
> probability and before n = 6 * 10^11 with 100% probability. (Of
> course, the first n for which [1] fails may be difficult to find.)

My heuristic argument is very simple. Note that the number of decimal
digits in n! is (obviously) given by
floor(logGamma(n + 1)/log(10)) + 1. The difference between
logGamma(n + 1)/log(10) and (log(2 pi n)/2 + n (log(n) - 1))/log(10)
is approximately 1/(12 log(10) n), the infinite sum of which goes to
infinity, albeit very slowly. (Harmonic series.) That is, the sum of
the lengths of the intervals
[(log(2 pi n)/2 + n (log(n) - 1))/log(10), logGamma(n + 1)/log(10)]
is infinite. If any such interval happens to contain an integer, then
floor of (log(2 pi n)/2 + n (log(n) - 1))/log(10) will be smaller by 1
than floor of logGamma(n + 1)/log(10), and so cause formula [1] to
fail.

But I have not found any n for which [1] fails. And so I am now led
to wonder if there is something happening which is not captured in the
heuristic argument -- number theoretically, maybe -- which causes
[1] to give the correct number of decial digits in n! for all n >= 2.
I would appreciate other people's thoughts on this matter.

> But there is a very simple modification to the formula which should,
> according to a heuristic argument, make it valid, with high
> probability, for all n >= 1:
>
> floor((log(2 pi n)/2 + n (log(n) - 1) + 1/(12 n))/log(10)) + 1 [2]

Using the same sort of heuristic reasoning:
The difference between
(log(2 pi n)/2 + n (log(n) - 1) + 1/(12 n))/log(10) and
logGamma(n + 1)/log(10) is approximately 1/(360 log(10) n^3).
The sum of the lengths of the intervals
[logGamma(n + 1)/log(10), (log(2 pi n)/2 + n (log(n) - 1) + 1/(12
n))/log(10)], from n = 1..oo, is zeta(3)/(360 log(10)) = 0.00145...,
making it unlikely that any such interval would happen to contain an
integer.

Investigations in other bases:

Of course, my formula [2] can be modified to give the number of base-b
digits in n! by simply replacing log(10) with log(b). It seems that,
so modified, [2] works correctly for all n >= 1 and for all bases b.

But, when similarly modified, formula [1] fails in some bases. (It
fails in all bases for n = 1, but that is of little interest.
Henceforth, we consider just n > 1.) It will fail at n in base b if
b = n! For example, it fails in base 6 for n = 3 because b = 6 = 3!
and it fails in base 2 for n = 2 because b = 2 = 2! Here are the
other failures I've found in various bases:

b n
5 47
6 2751
12 516
14 29
14 4799
14 11321

In all those cases, not surprisingly, n! in base b consists of 1
followed by some 0's. For example, in base 14, we have 29! =
1006bb5d6db12a825d9398ac0000

The fact that formula [1], when modified for base b, fails in some
bases strengthens my skepticism about its working in base 10 for
all n >= 2.

David W. Cantrell

David W. Cantrell

unread,
Feb 6, 2010, 1:56:19 AM2/6/10
to
"David W. Cantrell" <DWCan...@sigmaxi.net> wrote:

Oops. To avoid confusion I should have specified that my values of n are
given in _decimal_ notation.

> In all those cases, not surprisingly, n! in base b consists of 1
> followed by some 0's. For example, in base 14, we have 29! =
> 1006bb5d6db12a825d9398ac0000

That is, (29_10)! = 1006bb5d6db12a825d9398ac0000_14.

DWC

Richard Tobin

unread,
Feb 6, 2010, 11:05:05 AM2/6/10
to
In article <19ce633c-f872-4f7c...@z26g2000yqm.googlegroups.com>,

David W. Cantrell <DWCan...@sigmaxi.net> wrote:

>The fact that formula [1], when modified for base b, fails in some
>bases strengthens my skepticism about its working in base 10 for
>all n >= 2.

Could you determine whether it fails for 3121515!. I don't have the
facilities to hand to check it accurately.

-- Richard
--
Please remember to mention me / in tapes you leave behind.

David W. Cantrell

unread,
Feb 6, 2010, 11:55:25 AM2/6/10
to
ric...@cogsci.ed.ac.uk (Richard Tobin) wrote:
> In article
> <19ce633c-f872-4f7c...@z26g2000yqm.googlegroups.com>,
> David W. Cantrell <DWCan...@sigmaxi.net> wrote:
>
> >The fact that formula [1], when modified for base b, fails in some
> >bases strengthens my skepticism about its working in base 10 for
> >all n >= 2.
>
> Could you determine whether it fails for 3121515!. I don't have the
> facilities to hand to check it accurately.

Thanks for your suggestion.

3121515! = 1.000000695... * 10^18916606

and so it has 18916607 decimal digits. The formula given in the OEIS entry
works correctly in this case. For n = 3121515, the interval

[(log(2 pi n)/2 + n (log(n) - 1))/log(10), logGamma(n + 1)/log(10)] =

[18916606.000000290..., 18916606.000000301...]

and so does not contain an integer.

Any other suggestions?

David

Chip Eastham

unread,
Feb 6, 2010, 12:11:58 PM2/6/10
to
On Feb 6, 1:56 am, David W. Cantrell <DWCantr...@sigmaxi.net> wrote:

Hi, David:

I suspect your heuristic arguments are correct,
and we are simply frustrated by the horribly
slow divergence of the harmonic series.

Have the intervals of disagreement between [1]
and [2] been exhaustively checked up to some
point, e.g. n = 600000 as you supposed would
give a 50% chance of counterexample? It seems
this would be easy to program, and if none is
found, the data might suggest some systematic
reason for failure.

regards, chip

Richard Tobin

unread,
Feb 6, 2010, 12:25:23 PM2/6/10
to
In article <20100206120206.850$l...@newsreader.com>,

David W. Cantrell <DWCan...@sigmaxi.net> wrote:

>> Could you determine whether it fails for 3121515!. I don't have the
>> facilities to hand to check it accurately.

That was just the first case where 64-bit IEEE floating point numbers
give different results, calculating the factorial by adding the logs.
Of course, I couldn't tell whether it was the accumulated floating
point error or a real failure.

>Any other suggestions?

Using "long double", the first candidates are 48655817 and 91948153
(beyond that they occur often, suggesting that the floating point
error has become significant).

David W. Cantrell

unread,
Feb 6, 2010, 12:33:38 PM2/6/10
to
Chip Eastham <hard...@gmail.com> wrote:

> On Feb 6, 1:56=A0am, David W. Cantrell <DWCantr...@sigmaxi.net> wrote:
> > "David W. Cantrell" <DWCantr...@sigmaxi.net> wrote:
> >
> > > On Jan 26, 4:11 am, David W. Cantrell <DWCantr...@sigmaxi.net> wrote:
> > > > David W.Cantrell<DWCantr...@sigmaxi.net> wrote:
> >
> > > > > Also see the OEIS entry for the number of decimal digits in n!:

> > > > > <http://www.research.att.com/~njas/sequences/A034886>,
> > > > > in particular, the last item in the formula section. I wonder if
> > > > > the formula given in that item is correct.
> >
> > > The author of [1] recently informed me that, as indicated in the OEIS
> > > entry, it is merely based on Stirling's approximation, without proof
> > > that it is correct for all n >= 2. But see below...
> >
> > > > Now that I've had some time to think about it, I must say that I
> > > > doubt the correctness of that formula:
> >
> > > > floor((log(2 pi n)/2 + n (log(n) - 1))/log(10)) + 1 [1]
> >
> > > > Indeed, I now have a heuristic argument that [1] should fail to
> > > > give the number of decimal digits in n! correctly for infinitely

> > > > many n, with the first failure occurring before n = 600000 with
> > > > 50% probability and before n = 6 * 10^11 with 100% probability.
> > > > (Of course, the first n for which [1] fails may be difficult to
> > > > find.)
> >
> > > My heuristic argument is very simple. Note that the number of decimal
> > > digits in n! is (obviously) given by
> > > floor(logGamma(n + 1)/log(10)) + 1. The difference between
> > > logGamma(n + 1)/log(10) and (log(2 pi n)/2 + n (log(n) - 1))/log(10)
> > > is approximately 1/(12 log(10) n), the infinite sum of which goes to
> > > infinity, albeit very slowly. (Harmonic series.) That is, the sum of
> > > the lengths of the intervals
> > > [(log(2 pi n)/2 + n (log(n) - 1))/log(10), logGamma(n + 1)/log(10)]
> > > is infinite. If any such interval happens to contain an integer, then
> > > floor of (log(2 pi n)/2 + n (log(n) - 1))/log(10) will be smaller by
> > > 1 than floor of logGamma(n + 1)/log(10), and so cause formula [1] to
> > > fail.
> >
> > > But I have not found any n for which [1] fails. And so I am now led
> > > to wonder if there is something happening which is not captured in
> > > the heuristic argument -- number theoretically, maybe -- which causes
> > > [1] to give the correct number of decial digits in n! for all n >= 2

Thanks, Chip. That was my initial thinking too.

> Have the intervals of disagreement between [1]
> and [2] been exhaustively checked up to some
> point, e.g. n = 600000 as you supposed would
> give a 50% chance of counterexample? It seems
> this would be easy to program, and if none is
> found, the data might suggest some systematic
> reason for failure.

Yes, it's easy to program. Using Mathematica, apparently [1] works
correctly up to at least 10^7.

[And BTW, for the nondecimal bases, my investigations were exhaustive
up to n = 10^6 for binary through hexadecimal.]

Regards,
David

David W. Cantrell

unread,
Feb 6, 2010, 12:45:01 PM2/6/10
to
ric...@cogsci.ed.ac.uk (Richard Tobin) wrote:
> In article <20100206120206.850$l...@newsreader.com>,
> David W. Cantrell <DWCan...@sigmaxi.net> wrote:
>
> >> Could you determine whether it fails for 3121515!. I don't have the
> >> facilities to hand to check it accurately.
>
> That was just the first case where 64-bit IEEE floating point numbers
> give different results, calculating the factorial by adding the logs.
> Of course, I couldn't tell whether it was the accumulated floating
> point error or a real failure.
>
> >Any other suggestions?
>
> Using "long double", the first candidates are 48655817 and 91948153
> (beyond that they occur often, suggesting that the floating point
> error has become significant).

Thanks again. The formula in OEIS also works correctly for both of those
candidates.

David

David W. Cantrell

unread,
Feb 6, 2010, 4:52:51 PM2/6/10
to
David W. Cantrell <DWCan...@sigmaxi.net> wrote:

Oops, again. There was an important omission from my table of failures.

In some nondecimal bases, I searched beyond n = 10^6. In particular, in
binary, I had searched up to n = 2 * 10^8, but I had somehow overlooked a
failure which Mathematica had already found at n = 25463836:

Formula [1] {with log(10) replaced by log(2), of course} predicts that
25463836! should have 589723393 bits. But it actually has 589723394 bits.
The interval
[(log(2 pi n)/2 + n (log(n) - 1))/log(2), logGamma(n + 1)/log(2)] = [589723392.9999999961..., 589723393.00000000083...] contains an integer.
Expressing 25463836! in binary, between the leading 1 and the
next 1, there are about 30 0's, if my calculation is correct.

So, in binary, [1] fails for n = 1, 2 and 25463836. And I suppose
that there are other failures beyond n = 2 * 10^8.

David

Chip Eastham

unread,
Feb 8, 2010, 9:15:43 AM2/8/10
to

To summarize briefly, we seek integer n > 1 such that interval:

[ (log(2 pi n)/2 + n (log(n) - 1))/log(10),
(log(2 pi n)/2 + n (log(n) - 1) + 1/(12n))/log(10) ]

contains an integer. My question concerns the search method.
Are you iterating over n?

I'm curious whether defining:

f(n) = (log(2 pi n)/2 + n (log(n) - 1))/log(10)

and "solving" f(n) = M for increasing values of M will give
a more efficient search method. Am I off-base, or a day late
(and dollar short) with this idea?

regards, chip

David W. Cantrell

unread,
Feb 8, 2010, 11:34:32 AM2/8/10
to
Chip Eastham <hard...@gmail.com> wrote:
> To summarize briefly, we seek integer n > 1 such that interval:
>
> [ (log(2 pi n)/2 + n (log(n) - 1))/log(10),
> (log(2 pi n)/2 + n (log(n) - 1) + 1/(12n))/log(10) ]
>
> contains an integer. My question concerns the search method.
> Are you iterating over n?

Yes, regrettably. I have sought, in vain, a more efficient search method.

> I'm curious whether defining:
>
> f(n) = (log(2 pi n)/2 + n (log(n) - 1))/log(10)
>
> and "solving" f(n) = M for increasing values of M will give
> a more efficient search method. Am I off-base, or a day late
> (and dollar short) with this idea?

Your idea is equivalent, I think, to something I had thought about but then
rejected as being even less efficient than just iterating over n.

Picking an integer M and hoping that solving f(n) = M will have a
nearly-integer solution n is just about the same as solving n! = 10^M . The
basic idea in either case is that we want a factorial which is
approximately an integer power of the base.

Let me illustrate using an example which I had mentioned earlier:
base b = 12, n = 516. Of course, since I had already found that example,
I already knew that 516! is close to 12^1091, so I now use M = 1091:

In[79]:= FindRoot[(Log[2*Pi*n]/2 + n*(Log[n] - 1))/Log[12] == 1091,
{n, 516}, WorkingPrecision -> 12]
Out[79]= {n -> 516.000011528}


In[80]:= FindRoot[n! == 12^1091, {n, 516}, WorkingPrecision -> 12]
Out[80]= {n -> 515.999985676}

In originally finding my base-12 example, I had iterated over n up to 516.
To have found the same example by solving n! = 12^M (or f(n) = M), I would
have had to iterate over M up to 1091, and so it's even less efficient.

The crux of the matter is that n! eventually increases more rapidly than
does b^M. Thus, we "cover ground" more rapidly by iterating over n than by
iterating over M.

Regards,
David

David W. Cantrell

unread,
Feb 16, 2010, 8:08:45 AM2/16/10
to

For base-10, I can now say that [1] apparently works correctly up to
at least 10^9. I doubt I will search further for a base-10 failure.

> > [And BTW, for the nondecimal bases, my investigations were exhaustive
> > up to n = 10^6 for binary through hexadecimal.]
>
> Oops, again. There was an important omission from my table of failures.
>
> In some nondecimal bases, I searched beyond n = 10^6. In particular, in
> binary, I had searched up to n = 2 * 10^8, but I had somehow overlooked a
> failure which Mathematica had already found at n = 25463836:
>
> Formula [1] {with log(10) replaced by log(2), of course} predicts that
> 25463836! should have 589723393 bits. But it actually has 589723394 bits.
> The interval
> [(log(2 pi n)/2 + n (log(n) - 1))/log(2), logGamma(n + 1)/log(2)] =
> [589723392.9999999961..., 589723393.00000000083...] contains an integer.
> Expressing 25463836! in binary, between the leading 1 and the next 1,
> there are about 30 0's, if my calculation is correct.
>
> So, in binary, [1] fails for n = 1, 2 and 25463836. And I suppose
> that there are other failures beyond n = 2 * 10^8.

In my original programming, I simply looked for cases when [1] would fail.
But later, I decided it would also be interesting, in base-10, to find
cases when [1] is "close to failure". The closest near-failure (but it's
not really very close) is

252544447! = 1.0000000009... * 10^2012285104

that is, between its leading 1 and the next nonzero digit,
there are nine 0's.

I must still guess that [1] fails in base-10, and so I plan to add an
appropriate comment to the OEIS entry.

David

spudnik

unread,
Feb 16, 2010, 12:43:50 PM2/16/10
to
of course, there is a base-one;
what is it's digital counter, by induction on base-ten?

in factorial base, it has n digits; eh?

> In base 1, the factorial n! has n! digits.
> [OK I realize there's no "base 1"...]

thus:
sea-level is not rising, globally --
http://www.21stcenturysciencetech.com/Articles%202007/MornerInterview.pdf
-- and warming is mostly equatorial. however,
there is massive loss of soil, and that might change *relative* sea-
level,
in some locations, as well as dysplace some sea!

thus quoth:
Let’s take a look at the complexity of polar bear life. First, the
polar bear has been around for about 250,000 years, having survived
both an Ice Age, and the last Interglacial period (130,000 years ago),
when there was virtually no ice at the North Pole. Clearly, polar
bears have adapted to the changing environment, as evidenced by their
presence today.
(This fact alone makes the polar bear smarter than Al Gore and the
other global warming alarmists. Perhaps the polar bear survived the
last Interglacial because it did not have computer climate models that
said polar bears should not have survived!)
http://www.21stcenturysciencetech.com/Articles%202007/GW_polarbears.pdf
http://www.21stcenturysciencetech.com/Global_Warming.html


> It amuses me that the people who want "change" are so afraid of change--God
> forbid that the seal level rise two inches.

thus:
the photographic record that I saw,
in some rather eclectic compendium of Einsteinmania,
seemed to show quite a "bending" effect, I must say;
not that the usual interpretation is correct, though.

Nude Scientist said:
> > "Enter another piece of luck for Einstein. We now know that the light-
> > bending effect was actually too small for Eddington to have discerned

--Another Flower for Einstein:
http://www.21stcenturysciencetech.com/articles/spring01/Electrodynamics.html

--les OEuvres!
http://wlym.com

--Stop Cheeny, Rice & the ICC in Sudan;
no more Anglo-american quagmires!
http://larouchepub.com/pr/2010/100204rice

0 new messages