2
$\begingroup$

Let's say I have a list with the prime factors of a number $n$ and their corresponding exponents. Is there a formula to calculate how many multiplications with $k$ distinct factors of $n$ are possible? For example, $12 = 2^2 \times 3$, so $F(12,2)$ would be:

$[2,6][3,4][1,12] = 3$

Note that the order of each set of factors does not matter. I have a feeling I only need the exponents of each prime factor, but I'm at a loss on how to actually proceed.

$\endgroup$
  • 1
    $\begingroup$ Trying to cheat on project euler? :) $\endgroup$ – Elliot Gorokhovsky Oct 30 '15 at 17:32
  • $\begingroup$ @RenéG, after two days of thinking about it for hours I had no option but to give up, unfortunately, haha. $\endgroup$ – Nicolás Siplis Oct 30 '15 at 17:37
  • $\begingroup$ Ya that's how I found your question in the first place. $\endgroup$ – Elliot Gorokhovsky Oct 30 '15 at 17:37
2
$\begingroup$

This problem is a straightforward application of the Polya Enumeration Theorem and very similar to the problem discussed at this MSE link.

Remark. What follows below does permit repeated factors, the formula for distinct factors is at the end.

Recall the recurrence by Lovasz for the cycle index $Z(S_n)$ of the multiset operator $\def\textsc#1{\dosc#1\csod} \def\dosc#1#2\csod{{\rm #1{\small #2}}}\textsc{MSET}_{=n}$ on $n$ slots, which is $$Z(S_n) = \frac{1}{n} \sum_{l=1}^n a_l Z(S_{n-l}) \quad\text{where}\quad Z(S_0) = 1.$$

We have for example, $$Z(S_3) = 1/6\,{a_{{1}}}^{3}+1/2\,a_{{2}}a_{{1}}+1/3\,a_{{3}}$$ and $$Z(S_4) = 1/24\,{a_{{1}}}^{4}+1/4\,a_{{2}}{a_{{1}}}^{2} +1/3\,a_{{3}}a_{{1}}+1/8\,{a_{{2}}}^{2}+1/4\,a_{{4}}.$$

Suppose the prime factorization of $n$ is given by $$n=\prod_p p^v.$$

Applying PET it now follows almost by inspection that the desired count is given by the term $$F(n, k) = \left[\prod_p X_p^v\right] Z(S_k)\left(\prod_p \frac{1}{1-X_p}\right)$$ where the square bracket denotes coefficient extraction of formal power series.

As an example of what these substituted cycle indices look like consider $$Z(S_3)\left(\frac{1}{1-X_1}\frac{1}{1-X_2}\frac{1}{1-X_3}\right) \\ = 1/6\,{\frac {1}{ \left( 1-X_{{1}} \right) ^{3} \left( 1-X_{{2}} \right) ^{3}\left( 1-X_{{3}} \right) ^{3}}} \\+1/2\,{\frac {1}{ \left( -{X_{{1}}}^{2}+1 \right) \left( -{X_{{2}}}^{2}+1 \right) \left( -{X_{{3}}}^{2}+1 \right) \left( 1-X_{{1}} \right) \left( 1-X_{{2}} \right) \left( 1-X_{{3}} \right) }} \\+1/3\,{\frac {1}{ \left( -{X_{{1}}}^{3}+1 \right) \left( -{X_{{2 }}}^{3}+1 \right) \left( -{X_{{3}}}^{3}+1 \right) }}.$$

It should be clear that coefficient extraction here is fast and efficient using the Newton binomial series which says that $$[Q^{\mu k}] \left(\frac{1}{1-Q^\mu}\right)^\nu = {k+\nu-1\choose \nu-1}.$$

We get for $F(n,1)$ the sequence $$1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,\ldots$$ which is as it ought to be.

For $F(n,2)$ we get the sequence $$1, 1, 1, 2, 1, 2, 1, 2, 2, 2, 1, 3, 1, 2, 2, 3, 1, 3, 1, 3,\ldots$$ which points us to OEIS A038548.

