grhkm's blog

By grhkm, history, 3 months ago, In English

Hello as title, any math topic suggestions that have few resources on, and are not too hard (i.e. it's explainable in one blog, perhaps with high-school maths knowledge)? IMO my previous blogs on maths (lagrange interpolation and linear modulo inversions) are pretty pog and received decent feedback, and I think (not to boast or anything) I am a good explainer? So I want to contribute with my grade 10 maths knowledge and also practice my english and writing skills in general.

It's free too for the lazy people out there, so you save your time and get your credit while I get my contribution points. Win-win condition (win for the cf readers too).

tldr please read blog title and reply thx

Read more »

  • Vote: I like it
  • +20
  • Vote: I do not like it

By grhkm, history, 3 months ago, In English

I am just wondering, how does your place's (country?) OI team selection work?

I am from Hong Kong and to join the HKOI selection process, we must register through our school. Unfortunately, my teacher did not confirm through the three confirmation emails and thus did not sign me up for the heat event, therefore I miss an annual chance to join the HKOI. Therefore, I am wondering if it is the same anywhere else, that you must register through your school? That sounds stupid to me but HKOI (hi tony) simply replied with "rules are rules".

In addition, I have joined the mirror of China OI first round selection and did decent. I asked HKOI if they can make an exception for me, especially since it wasn't even my fault (again, the teacher did not check the confirmation emails), but they won't let me join the second round of HKOI with those results.

Hope this makes sense, since I am probably leaving out details of the HKOI selection process.

Read more »

  • Vote: I like it
  • +101
  • Vote: I do not like it

By grhkm, history, 3 months ago, In English

PROBLEM: (Source)

Given $$$n$$$ $$$d$$$-dimensional vectors $$$v_1, v_2, \cdots, v_n$$$, output $$$i > j$$$ such that $$$v_i \cdot v_j \equiv 0 \mod k$$$, or determine that such indexes do not exist.

Constraints: For 50% cases $$$n \leq 10^5, d \leq 30, 2 \leq k \leq 3$$$, for all others $$$n \leq 10^4, d \leq 100, 2 \leq k \leq 3$$$.

$$$k = 2$$$

The main idea is check for every $$$i$$$ if there exists a $$$j$$$ satisfying the requirement. Keep in mind that here, $$$i$$$ can be treated as a constant. Therefore, requirement is equivalent to negating

$$$\forall 1\leq j\leq i-1, v_i \cdot v_j \not\equiv 0\mod 2\\$$$
$$$\forall 1\leq j\leq i-1, \sum_{k=1}^d v_{ik}v_{jk} \equiv 1\mod 2\\$$$
$$$\sum_{j=1}^{i-1}\sum_{k=1}^d v_{ik}v_{jk} \equiv i — 1\mod 2\\$$$
$$$\sum_{k=1}^d v_{ik}\sum_{j=1}^{i-1}v_{jk} \equiv i — 1\mod 2\\$$$

The inner summation can be computed quickly using prefix sum. Therefore complexity is O($$$nd$$$) where the $$$n$$$ comes from looping through every $$$i$$$.

$$$k = 3$$$

Now we deal with $$$k = 3$$$. Same idea, but now we don't have the nice property of $$$\not\equiv 0 \implies \equiv 1$$$. However, we can do something clever:

$$$\forall 1\leq j\leq i-1, v_i \cdot v_j \not\equiv 0\mod 3\\$$$
$$$\forall 1\leq j\leq i-1, (\sum_{k=1}^d v_{ik}v_{jk})^2 \equiv 1\mod 3\\$$$
$$$\sum_{j=1}^{i-1}\Big(\sum_{k=1}^d v_{ik}v_{jk}\Big)^2 \equiv i — 1\mod 3\\$$$
$$$\sum_{j=1}^{i-1}\sum_{k=1}^d \sum_{l=1}^d v_{ik}v_{jk}v_{il}v_{jl} \equiv i — 1\mod 3\\$$$
$$$\sum_{k=1}^d \sum_{l=1}^d v_{ik}v_{il} \Big(\sum_{j=1}^{i-1}v_{jk}v_{jl}\Big) \equiv i — 1\mod 3\\$$$

The inner summation can again be stored from previous queries. Therefore the time complexity is $$$O(nd^2)$$$.


1) The logic above implies that if $$$\sum \not\equiv i-1 \mod 3$$$, there will be a pair of vectors with inner product divisible by $$$3$$$. HOWEVER, it is possible that such pair exists but still bypasses all tests. For example taking $$$d = 1, k = 2, v_1 = 0, v_i = 1 \forall 2\leq i\leq n$$$, the vectors will satisfy that test above. Therefore, it is necessary to shuffle the vectors before applying the check.

