LOTERY - Editorial

editorial
, , , ,
1

Problem Link

Practice

Contest

Difficulty

Medium-Hard

Pre-requisites

Math, persistent segment tree

Problem

You are given a number triangle, which is denoted as:

  • L(n, 1) = n
  • L(n, k) = (\frac{1}{L(n - 1, k - 1)}-\frac{1}{L(n, k-1)})^{-1}

You are given T on-line queries on finding the value of LCM(L(n, 1), L(n, 2), ..., L(n, k)) modulo 10^9+7.

How to get 30 points

The elements of the triangle, denoted by L(n, k) are actually the denominators of the terms of Leibniz Harmonic Triangle. So, the value of L(n, k) equals to (n \cdot \binom{n-1}{k-1} ). So, the answer for the query denoted by integers n and k will be equal to LCM(L(n, 1), L(n, 2), ..., L(n, k)) = n \cdot LCM(\binom{n-1}{0}, \binom{n-1}{1}, ..., \binom{n-1}{k-1}).

Letā€™s use the following equality: (n+1) \cdot LCM(\binom{n}{0}, \binom{n}{1}, ..., \binom{n}{k})=LCM(n-k+1, n-k+2, ..., n+1). Using it, itā€™s easy to see that the answer to the query denoted as (n, k) is LCM(n-k+1, ..., n).

So now we need to come up with a way to calculate the LCM of the consecutive numbers.

For doing it, first of all, we will need a way to factorize integers fast. According to the statement, there wonā€™t be any integer greater than 10^6, so we can use the modified sieve that not only determines whether the number prime or not, but also gives itsā€™ greater prime divisor. For doing that, any time we mark a number as a composite one, we can also remember the prime divisor that was the evidence of itsā€™ compositeness.

For better understanding, have a look at the following pseudocode. Modification is highlighted in bold:

    for i = 2 to n :
        isPrime[i] = true, largestPrimeDivisor[i] = 0
    for i = 2 to n :
        if isPrime[i] :
            <b>largestPrimeDivisor[i] = i</b>
            j = i + i
            while j <= n :
                isPrime[j] = false
                <b>largestPrimeDivisor[j] = i</b>
                j = j + i

Now, when we need to get a prime factorization of an integer, we can simply divide it by itsā€™ largest prime divisor, until we get one and write out these divisors.

The pseudocode that does that:

    primeDivisors = [empty list]
    while n > 1 :
        append largestPrimeDivisor[n] to the end of primeDivisors
        n = n / largestPrimeDivisor[n]

After that, primeDivisors[] will contain the factorization of an integer n. It isnā€™t hard to make another list that will store each prime divisor exactly once along with itsā€™ power in the decomposition. The runtime of this factorization is O(\log n).

Now letā€™s apply it to the calculation of the Least Common Multiple of the range of numbers. Consider the prime factorization of LCM(n-k+1, n-k+2, ..., n). Obviously it will contain only prime numbers not exceeding n, and for each of these primes, the power in the LCM will be equal to the maximal itsā€™ power in the integers between n-k+1 and n.

So, the solution will looks as follows:

  • Calculate the prime factorization for each number between n-k+1 and n.
  • Using it, find the maximal occurrence power for each prime from 2 to n, inclusively.
  • The answer will be equal to the product of primes from 2 to n, each in itsā€™ maximal occurrence power.

The factorization of a single integer runs in O(\log n), in a single query we will factorize no more than M integers. Then, we will have to calculate the powers of no more than M integers, each power can be calculated with binary exponentation in O(\log n). So, the final complexity will be O(TM \log M). This solution is enough for getting 30 points by solving subtasks #1 and #2, but it is too slow for solving two last subtasks.

How to get 100 points

Letā€™s maintain the array D_n[] with the following property: for each k \leq n it holds that LCM(k, k+1, ..., n) = D_n[k] \cdot D_n[k+1] \cdot ... \cdot D_n[n].

Letā€™s have a few examples:

  • D_1 = [1]
  • D_2 = [1, 2]
  • D_3 = [1, 2, 3]
  • D_4 = [1, 1, 3, 4]

Suppose that D_{n-1} is constructed correctly. How do we get D_n using it?

For the beginning, letā€™s take D_n = D_{n-1} with D_n[n] = n. Letā€™s note that at some k s, the required property wonā€™t be held. More precisely, that will be all k s such that there exists some l \geq k such that GCD(l, n)>1 i.e. there is some integer l greater or equal than k such that it shares some prime divisors with n. In this case, some prime divisors will be counted twice in the LCM, giving the wrong result.

Now, how do we get rid from these extra powers?

