Counting Abelian groups of a given size

Question

Background

Last time, we counted groups of a given size, which is a non-trivial problem.

This time, we'll only count Abelian groups, i.e., groups with a commutative operation. Formally, a group (G, ∗) is Abelian if x ∗ y = y ∗ x for for all x, y in G.

The problem becomes much simpler this way, so we're going to count them efficiently.

Task

Write a program or function that accepts a non-negative integer n as input and prints or returns the number of non-isomorphic Abelian groups of order n.

One way of calculating the number of groups – which we'll denote by A(n) – is by observing the following:

• A(0) = 0

• If p is a prime, A(p^k) is equal to the number of integer partitions of k. (cfr. OEIS A000041)

• If n = mk, and m and k are co-prime, A(n) = A(m)A(k).

You may use this or any other method of calculating A(n).

Test cases

Input               Output
0                   0
1                   1
2                   1
3                   1
4                   2
5                   1
6                   1
7                   1
8                   3
9                   2
10                  1
11                  1
12                  2
13                  1
14                  1
15                  1
16                  5
17                  1
18                  2
19                  1
20                  2
4611686018427387904 1300156
5587736968198167552 155232
9223371994482243049 2

(taken from OEIS A000688)

Additional rules

• Given enough time, RAM and a register size that can hold the input, your code should work (in theory) for arbitrarily large integers.

• Your code must work for all integers between 0 and 2⁶³ - 1 and finish in under 10 minutes on my machine (Intel i7-3770, 16 GiB RAM, Fedora 21).

  Please make sure you time your code for the last three test cases before submitting your answer.

• Built-ins that trivialize this task, such as Mathematica's FiniteAbelianGroupCount, are not allowed.

• Built-ins that return or count the integer partitions of a number or the partitions of a list are not allowed.

• Standard code-golf rules apply.

Answer 1

CJam (39 38 bytes)

qimF{~M\{_ee{~\)<1b}%1+a\+}*0=1be&}%:*

Online demo

This follows the suggested line of finding a prime factorisation (mF) and then calculating the partitions of each power and taking their product.

There are two special cases for mF: it factors 0 as 0^1 and 1 as 1^1. The latter doesn't require special treatment: there is one Abelian group of size 1, and one partition of 1. However, zero does require a special case.

The partition counting uses a recurrence for A008284(n, k), the number of partitions of n into k parts. In OEIS it is given as

T(n, k) = Sum_{i=1..k} T(n-k, i), for 1<=k<=n-1; T(n, n) = 1 for n >= 1.

but I think it's more useful to think of the sum as ranging from 1 to min(k, n-k).

Dissection

q~              e# Parse input into an integer
mF              e# Factorise it
{               e# For each factor p^a
  ~             e#   Split the array [p a]
                e#   The following lines count partitions of a
                e#   (Note: they would be buggy if a were ever 0, but it isn't)
  M\{           e#   Starting with a table of zero rows, repeat a times
    _ee         e#     Copy table and pair each row with its index
    {~\)<1b}%   e#     Extract that prepended index and use it to sum for each j
                e#     the first jth items of row j
    1+          e#     Append a 1 for P(i, i)
    a\+         e#     Prepend the new row to the table (which is stored in reverse)
  }*
  0=1b          e#   Sum the elements in the latest (first) row

  e&            e#   If p was 0 then replace with 0
}%
:*              e# Take the product

Answer 2

CJam, 50 49 47 43 bytes

ri_mF{1=_L{1$0>{,f{):X-Xj}:+}{;!}?}2j}%:*e&

Uses CJam's builtin mF factorisation and a memoised port of this Python partition number function:

p=lambda n,x:n==0 or n>0and sum(p(n+~a,a+1)for a in range(x))

or ungolfed:

def p(n, x): # Call like p(n, n). n is number remaining, x is max part size
  if n > 0:
    return sum(p(n-a-1,a+1)for a in range(x))
  else:
    return (n==0)

Like @RetoKoradi's answer, the last case takes about 17 seconds on the offline interpreter because that's how long it takes CJam to factorise the number. Hence I've left it out of this online test suite.

Full explanation

[Main body]
ri                                Read input and convert to int
  _          e&                   Logical AND input with final result to special case 0 
   mF                             Factorise input into [base, exponent] pairs
     {...}%                       Map, converting each pair to a partition number
           :*                     Take product

[Pair -> partition]
1=_                               Get exponent and copy (n,x in above Python)
   L                              Initialise empty cache
    {                       }2j   Memoise with 2 arguments
     1$0>                         Check if n > 0
         {            }{  }?      Execute first block if yes, else second block
                        ;!        Return (n == 0)
          ,f{      }              For each a in range(x) ...
             ):X-Xj               Call p(n-a-1,a+1) recursively
                    :+            Sum the results

Answer 3

Mathematica, 96 94 88 bytes

f=1##&@@#&;f[SeriesCoefficient[1/f[1-x^Range@#],{x,0,#}]&/@Last/@FactorInteger@#]Sign@#&

I'm not that proficient with Mathematica, but I thought I'd give it a try. Thanks to @MartinBüttner for -6 bytes.

This uses the generating function formula for integer partitions.

Answer 4

CJam, 58 bytes

li_mF{1=_L{_1>{_2$<{\;_j}{\,f{)_@\-j}:+}?}{;;1}?}2j}%:*\g*

Try it online

The very last test example takes forever (or at least longer than I was willing to wait) in the online interpreter, but finishes in 17 seconds with the offline version of CJam on my laptop. All other test examples are pretty much instantaneous.

This uses the CJam mF operator, which gives the prime factorization with exponents. The result is then the product of the partition counts for each exponent.

The main part of the code is calculating the partition counts. I implemented the recursive algorithm on the wikipedia page, using the j operator that supports recursion with memoization.

Explanation:

li    Get input and convert to int.
_     Make a copy to handle 0 special case at the end.
mF    Factorization with exponents.
{     Loop over factors.
  1=    Take exponent from [factor exponent] pair.
  _     Repeat it, recursive calls are initiated with p(n, n).
  L     Empty list as start point of memoization state.
  {     Start recursive block. Argument order is (m, n), opposite of Wikipedia.
    _1>   Check for n > 1.
    {     Start n > 1 case.
      _2$   Copy both m and n.
      <     Check for n < m.
      {     n < m case.
        \;    Pop m.
        _     Copy n.
        j     Make the p(n, n) recursive call.
      }     End n < m case.
      {     Main part of algorithm that makes recursive calls in loop.
        \,    Generate [0 1 ... m-1] range for k.
        f{    Start loop over k.
          )     Increment, since k goes from 1 to m.
          _     Copy k.
          @\    Rotate n to top, and swap. Now have k n k at top of stack.
          -     Subtract, now have k n-k at top of stack.
          j     Make the p(n-k, k) recursive call.
        }     End loop over k.
        :+    Sum up all the values.
      }?    Ternaray operator for n < m condition.
    }     End n > 1 case.
    {     n <= 1 case.
      ;;1   Pop m, n values, and produce 1 as result.
    }?    Ternary operator for n > 1 condition.
  }2j   Recursive call with memoization, using 2 values.
}%    End loop over factors.
:*    Multiply all values.
\     Swap original input to top.
g     Signum.
*     Multiply to get 0 output for 0 input.