You are not logged in. Please login at www.codechef.com to post your questions!

×

TUPLES - Editorial

3
1

Practice
Contest

Author: Kevin Charles Atienza
Tester: Istvan Nagy
Editorialist: Misha Chorniy

Difficulty:

Hard

Pre-Requisites:

fft, number theory, inclusion-exclusion

Problem Statement

You are given array $X$ consisting $N$ elements. All the elements in the input are integers. 2-dimensional array $Y$ is defined in next way, $Y_{i,j}$ = $X_{i} * X_{j}$ mod $359999$. Find number of integer tuples $a, b, c, d, e, f$, where gcd($Y_{a,b}$, $Y_{c,d}$, $Y_{e,f}$) = 1 modulo $10^9+7$. Denote that $gcd(0, 0) = 0$.

Explanation

The first part of the solution is using not 2-dimensional array $Y$, instead of that calculate array $Z$ where $Z_{i}$ denotes the frequency of element $i$ in array $Y$, $i$ will be in a range between 0 and 359998. After that problem looks in a next way: we need to find the sum of all $Z_{i} * Z_{j} * Z_{k}$, where $gcd(i, j, k) = 1$.

Subtask 1

$N$ is less or equal than 1000. Let's create array $Y$ in a naive way in $O(N^{2})$. Also, simply calculate array $Z$ in $O(N^{2})$.

Subtask 2

A problem of the solution for subtask 1 in a too slow calculation of array $Y$, how to improve it? Let's consider number $359999$ in more detail, $359999 = 599 * 601$, both numbers $599$ and $601$ are primes. What can this fact give to us? Key observation: let we have a prime number $P$, from number theory follows next fact, any prime $P$ has primitive root $G$, in other words: For every $x$ = $1..P-1$, exists y, such $G^y$ = $x$, $0 <= y < (P-1)$.

Assume that we have the similar problem but with modulo $P$(not 359999), where $P$ is a prime number. Let's find the primitve root for P, call it G. Problem can be divided into two subproblems: where $X_{i}$ is zero or $X_{i}$ isn't zero. Because zero can't be represented using the primitive root of number $P$. $f(X)$ denoting such number that $G^{f(X)} = X$ and $ 0 <= f(X) < (P-1)$

$Y_{i,j}$ = $X_{i}$ * $X_{j}$ = $G^{f(X_{i})}$*$G^{f(X_{j})}$ = $G^{(f(X_{i}) + f(X_{j}))}$, in case if $X_{i} != 0$ and $X_{j} != 0$.