For each prime number (say, p) letā€™s store the stack of the positions of such indexes of D_n[], at which the corresponding element of D_n[] is currently divisible by p. For each position, we can also calculate the power of p with which it will appear - this is basically the maximal number of times that the position number can be divided by p. For each prime number (say, p, again), itsā€™ stack will contain the positions in the decreasing order of the power of p that occurs at this position, having the position with the smallest power of p at the top.

So in order to get rid of counting the same factor twice in the array D_n[], we first decompose n into a product of primes and for each prime (let p be the corresponding prime number) we do the following:

  • First, we process all the positions from the corresponding stack that have the power of p less or equal to the power of p in n (i.e. the number of times that n can be divided by p without remainder). We divide these positions by p in their respective powers and pop them from stack. For example, if we have n=4, we will firstly have D_4=[1, 2, 3, 4]. But the stack of prime number 2 will have the information that at position 2, D_4[2] is divisible by 2. The power of 2 in 2 is obviously 1 and it is less than the power of 2 in the integer 4. So we divide D_4[2] by the power of 2 in 2 and remove the information that D_4[2] is divisible by 2 (because it is no longer divisible by 2). Finally, we get D_4 = [1, 1, 3, 4].
  • Then, there can be an entry in the corresponding stack with the power of p greater that the power of p in n. In this case, we divide the value of D_n[] at that position by p in the power it occurs in n. This way, the sought invariant will be held.
  • Finally, we add a new entry to the corresponding stack, containing the position n and the maximal power p in n.

The divisions can be made with the use of modular inverse, because the modulo 10^9+7 is a prime number.

It is very convenient to use persistent segment tree for calculating D_n[]. The segment tree should be capable of multiplying a single element and finding the product of the segment. Anytime, the wonā€™t be more than O(log M) elements in each stack, and each prime decomposition wonā€™t contain more than O(\log M) distinct prime numbers. So, it will take O(M \log^2 M) to precalculate D_1[], D_2[], ā€¦, D_M[].

Using the property of the array D_n[] you can find the LCM of the consecutive range of integers with a single query to the segment tree. So it will take O(\log M) to process a single query. So, the final complexity turns out to be O(T \log M + M \log^2 M).

Setterā€™s Solution

Can be found here

Testerā€™s Solution

Can be found here

1 Like
2

Could not understand the 100 point soln .
Can someone explain it properly or another way of solving thsi problem ?

3 
Well, There was another approach I used for solving this problem.

We can easily reduce the above question to following :
For each query of n and k, we need LCM(nāˆ’k+1,nāˆ’k+2,...,n).

Now, consider an example of LCM(4,5,6,7,8). Here, if we find highest powers of all the primes which divide any of the numbers given and take product, we can get LCM.
That is, for 2, we get 8; for 3, we get 3 ( 6 is divisible by 3 but no number is divisible by 9); for 5, we get 5; for 7 we get 7; for rest of the primes, we get 1.
Thus, LCM(4,5,6,7,8) = 8*3*5*7 = 840.
So, we build our solution on the similar basis.

Say we need to find LCM(A,A+1,A+2,....,B-1,B)
BRIEF SOLUTION:
 1. For all the primes <= sqrt(b) : we find highest power of that prime that divides a number which lies in [A,B].
 2. We find all the primes which divide a number in range [A,B]

EXPLANATION:

PRECALCULATIONS:
1. We can get all primes upto N in O(N) (Sieve of Eratosthenes).
2. We can also calculate PRIMEPRODUCT[N], where PRIMEPRODUCT[i] will give product of all the prime numbers <= i (mod 10^9 + 7) ( done in O(N) )
3. We can calculate Multiplicative Modular Inverse of all PRIMEPRODUCT[n] and store in MMI[n] ( done in O(No. of Primes upto N * log(10^9 + 7))

Now for each query (A,B), we do the following steps : 
 1. if( A == B)     return A;
 2. For each prime <= sqrt(B) : Get Highest power of prime that divides a number in [A,B]. We do this via following snippet :

   num = primes[i];
   highestpower = 1;
   while (num <= b) {
	     diva = ceil(a/ num);
	     divb = floor(b/num);
	     if(diva<=divb){
			highestpower = num;
             num *= primes[i];
             }else{
	             break;
             }
    }

	Correctness Proof : 

	say num * x = X;
	for X to lie in range [A,B], A<=X<=B
	=> A <= num * x <= B
	=> A/num <= x <= B/num<
	=> There should exist an integer x, s.t. A/num <=x<=B/num.
	   We are only concerened with integral values >= A/num and <=B/num
	=> ceil(A/num) <=x <= floor(B/num)
	   Thus, if diva <= divb; there exists a multiple of num in range [A,B].

 3. Now we are done with all primes <= sqrt(b). Thus, for rest of the primes, we know their hightest powers can be the numbers themselves (For query [50,100] , all the primes greater than sqrt(100) can't have a square of prime as a divisor of any number as they will exceed b.
 We can Calculate this by finding all primes in ranges [A/1 , B/1] , [A/2,B/2], .... [A/sqrt(b) , B/sqrt(b)] and avoiding any repitition. But this will give TLE as for this, the time complexity for each query can go upto sqrt(N). So, we use some little trick. We keep a two checks, lowcheck and highcheck. lowcheck tells that we have calculated highest powers for all primes <= lowcheck and highcheck tells that we have calculated highest powers for all primes >=highcheck.
 We initialize highcheck by infinity and lowcheck by sqrt(b)
 Now , study the following snippet : 

 int divisor = 1;
 while(highcheck > lowcheck ){
	high = floor(b/divisor);
	low = ceil(a/divisor);
	if (high >= highcheck) {
		high = highcheck-1;
	}
	if (low <= lowcheck) {
		low = lowcheck+1;
		highcheck = -INFINITY;
	}
	product = (PRIMEPRODUCT[high] * MMI[low-1]) mod 10^9 + 7;
        ans = (ans * product) mod 10^9 + 7
	highcheck = min(highcheck, low);
	divisor = max(divisor+1 , floor(b/highcheck));
 }

 Here, low and high gives the range for which we haven't calculated the primes till now. We get product for this by getting product for all primes <= high and dividing it by product of all primes < low.
 Initially we take divisor =1. We update highcheck as min(highcheck, low) because for each range [A/1 , B/1] , [A/2,B/2], .... [A/sqrt(b) , B/sqrt(b)], high of former range is less than high of later range. So, we are sure that for all primes p s.t. high < p < highcheck, contribution will be of none. we calculate from low to high. thus, we have checked all the primes >= low. 

 Now, we only need to see for high < highcheck as for others, we have already calculated. high = floor(b/divisor). Thus, divisor = max(divisor+1 , b/highcheck); (divisor + 1, to make sure it does not check for same divisor again and again).

 Hence, We can get the answer.

COMPLEXITY:
	PRECALCULATIONS: O(N) + O(N) + O(1000*log(10^9 + 7)) [1000 is the number of primes <= 10^5].
	FOR EACH QUERY : O(65*log(N)) [65 is number of primes upto sqrt(N)] + O(sqrt(N)).
THUS, overall complexity = O(N + 1000*log(10^9 + 7) + 65*Q*log(N) + Q*sqrt(N))

SOLUTION LINK

6 Likes
4

basically the idea is that for all primes <= sqrt(100000) , we can bruteforce them and for primes larger than that , their highest power will be just 1 , so we have to answer a query of the form (product of all distinct primes from l to r) where l = n - k + 1 and r = n , we this we can use a persistent segment tree , for more details about queries regarding distinct numbers , refer to august long challenge -> simple queries editorial . here is my implementation https://www.codechef.com/viewsolution/8442793

1 Like
5

I used a similar computation as @uhateme but somewhat simpler in the second part. I used a static decomposition in primes lower than sqrt(N) and higher than sqrt(N). For the ā€œsmallā€ primes i used the same method as @uhatme. for ā€œlargerā€ primes I can be shure that each prime occurs in each interval queried at most once.

So I precomputed a prefix-product array Ā§ for an array with the value of the prime at the prime position an d 1 otherwise. The product of large primes in an interval [a,b] can then be easily computed by

p[b]*p[a-1]^-1 mod 1000000007

2 Likes
6

Even i used kind of similar approach, for below sqrt(n) precomputed and for higer did in a loop.

7

This is the way I used for factorizing integers from 1 to m. This does the job with almost no extra cost in terms of time than the Sieve of Eratosthenes. It exploits the fact that the largest exponent of a prime P in an integer N is equal to ( 1 + largest exponent of the prime in N/P ).

  0 ## FAST PRIME FACTORIZATION OF ALL INTEGERS UPTO N
  1 def initSeive(m):
  2     global sv
  3     sv = [[True, {}, i] for i in range(m+1)]   # sv[i][0] = is prime (boolean), 
  4                                                # sv[i][1] = factorization representation (dict) 
  5                                                # sv[i][2] = residue of division by primes
  4     sv[1] = [False, {}, 1]
  5     sv[2] = [True, {2:1}, 1]
  6     lim = int(m**0.5)+2
  7     for i in range(2, lim):
  8         if sv[i][0]:
  9             sv[i][1][i]=1
 10             sv[i][2] = 1
 11             for j in range(2*i, m+1, i):
 12                 sv[j][0] = False
 13                 quotient = sv[j][2]//i
 14                 sv[j][1][i] = 1 + (sv[quotient][1][i] if i in sv[quotient][1] else 0)
 15                 sv[j][2] //= i**(sv[j][1][i])
 16 
 17     for i in range(2,m+1):
 18         if sv[i][2] != 1:
 19             sv[i][1][sv[i][2]] = 1
 20             sv[i][2] = 1
8

Can someone explain why this would be true?

(n+1)ā‹…LCM((^n_0),(^n_1),...,(^n_k))=LCM(nāˆ’k+1,nāˆ’k+2,...,n+1)

1 Like
9

To understand it better, lets move factor or (n+1) inside, so we get
LCM(C(n,0)(n+1), C(n,1)(n+1), C(n,2)(n+1),ā€¦,C(n,k)(n+1))
Now expand the terms according to formula C(n,r) = [n*(n-1)*(n-2)ā€¦(n-r+1)]/r!.
we get => [(n+1), (n+1)*n, (n+1)n(n-1)/2,ā€¦] Observe that for some denom 1.2.3ā€¦d, we have terms (n+1)n(n-1)ā€¦(n-d+1), i.e d+1 consecutive terms, which means that if we only consider numerator from (n+1).n.(n-1)ā€¦(n-d), there will be some number divisible by d, which means we donā€™t need to disturb the new term i.e (n-d+1), so these new terms (n+1),n,(n-1)ā€¦ contribute in LCM.

10

they contribute in LCM on the RHS, but what about the denominators on the LHS? I donā€™t understand your proof. And after googling a bit, I found this link: http://math.stackexchange.com/questions/1442/is-there-a-direct-proof-of-this-lcm-identity

11

What i meant is all terms in denom will cancel out by numerator, for that it is sufficient to use only first d terms, no the newly added term need not take part in division, means we can use it for lcm as it differs from all previous added numbers. It can also be shown using Proof by Induction.

12

How many elements does each D_i[] array have?
I think we also need to keep track of the updated segment trees values but what would be the time complexity? Is it included in the O(Mlog^2 M) of precalculation?

13

How can someone come up with the fact that the F(N,K) of the problem was related to Leibniz Harmonic Triangle?

14

Hi all, I have just a little a supplement.

  • Sorry or the delay, the proofs of the main equations is provided here
    Also, notice that it was not intended you to know about Leibniz Harmonic Triangle(but itā€™s good if you know). The relation between given recurrent formula and binomial coefficients can be found if you just write down some initial values of function, and guess about the pattern.

  • Soon after contest started I noticed that time complexity can be also improved. Iā€™d like to share this details with you(here). And now total time complexity is O(MlogM + TlogM).

  • Also, would like to say, that I didnā€™t want O(X*sqrt(X)) solution to be accepted(and my such a solution didnā€™t pass), but it appeared a lot of solutions(seems like at least more that half of all AC > 30) with such a complexity which passed(some of them are very fast, actually).

  • And the last one, for those how didnā€™t understand explanation below, let me give you some clarifications(and example of calculating D[]):

    D_i[] contains exactly i elements.

    stack[prime] contains all indexes of elements of D[] which are divisible by prime in decreasing order.

    for each n, i, D_n[i] <= i

    Consider array D_11[]: D = {1,1,1,1,1,1,7,4,9,10,11} // 1-based

    stack2 = {10, 8}
    stack[3] = {9}, stack[5] = {10}, stack[7] = {7}, stack[11] = {11}

    Adding 12:

      D[12] = 12

  p = 4:

  p /= 2, D[stack[2 ].top] /= 2 (and then remove from stack, because it's equal one)

  p /= 2, D[stack[2 ].top] = D[8] /= 2 -> D[8] = 2

  p = 3:

  p /= 3, D[9] /= 3 -> D[9] = 3

    Array D_12[] = {1,1,1,1,1,1,7,2,3,5,11,12}

    So LCM(8,9,10,11,12) = 2 * 3 * 5 * 11 * 12

1 Like
15

@ceilks can you a suggest a way to calculate p[b]*p[a-1]^-1 mod 1000000007 !

16

I think the correct time complexity for the authorā€™s solution is \mathcal{O}\left(N \log N \log M \sum_{p}\frac{1}{p\log p}\right) where p are primes are less than M. But since \sum_{p}\frac{1}{p\log p} \approx 1.55 for M = 10^5, we can say that the complexity is equivalent to \mathcal{O}(N\log N \log M)