2) In the original problem, we also have to restore the index. If we know a given $$$i$$$ fails, we can simply check all $$$j<i$$$ in $$$O(nd)$$$ time.

My implementation can be found here

Read more »

  • Vote: I like it
  • +26
  • Vote: I do not like it

By grhkm, history, 6 months ago, In English

I am messing around in (compiler), and I tried to compile this code with -O2

unsigned int func(unsigned int x) {  return x / 7; } 

Surprisingly, this is what I get

func:  mov eax, edi # edi is argument (x). Move x to eax (register)  imul rax, rax, 613566757 # eax *= 613566757 (unsigned 64-bit multiplication, rax is 64-bit version of eax)  shr rax, 32 # eax >>= 32 (gets upper 32 bit of product)  sub edi, eax # x -= eax  shr edi # x >>= 1  add eax, edi # eax += x  shr eax, 2 # x >>= 2  ret # return x 

What is imul rax, rax, 613566757 doing? Where does the constant come from?

Also, this translated code works as intended:

#include<bits/stdc++.h> using namespace std; typedef long long ll;

int main() { int i = 137;

long long edi = i;
    long long eax = i;
    eax *= 613566757LL;
    eax &gt;&gt;= 32;
    edi -= eax;
    edi &gt;&gt;= 1;
    eax += edi;
    eax &gt;&gt;= 2;
    cout &lt;&lt; i &lt;&lt; &quot; &quot; &lt;&lt; i / 7 &lt;&lt; &quot; &quot; &lt;&lt; eax &lt;&lt; endl;


Thank you very much.

(Originally I had int instead of unsigned int, but it's more complicated)

Read more »

  • Vote: I like it
  • +41
  • Vote: I do not like it

By grhkm, history, 6 months ago, In English

Thank you


We're going to learn how to find inverses mod p today (efficiently).

Pre-req: Know how to find inverses

Main results:

As we know, finding the inverse of n numbers is $$$O(n\log{p})$$$. That is too slow, especially when time limit is tight. Therefore, we want a faster way. I present:

  • Find inverse of all numbers between 1 and n in $$$O(n)$$$

  • Find inverse of n numbers in $$$O(n + \log{p})$$$

Idea 1:

Notice that $$$p = i \cdot \lfloor \frac{p}{i} \rfloor + (p \% i)$$$. (If you don't know why, compare it to $$$13=3\cdot 4+1$$$)

$$$p = i \cdot \lfloor \frac{p}{i} \rfloor + (p \% i)$$$
$$$\implies i \cdot \lfloor \frac{p}{i} \rfloor + (p \% i) \equiv 0 \mod p$$$
$$$p \% i \equiv -i \cdot \lfloor \frac{p}{i} \rfloor \mod p$$$

Now NOTICE that $$$p \% i$$$ is actually less than $$$i$$$. (Big brain mode)

So if we've calculated the inverses of $$$1-(i-1)$$$ already, we can find inverse of $$$(p \% i)$$$. This is what we're going to use.

$$$1 \equiv i \cdot (-\lfloor\frac{p}{i}\rfloor \cdot (p\% i)^{-1}) \mod p$$$

And therefore, by definition of inverse:

$$$i^{-1} \equiv (-\lfloor\frac{p}{i}\rfloor \cdot (p\% i)^{-1}) \mod p$$$

Or to make it easier to read: (Integer division)

$$$inv[i] \equiv (-(p/i) \cdot inv[p\% i]) \mod p$$$

So we can calculate inverse of $$$1~n$$$ in $$$O(n)$$$.

Code: (Be aware of overflow) (Also, -(p/i) mod p = (p-p/i) mod p)

int n = 10, p = 1000000007; int inv[n + 1]; inv[1] = 1; for (int i = 2; i <= n; i ++) inv[i] = 1LL * (p — p / i) * inv[p % i] % p; 

Idea 2:

Our target is $$$O(n+\log p)$$$, so we shall only take inverse once. How?

Let's look at what we want to find:

$$$\frac{1}{a_1}, \frac{1}{a_2}, \cdots, \frac{1}{a_n} \mod p$$$

Let's try to make them have the same denominator!

$$$\frac{a_2a_3\cdots a_n}{a_1a_2\cdots a_n}, \frac{a_1a_3a_4\cdots a_n}{a_1a_2\cdots a_n}, \cdots, \frac{a_1a_2\cdots a_{n-1}}{a_1a_2\cdots a_n} \mod p$$$

As we can see, the numerator is a prefix times a suffix, and the denominator is constant. Therefore, we can precompute those and then take inverse of denominator once.


#include "bits/stdc++.h" using namespace std;

const long long P = 1000000007; long long qpow(long long a, long long b) { long long ans = 1; while (b) { if (b & 1) ans = ans * a % P; a = a * a % P; b >>= 1; } return ans; }

long long n, a[1000010], pre[1000010], suf[1000010], pr = 1; // prefix, suffix, product int main() { cin >> n; pre[0] = suf[n + 1] = 1; for (int i = 1; i <= n; i ++) cin >> a[i]; for (int i = 1; i <= n; i ++) { pre[i] = pre[i — 1] * a[i] % P; suf[n + 1 — i] = suf[n + 2 — i] * a[n + 1 — i] % P; pr = pr * a[i] % P; } pr = qpow(pr, P — 2); for (int i = 1; i <= n; i ++) { cout << (pre[i — 1] * suf[i + 1] % P) * pr % P << endl; } }

Thank you for reading. If you have any suggestions or any topic you want to learn about I guess I can try writing! If this is helpful vote pls thx


Read more »

  • Vote: I like it
  • +139
  • Vote: I do not like it

By grhkm, history, 7 months ago, In English

Lagrange Interpolation by grhkm

No polynomial stuff because I don't know $$$O(n\log{n})$$$ polynomial multiplication :(

Main results:

Given $$$f(x_1)=y_1, f(x_2)=y_2, ..., f(x_n)=y_n$$$ we can calculate $$$f(X)$$$ where $$$f$$$ is the unique $$$(n-1)$$$-degree polynomial

For general {$$$x_i$$$} we can do so in $$$O(n^2)$$$

For consecutive {$$$x_i$$$} i.e. $$$x_j=x_1+(j-1)$$$, we can do so in $$$O(n)$$$ excluding binary exponentiation

Why useful?

Calculate $$$1^n+2^n+...+X^n$$$ in $$$O(n)$$$

DP optimization

With polynomials you can do cool things but not gonna cover this blog

Let's get started :)

Assume all polynomials below are $$$(n-1)$$$ degree unless specified.

Why (n-1) degree?


To start let's deal with first problem:

Given $$$f(x_1)=y_1, f(x_2)=y_2, ..., f(x_n)=y_n$$$ and an integer $$$X$$$ we can calculate $$$f(X)$$$ where $$$f$$$ is the unique $$$(n-1)$$$-degree polynomial

Let's consider an easier version:

Find polynomial satisfying $$$f(x_1)=y_i, f(x_2)=0, ..., f(x_n)=0$$$

As we learn in school, $$$f(r)=0 \implies (x-r)$$$ is divisor

We can write this down:

$$$f(x)=(x-x_2)(x-x_3)\cdots (x-x_n)$$$

Now notice that $$$f(x_1)=\prod_{i=2}^n (x_1-x_i)$$$ (a constant), but we want it to be $$$y_1$$$

We can simply scale the polynomial by a correct factor i.e.

$$$f_1(x)= \prod_{i=2}^n (x-x_i) \cdot \frac{y_1}{\prod_{i=2}^n (x_1-x_i)}=y_1\prod_{i=2}^n \frac{x-x_i}{x_1-x_i}$$$

Similarly, we can find a similar polynomial such that $$$f_i(x_i)=y_i$$$ and $$$f_i(x_j)=0 \forall j \neq i$$$

$$$f_i(x) = y_i \prod_{j=1,j\neq i}^n \frac{x-x_i}{x_i-x_i}$$$

Okay now finally to answer the original question, we can notice that $$$f(x)=\sum_{i=1}^n f_i(x)$$$ satisfies the constraints. So there we have it! Lagrange interpolation:

$$$f(x)=\sum_{i=1}^n y_i \prod_{j=1,j\neq i}^n \frac{x-x_j}{x_i-x_j}$$$

$$$O(n^2)$$$ excluding inverse for division


Let's try to optimise it now! Assume that we're given $$$f(1), f(2), \cdots, f(n)$$$, how can we calculate $$$f(X)$$$ quickly?

Notice that here, $$$x_i=i$$$ and $$$y_i=f(i)$$$. Substitute!

$$$f(x) = \sum_{i=1}^n f(i) \prod_{j=1,j\neq i}^n \frac{x-j}{i-j} = \sum_{i=1}^n f(i) \frac{\prod_{j=1,j\neq i}^n x-j}{\prod_{j=1,j\neq i}^n i-j}$$$

Let's consider the numerator and denominator separately


$$$\prod_{j=1,j\neq i}^n x-j = [(x-1)(x-2)\cdots (x-(i-1))] \cdot [(x-(i+1))(x-(i+2))\cdots (x-n)]$$$

We can pre-calculate prefix and suffix product of x, then we can calculate the numerator (for each $$$i$$$) in $$$O(1)$$$


$$$\prod_{j=1,j\neq i}^n i-j = [(i-1)(i-2)(i-3)\cdots(i-(i-1))] \cdot [(i-(i+1))(i-(i+2))\cdots (i-n)]$$$
$$$=(i-1)! \cdot (-1)^{n-(i+1)+1} \cdot (n-i)!$$$
$$$=(-1)^{n-i} (n-i)! (i-1)!$$$

So we can precompute factorials (preferably their inverse directly) then $$$O(1)$$$ calculation

So now we can calculate $$$f(X)$$$ in $$$O(n)$$$

Example 0

Find_ $$$\sum_{i=1}^N i^k$$$, $$$k\leq 1e6, N\leq 1e12$$$

Solution: Notice that the required sum will be a degree $$$(k+1)$$$ polynomial, and thus we can interpolate the answer with $$$(k+2)$$$ data points (Remember, we need $$$deg(f)+1$$$ points)

To find the data points simply calculate $$$f(0)=0, f(x)=f(x-1)+x^k$$$

Code can be found below

Example 1

(BZOJ 3453, judge dead): Find_ $$$\sum_{i=0}^n \sum_{j=1}^{a+id} \sum_{k=1}^j k^x \mod 1234567891, x \leq 123, a, n, d \leq 123456789$$$


Example 2 (ADVANCED)

(Luogu P4463 submit here): Given $$$n \leq 500, k \leq 1e9$$$, we call a sequence $$${a_n}$$$ nice if $$$1\leq a_i \leq k \forall i$$$ and $$$(a_i\neq a_j) \forall i\neq j$$$. Find the sum of (products of elements in the sequence) over all nice sequences. MY code below


Example 3 (ADVANCED TOO)

Problem : Find $$$\sum_{i=1,gcd(i,N)=1}^N i^k, N \leq 1e12, k \leq 256$$$. My code below


If you have any questions feel free to ask below

Practice question:

CF 622F (Example 0)

Luogu P4463 (Example 2)

Codechef XETF (Example 3)

Not in order of difficulty:


SPOJ ASUMHARD which is harder than EXTREME

Codechef COUNTIT

Atcoder ARC33D (Japanese)

Atcoder S8PC3G

  • Unfortunately you can't get full score, but you can get up to subtask 4

  • Still not trivial, do first 3 subtasks would give insight

  • For full solution please wait for my next blog (next year)

HDU 4059

SPOJ SOMESUMS very interesting

PE 487 Direct application

Read more »

  • Vote: I like it
  • +229
  • Vote: I do not like it