There are 4 subproblems:

  • $0 * 0 = 0$
  • $0 * B = 0$, $B mod P > 0$
  • $A * 0 = 0$, $A mod P > 0$
  • $A * B = G^{f(A)} * G^{f(B)} = G^{f(A)+f(B)}$, $A mod P > 0, B mod P > 0$
  • zeroes = 0 for i = 1..N if X[i] % P == 0 zeroes += 1 Z[0] = zeroes * zeroes + zeroes * (N - zeroes) + (N - zeroes) * zeroes for i = 1..N if X[i] mod P == 0 continue for j=1..N if X[j] mod P == 0 continue da = f(X[i]) // G^da = X[i], 0<=da<(P-1) db = f(X[j]) // G^db = X[j], 0<=db<(P-1) dab = da + db // X[i]X[j] = G^(da+db) Z[(G^dab)%P] += 1

    How to calculate function $f$ in $O(1)$. Let's make some precomputations, $ff$[$G^i mod P$] = i, for i in the range between $0$ and $P-2$. $f(i)=ff[i]$. Let's rewrite this code a bit, namely precalculate $C_{i}$ - frequency of the number $i$ in the array $X$.

    for i=1..N if X[i]%P==0 C[X[i] % P] += 1 else C[f(X[i])] += 1 Z[0] = C[0] * C[0] + 2 * C[0] * (N - C[0]) for i = 0..P-2 for j = 0..P-2 Z[(G^(i+j))%P] += C[i] * C[j]

    What we see, last code is very similar to multiplication of two polynomials

    for i = 1..N for j = 1..M C[i + j] += A[i] * B[j] //Standard multiplication of polynomials

    Try to change the code above a bit

    for i = 0..P-2 for j = 0..P-2 F[i + j] += C[i] * C[j] for i = 0..2
    (P-2) Z[(G^i)%P] += F[i]

    How to multiplicate two polynomials in way faster than $O(N^2)$, there are many algorithms, let's choose Fast Fourier Transform. Below pseudocode of multiplication using Fast Fourier Transform.

    D = FFT(C) for i = 0..2 * (P-2) E[i] = D[i] * D[i] F = inverseFFT(E) for i = 0..2*(P-2) Z[(G^i)%P] += F[i]

    But we have two primes 599 and 601, what to do in this case? This solution can be modified for the case with two primes, find primitive roots for 599 and 601, let's call them $G$ and $H$.

    Now every number between 0 and 359998 can be in one from the four types:

  • (0, 0), only 0 satisfied this condition
  • ($G^i mod 599$, 0), only 598 numbers satisfied this condition
  • (0, $H^i mod 601$), only 600 numbers satisfied this condition
  • ($G^i mod 599$, $H^j mod 601$)
  • How to multiply two numbers from these types, we can multiply first and second components independently from each other. In other words: $(A,B)*(C,D)$=$(A*C mod 599, B*D mod 601)$ from Chinese Remainder Theorem Consider all possible multiplications from types:

  • $(0, 0) = (0, 0) * (A, B)$ OR $(A, B) * (0, 0)$, $0 <= A <= 598, 0 <= B <= 600$
  • $(0, 0) = (0, A) * (B, 0)$ OR $(C, 0) * (0, D)$, $0 < A, D <= 598, 0 < B, C <= 600$
  • $(0, E) = (0, A) * (0, B)$ OR $(0, A) * (C, B)$ OR $(C, A) * (0, B)$, $0 < C <= 598, 0 < A, B, E <= 600$
  • $(E, 0) = (A, 0) * (B, 0)$ OR $(A, 0) * (B, C)$ OR $(A, C) * (B, 0)$, $0 < A, B, E <= 598$, $0 < C <= 600$
  • $(X, Y) = (A, B) * (C, D)$, $0 < A, C, X <= 598, 0 < B, D, Y <= 600$
  • First four convolutions can be processed in time $O(599 * 601) = O(359999)$, consider last convolution(it will be between two multipliers of fourth type):

    $(X, Y) = (A, B) * (C, D)$ = $(G^{i1} mod 599, H^{j1} mod 601) * (G^{i2} mod 599, H^{j2} mod 601)$ = $(G^{i1+i2} mod 599, H^{j1+j2} mod 601)$

    Let's create polynomial $P$, where $P$[$i * 2 * 601 + j$] is the number of numbers $G^{i}*H^{j} mod 359999$, if we'll multiply it by itself with FFT, we can take data of $(X, Y)$

    $(i1 * 2 * 601 + j1)$ * $(i2 * 2 * 601 + j2)$ = $(i1 + i2) * 2 * 601 + (j1 + j2)$, $0 <= (j1 + j2) < 2 * 601$

    for i = 0..598 for j = 0..600 P[i * 2 * 601 + j] += A[i][j] // number of numbers which represented in type G^i * H^j % 359999 fft(P) for i = 0..(1<<21)-1 D[i] = P[i] * P[i] //Multiplying with FFT Q[i] = inverseFFT(D) for i = 0..4*359999 i12 = i / (2 * 601) j12 = i % (2 * 601) Z[G^(i12)%599 * H^(j12)%601] += Q[i]

    After we'll find array $Z$ in fast way. How to speed up the following code:

    ans=0 for i = 0..359999-1 for j = 0..359999-1 for k = 0..359999-1 if gcd(i, j, k) = 1 ans += Z[i] * Z[j] * Z[k] ans %= 1000 * 1000 * 1000 + 7 Let's use some inclusion-exclusion, more precisely mobius-function, you can see this tutorial for similar problem: But in this problem, we don't need to order $i, j, k$, therefore for i = 1..359999 for j = i;j <= 359999; j+=i D[i] += Z[j]; ans += mu[i] * D[i] * D[i] * D[i] //number of positive integer triplets, where each number is divisible by i multiplied by mobius-function from i We forget about zeroes. How to count it:
  • $i = j = k = 0$, $gcd(i, j, k) = 0$, skip it
  • $i * j = 0$ OR $j * k = 0$ OR $i * k = 0$, we need to add $Z_{1} * Z_{0} * Z_{0}$
  • $i = 0$ OR $j = 0$ OR $k = 0$, we need to find the number of pairs $(a, b)$, where $gcd(Z_{a}, Z_{b}) = 1$ and multiply it by value of $Z_{0}$
  • ans = 3 * Z[1] * Z[0] * Z[0] for i=1..359999 for j=i;j<=359999; j+=i D[i]+=Z[j]; ans += mu[i] * D[i] * D[i] * D[i] //number of positive integer triplets, where each number is divisible by i multiplied by mobius-function from i ans += 3 * mu[i] * D[i] * D[i] //Don't forget about multipling by 3

    Overall time complexity of this approach is $O(2^{M} * M + X log X)$, where M = 21, X = 359999.

    Solution:

    Setter's solution can be found here
    Tester's solution can be found here

    Please feel free to post comments if anything is not clear to you.

    This question is marked "community wiki".

    asked 18 Jan, 07:56

    mgch's gravatar image

    7★mgch
    124214
    accept rate: 7%

    edited 18 Jan, 20:13

    Very nice Editorial.

    (19 Jan, 00:39) likecs6★
    1

    A very good problem and editorial! Though i could not solve it during the contest, I can definitely solve it now on my own :D

    (19 Jan, 02:19) swetankmodi5★

    There's a slight difference in my solution from the one described in this editorial, particularly in the part where you have to count the occurrences of $Y_{i,j}$ coprime to $m = 359999 = 599\cdot 601$.

    I'll introduce the notation I'll use. Let $x_v$ be the number of indices $i \in [1, N]$ such that $X_i = v$. Let $y_v$ be the number of index pairs $(i,j) \in [1, N]^2$ such that $Y_{i,j} = v$. Our goal is to compute the values $y_v$ for all $v \in \left(\mathbb{Z}/m\mathbb{Z}\right)^\times$. Here, $\left(\mathbb{Z}/m\mathbb{Z}\right)^\times$ denotes the multiplicative group of integers mod $m$ and consists of numbers coprime to $m$. This group has $\phi(m)$ elements.

    I'll also mention that this explanation requires some standard abstract algebra facts. Nothing too deep, though.


    If there were a number $g$ that can generate all numbers in $\left(\mathbb{Z}/m\mathbb{Z}\right)^\times$, then we can use FFT to compute the $y_v$ by rearranging the elements of $\left(\mathbb{Z}/m\mathbb{Z}\right)^\times$ according to their logarithms base $g$. Specifically, every element of $\left(\mathbb{Z}/m\mathbb{Z}\right)^\times$ is of the form $g^i$ for $i \in [0, \phi(m))$, and we have the equality:

    $$ y_{g^i} = \sum_{j=0}^{\phi(m)-1} x_{g^j} \cdot x_{g^{i-j}} $$

    for all $i \in [0, \phi(m))$. This is a convolution and can be computed in $O(\phi(m) \log \phi(m))$ time, and the details of this convolution has already been described in the editorial.

    But the generator $g$ doesn't exist. In fact, it only exists if $m$ is one of the forms $2$, $4$, $p^k$, or $2p^k$, for $k > 0$ and $p$ an odd prime.

    Our case $m = 359999 = 599\cdot 601$ isn't of that form. But we can modify the algorithm above to work on our $m$ by using a near generator instead. To explain what I mean by this, note that the most that any number $g$ can generate mod $m = 359999$ is $\phi(m)/2$, or exactly half the numbers. Let's call such a number a near generator. An example of a near generator is $g = 7$. This number generates exactly half of the numbers in $\left(\mathbb{Z}/m\mathbb{Z}\right)^\times$. Note that $g^{\phi(m)/2} \equiv 1 \pmod{m}$.

    An example of a number that $g = 7$ cannot generate is $17$. Let's define $h = 17$.

    Let's denote the set of numbers that $g$ can generate as $(g)$. Notice that:

    • $(g)$ is a subgroup of $\left(\mathbb{Z}/m\mathbb{Z}\right)^\times$, so if we form the cosets of $(g)$, we get two sets: $(g)$ itself and $h\cdot (g)$. This implies that every element of $\left(\mathbb{Z}/m\mathbb{Z}\right)^\times$ is exactly one of the following forms: $g^i$ or $hg^i$, for $i \in [0, \phi(m)/2)$.
    • $(g)$ is a normal subgroup of $\left(\mathbb{Z}/m\mathbb{Z}\right)^\times$, which means we can form the quotient $\left(\mathbb{Z}/m\mathbb{Z}\right)^\times/(g)$. Furthermore, it must be isomorphic to $\mathbb{Z}/2\mathbb{Z}$, the only group of size 2 (up to isomorphism). This implies that $(h\cdot (g))\cdot (h\cdot (g)) = (g)$. This is the same as saying $h^2 \in (g)$. This is the same as saying $h^2 = g^H$ for some $H \in [0, \phi(m)/2)$. (In fact, $H = 21370$.)

    Using these facts, we can now compute $y_v$ for various kinds of $v$ using the following:

    $$\begin{align*} y_{hg^i} &= \sum_{j=0}^{\phi(m)/2-1} x_{g^j} \cdot x_{hg^{i-j}} + \sum_{j=0}^{\phi(m)/2-1} x_{hg^j} \cdot x_{g^{i-j}} \\\ y_{g^i} &= \sum_{j=0}^{\phi(m)/2-1} x_{g^j} \cdot x_{g^{i-j}} + \sum_{j=0}^{\phi(m)/2-1} x_{hg^j} \cdot x_{hg^{i-j-H}} \end{align*}$$

    for all $i \in [0, \phi(m)/2)$. Now, we get 4 convolutions which can all be computed with FFT!

    link

    answered 22 Jan, 04:58

    kevinsogo's gravatar image

    7★kevinsogo
    1.7k471134
    accept rate: 11%

    Please upload the full editorial.

    link

    answered 18 Jan, 18:05

    likecs's gravatar image

    6★likecs
    3.1k633
    accept rate: 9%

    toggle preview
    Preview

    Follow this question

    By Email:

    Once you sign in you will be able to subscribe for any updates here

    By RSS:

    Answers

    Answers and Comments

    Markdown Basics

    • *italic* or _italic_
    • **bold** or __bold__
    • link:[text](http://url.com/ "title")
    • image?![alt text](/path/img.jpg "title")
    • numbered list: 1. Foo 2. Bar
    • to add a line break simply add two spaces to where you would like the new line to be.
    • basic HTML tags are also supported
    • mathemetical formulas in Latex between $ symbol

    Question tags:

    ×11,354
    ×779
    ×442
    ×194
    ×123
    ×64
    ×12
    ×10

    question asked: 18 Jan, 07:56

    question was seen: 1,361 times

    last updated: 22 Jan, 04:58