As another example consider $F(n,4)$ which yields $$1, 1, 1, 2, 1, 2, 1, 3, 2, 2, 1, 4, 1, 2, 2, 5, 1, 4, 1, 4,\ldots$$ which points us to OEIS A21320 where the above calculation is confirmed.

As a last example consider $F(n,5)$ which yields starting at $n=20$ $$4, 2, 2, 1, 7, 2, 2, 3, 4, 1, 5, 1, 7, 2, 2, 2, 9, 1, 2, 2,\ldots$$ which points us to OEIS A252665.

The following Maple code implements this cycle index computation.

with(combinat);
with(numtheory);

pet_cycleind_symm :=
proc(n)
option remember;

    if n=0 then return 1; fi;

    expand(1/n*add(a[l]*pet_cycleind_symm(n-l), l=1..n));
end;

pet_varinto_cind :=
proc(poly, ind)
local subs1, subs2, polyvars, indvars, v, pot, res;

    res := ind;

    polyvars := indets(poly);
    indvars := indets(ind);

    for v in indvars do
        pot := op(1, v);

        subs1 :=
        [seq(polyvars[k]=polyvars[k]^pot,
             k=1..nops(polyvars))];

        subs2 := [v=subs(subs1, poly)];

        res := subs(subs2, res);
    od;

    res;
end;

F :=
proc(n, k)
    option remember;
    local v, q, sind, res;

    v := op(2, ifactors(n));

    sind :=
    pet_varinto_cind(mul(1/(1-cat(`X`, q)), q=1..nops(v)),
                     pet_cycleind_symm(k));

    res := sind;
    for q to nops(v) do
        res := coeftayl(res, cat(`X`, q)=0, v[q][2]);
    od;

    res;
end;

Addendum. For distinct factors the formula is $$G(n, k) = \left[\prod_p X_p^v\right] Z(P_k)\left(\prod_p \frac{1}{1-X_p}\right)$$ where the square bracket denotes coefficient extraction of formal power series and $Z(P_n)$ is the cycle index of the set operator $\def\textsc#1{\dosc#1\csod} \def\dosc#1#2\csod{{\rm #1{\small #2}}}\textsc{SET}_{=n}$ which was also used in the linked-to computation from above.

Starting at $n=20$ we get for $G(n,3)$ the sequence $$2, 1, 1, 0, 4, 0, 1, 1, 2, 0, 4, 0, 2, 1, 1, 1, 4, 0, 1, 1,\ldots$$ which points us to OEIS A088434.

Starting at $n=120$ we get for $G(n,4)$ the sequence $$8, 0, 0, 0, 0, 0, 3, 0, 1, 0, 1, 0, 3, 0, 0, 1, 1, 0, 1, 0,\ldots$$

The Maple code for $G(n,k)$ is as follows.

with(combinat);
with(numtheory);

pet_cycleind_set :=
proc(n)
option remember;

    if n=0 then return 1; fi;

    expand(1/n*add((-1)^(l-1)*a[l]*pet_cycleind_set(n-l), l=1..n));
end;

# duplicate code omitted

G :=
proc(n, k)
    option remember;
    local v, q, sind, res;

    v := op(2, ifactors(n));

    sind :=
    pet_varinto_cind(mul(1/(1-cat(`X`, q)), q=1..nops(v)),
                     pet_cycleind_set(k));

    res := sind;
    for q to nops(v) do
        res := coeftayl(res, cat(`X`, q)=0, v[q][2]);
    od;

    res;
end;

Addendum. Observe that we can easily compute $$G(n,5)$$ where $$n=2^{50}\times 3^{20} \times 5^{10}$$ which has $11781$ divisors so that the brute force method would have to check $$1889560697901637761$$ five-tuples, which is not feasible. The answer is $$27879304592.$$

