[{"id":247553,"date":"2026-08-06T08:23:43","date_gmt":"2026-08-06T13:23:43","guid":{"rendered":"https:\/\/www.johndcook.com\/blog\/?p=247553"},"modified":"2026-08-06T09:22:07","modified_gmt":"2026-08-06T14:22:07","slug":"log1000","status":"publish","type":"post","link":"https:\/\/www.johndcook.com\/blog\/2026\/08\/06\/log1000\/","title":{"rendered":"Calculating log(1000!)"},"content":{"rendered":"

The previous post<\/a> pointed out that the following code such as the following unexpectedly works.<\/p>\n

>>> from math import log, factorial\r\n>>> log(factorial(1000))\r\n5912.128178488163\r\n<\/pre>\n
If you don’t find this unexpected, note that if you replace math.log<\/code> with numpy.log<\/code> the code will fail [1]. Functions like natural logarithm operate on real numbers. Real numbers are represented as floating point numbers in programming languages, and 1000! factorial is too large to represent as a standard floating point number. (More on that here<\/a>.)<\/p>\n
In this post I’d like to look at how you might calculate log(1000!) with less capable software, and even without software.<\/p>\n
One approach would be to sum the logarithms of the numbers 1 through 1000. This will give essentially the same result as above, with a little difference in the last couple decimal places due to rounding error.<\/p>\n
If you have a way to calculate 1000! but not a way to cast it to a floating point number, you could do this manually.<\/p>\n
>>> s = str(factorial(1000))\r\n>>> s[:16]\r\n'4023872600770937'\r\n>>> len(s)\r\n2568\r\n<\/pre>\n
This tells us 1000! = 4.023872600770937 \u00d7 102567<\/sup>. Therefore<\/p>\n
log(1000!) = log(4.023872600770937) + 2567 log(10)<\/p>\n
which only requires working with numbers of modest size.<\/p>\n
Calculating by hand<\/h2>\n
Now suppose it’s 1964. You don’t have a computer, or even a calculator, but you do have a copy of the recently published Handbook of Mathematical Functions by Abramowitz and Stegun (A&S). You turn to Table 6.6 “Factorials for large arguments.” This has values of factorial for 100, 200, 300, \u2026, 1000, so you can simply look up your answer to 20 decimal places.<\/p>\n
That was too easy; I didn’t expect that to be there when I started writing this post. If you wanted to compute log(950!), for example, you’d have to work harder. You could find A&S equation 6.1.41 (Stirling’s series) which says<\/p>\n
\"\\begin{align*}<\/p>\n
So how would you use this formula to calculate log(1000!)? Since n<\/em>! = \u0393(n<\/em> + 1), you set z<\/em> = 1001.<\/p>\n
You’d need to decide how many terms you need to use. Assuming the error is on the order of the first term you leave out, you’d reason that you could probably stop with the 1\/12z<\/em> term because the next term is between 10\u221211<\/sup> and 10\u221212<\/sup>.<\/p>\n
You find Table 4.2 has natural logarithms, but not for 1001. You can look up log(1.001), however, and at the bottom of the same page is log(10) to 16 decimal places, and you can find log(10) to 24 decimal places in Table 1.1. So you calculate<\/p>\n
log(1001) = log(1.001 \u00d7 10\u00b3) = log(1.001) + 3 log(10).<\/p>\n
You can find log(2) and log(\u03c0) in Table 1.1, and average them to find \u00bd log(2\u03c0).<\/p>\n
Here’s Python code to simulate the hand calculations.<\/p>\n
log2     = 0.6931_47180_55994_53094_172321 # Table 1.1\r\nlog10    = 2.3025_85092_99404_56840_179915 # Table 1.1\r\nlogpi    = 1.1447_29885_84940_01741_43427  # Table 1.1\r\nlog1_001 = 0.00099_95003_330835            # Table 4.2\r\n\r\nz = 1001\r\nlogz = log1_001 + 3*log10\r\ns = (z - 0.5)*logz - z + (log2 + logpi)\/2 + 1\/(12*z)\r\n\r\nprint(s)\r\n<\/pre>\n
This result differs from the one at the top of the post only in the last decimal place.<\/p>\n
Related posts<\/h2>\n
Doing calculations with tables is not as simple as “just look it up.” It takes a bit of skill.<\/p>\n