$\endgroup$
  • $\begingroup$ Slightly overwhelming but I will study the code to see if I can translate it into another language, a million thanks! $\endgroup$ – Nicolás Siplis Apr 24 '15 at 23:43
  • $\begingroup$ Being a Linux user and programmer I strongly encourage translating these MAPLE routines for PET computations and making them available for anyone to use. (The link has contributions from a number of users.) $\endgroup$ – Marko Riedel Apr 24 '15 at 23:50
  • $\begingroup$ I will certainly upvote a functional contribution in another language that would then confirm the results from the above post assuming I have not made any mistakes. $\endgroup$ – Marko Riedel Apr 24 '15 at 23:58
  • $\begingroup$ Hmm, interesting. It seems this algorithm can compute large values of n but it becomes inefficient when k increases. For example, (10!,30) hogs a lot of memory and takes too long to finish. Are there any particular optimizations that could be made to calculate when k = 30 and n is a factorial? I'm assuming that, since calculating the divisors of a factorial is trivial, there's a significant speed improvement that can be made. $\endgroup$ – Nicolás Siplis Apr 25 '15 at 16:00
  • $\begingroup$ The factorization is not a determining factor here, rather it is the coefficient extraction from the $5604$ terms in the cycle index $Z(P_{30})$. Fine-tuning this would appear to be a thesis type project involving a considerable effort of optimization and innovative coding. $\endgroup$ – Marko Riedel Apr 25 '15 at 20:35
-1
$\begingroup$

Any number in prime factorization would look like:

N = X_1^a * X_2^b * X_3^c * ... where X_1, X_2, X_3, ... are its prime factors

Notice that all powers of X_1 will be factors of N. So we can vary a from 0 through 1, 2, 3, ... to a. i.e. a can have a+1 values and for all these, X_1^a will be a factor. Similarly X_2 and X_3 and so on will be factors from 0 to b+1 and from 0 to c+1 etc. Also, we can choose X_1, X_2, X_3 etc. simultaneously. So, the total factors of N should be (a+1)(b+1)(c+1)*....

Since you want only pairs and any one factor of N can pair up with exactly one other factor to give N on multiplication, you simply divide number of factors by 2 and you have

[(a+1)(b+1)(c+1)....]/2

$\endgroup$
  • $\begingroup$ I think I understand the basic idea so maybe I'm applying it wrong, but wouldn't this approach for higher of k? For example, $144 = 2^4 \times 3^2$, so $F(144,4) = [(4+1)\times(2+1)] / 4 = 15 / 4 = 4,25$, which obviously since $F(144,4) = [1,2,4,18] [1,2,8,9] [1,2,3,24] [1,2,6,12] [1,3,4,12] [1,3,6,8] [2,3,4,6] = 7$ $\endgroup$ – Nicolás Siplis Apr 24 '15 at 21:29
  • $\begingroup$ Sorry, couldn't understand. What do you mean by 'would it work for higher values of k'? What is k? $\endgroup$ – Soham Apr 24 '15 at 21:33
  • $\begingroup$ I think I expressed myself wrong in the title. F(n,k) is the amount of ways in which you can express n as the product of k distinct numbers. $\endgroup$ – Nicolás Siplis Apr 24 '15 at 21:35
-2
$\begingroup$

Yes! From a number theoretical perspective, this is known as the divisor function. Or, even more generally, $\sigma_k (n)$ where the sigma function denotes the sum of all divisors of an integer to the power $k$. If we set $k=0$ we see that this is exactly the same as counting all possible divisors of this integer.

The general form of the divisor function (sometimes denoted $d(n)$) is as follows:

$$\prod_{p^\alpha || n} (\alpha +1)$$ where $p^\alpha || n$ means that $p^\alpha |n$ but $p^{\alpha +1}$ does not. Also, p denotes a prime divisor. Also, if the notation is getting to you, just thinking of it as factorizing an integer, and then taking the product of 1 more than all the exponents of each divisor. (Do an example!)

$\endgroup$
  • $\begingroup$ I think I didn't express myself as clear as I should have. The formula you gave, if I'm correct, shows all the ways in which you can express a number as a product of its divisors, right? However, what I'm looking for is in how many ways can a number be represented as the product of k distinct divisors. So F(144,4) would be the same as asking "How many unique sets (a,b,c,d) where a * b * c * d = 144 are there?" Remembering that permutations of the same set do not count and each element only appears once in each set. $\endgroup$ – Nicolás Siplis Apr 24 '15 at 21:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.