Divisibility

题目描述

Inspired by Stephen Graham, the King of Berland started to study algorithms on strings. He was working days and nights, having a feeling that the full potential in this area is still to be unlocked. And he was right!
One day, all the sudden, he made a huge breakthrough by discovering the fact that strings can be magically transformed into integer numbers. It was so simple! You just have to map different letters to different digits and be careful enough not to introduce any leading zeroes.
Here is what he wrote in his textbook about the string ‘lalala’:

  • it can be transformed to an $282828$ by mapping ‘l’ to $2$, and ‘a’ to $8$,
  • it can also be transformed to $909090$ by mapping ‘l’ to $9$, and ‘a’ to $0$,
  • acouple of examples of invalid transformations are $050505$ (the resulting number has a leading zero), $333333$ (different letters are mapped to the same digit), $123456$ (no mapping to the original letters at all).

But then things started to become more interesting. Obviously, it was known from very beginning that a single string can potentially be mapped to a variety of different integer numbers. But the King couldn’t even imagine that all numbers produced by the same string pattern might have common properties!
For example, every single number that can be produced from string ‘lalala’ is always divisible by $259$, irrespective of the letter-to-digit mapping you choose. Fascinating!
So the King ended up with the following problem. For any given string, he wanted to come up with an algorithm to calculate the set of its divisors. A number is called a divisor of the given string if all positive integers, that could possibly be produced from the given string, are divisible by it.
As usual, the King desperately wants you to help him, so stop thinking and start acting!

题意概述

给定一个字符集大小不超过$10$的字符串$s$。定义一个数字是合法的,当且仅当它没有前导零,位数等于$|s|$,且它从高到低的第$i$位和第$j$位相等当且仅当$s_i=s_j$。求出所有能整除每一个合法数字的数。有$n$组数据。
数据范围:$1 \le n \le 100, \; 1 \le |s| \le 14$。

算法分析

如果我们求出了所有合法数字的GCD,那么答案就是GCD的所有约数。
为了方便,我们设字符串的下标从$1$开始。
令$g$等于所有合法数字的GCD。首先考虑它的下界。假设字符$i$在字符串中出现的次数为$a_i$,出现的位置是$p_{i, 1}, p_{i, 2}, \ldots, p_{i, a_i}$。定义字符$i$的价值$v_i=\sum_{j=1}^{a_i} 10^{|s|-p_{i, j}}$。显然,所有合法的数字都可以被表示为$\sum k_iv_i$,其中$k_i$表示字符$i$变成的数字。所以,$g$的下界为$(v_a, v_b, \ldots, v_z)$。
考虑$g$的上界。分两种情况讨论。

  • 引理 $\forall a \le b, \; (a, b) \mid b-a$。

当字符集大小小于$10$时,对于每一个$i$,总能找到一个合法的数$x$使得$x+v_i$也合法。这是因为总有一个合法的数包含$k_i$但不包含$k_i+1$。根据引理,$(x, x+v_i) \mid v_i$。由于$x$和$x+v_i$均合法,因此$g \mid (x, x+v_i)$,即$g \mid v_i$。所以,$g$的上界为$(v_a, v_b, \ldots, v_z)$。结合下界,$g=(v_a, v_b, \ldots, v_z)$。
当字符集大小等于$10$时,要将某个$k_i$加一就必须将另一个$k_j$减一。如果对于不同的$i, j$,有$v_i, v_j \gt 0 \land v_i \lt v_j$,那么$g \mid v_j-v_i$。令$S=\{v_j-v_i \mid v_i, v_j \gt 0 \land v_i \lt v_j\}$。对于一个合法的数,它加上/减去$S$中的某些数后仍是一个合法的数;对于$S$中的数,存在一个合法的数加上/减去它之后仍是一个合法的数。即一个合法的数$x$与$S$中的数加减可以得到所有合法的数。令$h$等于$S$中所有数的GCD,那么显然,$g$的下界为$(x, h)$。又因为$g \mid x$且$g$整除$S$中的所有数,所以$g$的上界也是$(x, h)$。即$g=(x, h)$。
求$g$的所有约数可以先分解质因数再DFS。

代码

#include <algorithm>
#include <cstdio>
#include <cstring>

static int const L = 15;
static int const M = 10000005;
char s[L];
int top, cnt, cc, p[M];
long long v[26], dv[L], dc[L], ans[M];

long long get_gcd(long long a, long long b) {
  for (long long c; b; c = a % b, a = b, b = c)
    ;
  return a;
}

void init() {
  top = 0;
  for (int i = 2; i < M; ++i) {
    if (!p[i])
      p[top++] = i;
    for (int j = 0; j < top; ++j) {
      int k = i * p[j];
      if (k >= M)
        break;
      p[k] = 1;
      if (!(i % p[j]))
        break;
    }
  }
}

void dfs(int t, long long v) {
  if (t == cnt)
    return void(ans[cc++] = v);
  for (int i = 0; i <= dc[t]; ++i, v *= dv[t])
    dfs(t + 1, v);
}

int main() {
  int n;
  scanf("%d", &n), init();
  for (int _ = 1; _ <= n; ++_) {
    printf("Case %d:", _), scanf("%s", s);
    int len = strlen(s), st = 0;
    for (int c = 0; c < 26; ++c) {
      v[c] = 0;
      for (int i = 0; i < len; ++i) {
        v[c] *= 10;
        if (s[i] == c + 'a')
          ++v[c];
      }
      st += v[c] > 0;
    }
    long long gcd = 0;
    if (st < 10)
      for (int i = 0; i < 26; ++i)
        gcd = get_gcd(gcd, v[i]);
    else {
      int k = 0;
      for (int i = 0; i < 26; ++i)
        if (v[i])
          gcd += k++ * v[i];
      for (int i = 0; i < 26; ++i)
        for (int j = 0; j < 26; ++j)
          if (v[i] && v[j] && v[i] < v[j])
            gcd = get_gcd(gcd, v[j] - v[i]);
    }
    cnt = cc = 0;
    for (int i = 0; i < top && 1ll * p[i] * p[i] <= gcd; ++i)
      if (!(gcd % p[i])) {
        dv[cnt] = p[i], dc[cnt] = 0;
        for (; !(gcd % p[i]);)
          gcd /= p[i], ++dc[cnt];
        ++cnt;
      }
    if (gcd > 1)
      dv[cnt] = gcd, dc[cnt] = 1, ++cnt;
    dfs(0, 1), std::sort(ans, ans + cc);
    for (int i = 0; i < cc; ++i)
      printf(" %lld", ans[i]);
    puts("");
  }
  return 0;
}

Clever Y

题目描述

Little Y finds there is a very interesting formula in mathematics:
$$
X^Y \bmod Z=K
$$
Given $X, Y, Z$, we all know how to figure out $K$ fast. However, given $X, Z, K$, could you figure out $Y$ fast?

题意概述

给定$X, Z, K$,求满足$X^Y \bmod Z=K$的$Y$的最小值。若不存在$0 \le Y \lt Z$满足条件,则输出无解。
数据范围:$0 \le X, Z, K \le 10^9$。

算法分析

显然$K \ge Z$时无解。在其他情况下,条件可以写成
$$
X^Y \equiv K \pmod Z
$$
特殊处理两种情况:$X, K$中至少有一个为$0$时,$Y=1$或无解;$K=1$时,$Y=0$。
先考虑$(X, Z)=1$的情况。令$M=\lceil \sqrt{Z} \rceil, \; Y=aM-b \; (0 \lt a, b \le M)$,那么
$$
\begin{align}
X^{aM-b} &\equiv K \pmod Z \\
X^{aM} &\equiv K \cdot X^b \pmod Z
\end{align}
$$
可以在$O(M)$的时间内预处理出$K \cdot X^b$,存在哈希表中,接着$O(M)$枚举$a$,判断哈希表中是否存在$X^{aM}$,若存在,则找到解$Y=aM-b$。由于要求最小值,因此哈希表中有重复时应记录较大的$b$。
下面考虑$(X, Z) \gt 1$的情况。令$D=(X, Z)$,易知若$K \nmid D$则无解。设$X=Dx, \; K=Dk, \; Z=Dz$,据同余的性质可得
$$
\begin{align}
(Dx)^Y &\equiv Dk \pmod{Dz} \\
x \cdot (Dx)^{Y-1} &\equiv k \pmod z
\end{align}
$$
因为$D=(X, Z)$,所以$(x, z)=1$,即$x$在模$z$意义下有逆元,因此
$$
(Dx)^{Y-1} \equiv k \cdot x^{-1} \pmod z
$$
令$X’=Dx=X, \; Y’=Y-1, \; K’=k \cdot x^{-1}, \; Z’=z$,则
$$
(X’)^{Y’} \equiv K’ \pmod{Z’}
$$
这和初始时的方程具有一样的形式。若$(X’, Z’) \gt 1$,则重复以上过程,否则可以参照$(X, Z)=1$的情况。

代码

/*
 * You will be singled out for promotion in your work.
 */

#include <algorithm>
#include <cmath>
#include <cstdio>
#include <cstring>

class HashTable {
private:
  static int const MOD = 3214567;
  int numr, h[MOD], id[MOD], val[MOD], nxt[MOD];

public:
  void clear() { memset(h, numr = 0, sizeof h); }

int count(int const &n) {
    for (int i = h[n % MOD]; i; i = i[nxt])
      if (i[id] == n)
        return 1;
    return 0;
  }

int &operator[](int const &n) {
    for (int i = h[n % MOD]; i; i = i[nxt])
      if (i[id] == n)
        return i[val];
    ++numr, numr[id] = n, numr[nxt] = h[n % MOD], h[n % MOD] = numr;
    return numr[val];
  }
} mp;

int ex_gcd(int a, int b, int &x, int &y) {
  if (!b) {
    x = 1, y = 0;
    return a;
  }
  int ret = ex_gcd(b, a % b, y, x);
  y -= a / b * x;
  return ret;
}

int get_gcd(int a, int b) {
  int x, y;
  return ex_gcd(a, b, x, y);
}

int get_inv(int a, int b) {
  int x, y;
  return ex_gcd(a, b, x, y), (x + b) % b;
}

int bsgs(int a, int b, int c) {
  a %= c;
  if (a == 0)
    return b == 0 ? 1 : -1;
  if (b == 0)
    return -1;
  if (b == 1)
    return 0;
  int k = 0;
  for (int t; (t = get_gcd(a, c)) > 1; ++k) {
    if (b % t)
      return -1;
    c /= t, b /= t, b = 1ll * b * get_inv(a / t, c) % c;
  }
  mp.clear();
  int d = a, m = sqrt(c) + 1;
  for (int i = 1; i <= m; ++i, d = 1ll * d * a % c)
    mp[1ll * d * b % c] = i;
  d = 1ll * d * get_inv(a, c) % c;
  int e = d;
  for (int i = 1; i <= m; ++i, e = 1ll * e * d % c)
    if (mp.count(e))
      return i * m - mp[e] + k;
  return -1;
}

int main() {
  for (int a, b, c; scanf("%d%d%d", &a, &c, &b), a || b || c;) {
    int ans = bsgs(a, b, c);
    if (~ans)
      printf("%d\n", ans);
    else
      puts("No Solution");
  }
  return 0;
}

Ceizenpok’s Formula

题目描述

Dr. Ceizenp’ok from planet i1c5l became famous across the whole Universe thanks to his recent discovery – the Ceizenpok’s formula. This formula has only three arguments: $n, k$ and $m$, and its value is a number of $k$-combinations of a set of $n$ modulo $m$.
While the whole Universe is trying to guess what the formula is useful for, we need to automate its calculation.

题意概述

求${n \choose k} \bmod m$。
数据范围:$1 \le n \le 10^{18}, \; 0 \le k \le n, \; 2 \le m \le 10^6$。

算法分析


$$
x={n \choose k}, \; m=p_1^{a_1}p_2^{a_2} \cdots p_q^{a_q}
$$
可以列出方程组
$$
\left\{
\begin{array}{c}
x \equiv r_1 \pmod{p_1^{a_1}} \\
x \equiv r_2 \pmod{p_2^{a_2}} \\
\cdots \\
x \equiv r_q \pmod{p_q^{a_q}}
\end{array}
\right.
$$
由于模数两两互质,所以该方程组在模$m$意义下有唯一解。
考虑如何求$r_i$。实际上,我们要求的就是${n \choose k} \bmod p_i^{a_i}$。我们知道
$$
{n \choose k}={n! \over k!(n-k)!}
$$
那么只要求出$n! \bmod p_i^{a_i}, \; k! \bmod p_i^{a_i}, \; (n-k)! \bmod p_i^{a_i}$的值,就可以用逆元求出${n \choose k} \bmod p_i^{a_i}$。
对于如何求$n! \bmod p_i^{a_i}$,令
$$
f(n)=n! \bmod p_i^{a_i}
$$
由$x \equiv x+p_i^{a_i} \pmod{p_i^{a_i}}$,可得
$$
f(n)=\left(f\left(\left\lfloor{n \over p_i}\right\rfloor\right)\cdot p_i^{\lfloor n/p_i \rfloor}\cdot\left(\prod_{i \in [1, p_i^{a_i}], \; p_i\not\mid i}i\right)^{\lfloor n/p_i^{a_i} \rfloor}\cdot\prod_{i \in [1, n \bmod p_i^{a_i}], \; p_i\not\mid i}i\right)\bmod p_i^{a_i}
$$
但是$k!, \; (n-k)!$在模$p_i^{a_i}$意义下可能不存在逆元,因此需要将$n!, \; k!, \; (n-k)!$中的$p_i$因子提取出来,求出逆元后再乘回去。
这样就得到了所有$r_i$。用中国剩余定理求解方程组即可。

代码

/*
 * Your lucky number has been disconnected.
 */

#include <algorithm>
#include <cstdio>
#include <cstring>

#define int long long

int power(int a, int b, int m) {
  int ret = 1;
  for (; b; b >>= 1)
    b & 1 && ((ret *= a) %= m), (a *= a) %= m;
  return ret;
}

void ex_gcd(int a, int b, int &x, int &y) {
  if (!b) {
    x = 1, y = 0;
    return;
  }
  ex_gcd(b, a % b, y, x), y -= a / b * x;
}

int get_inv(int a, int b) {
  int x, y;
  ex_gcd(a, b, x, y);
  return (x + b) % b;
}

int get_fac(int n, int p, int k) {
  int m = power(p, k, 1e9), ret = 1;
  for (; n; n /= p) {
    if (n / m) {
      int rec = 1;
      for (int i = 2; i < m; ++i)
        if (i % p)
          (rec *= i) %= m;
      (ret *= power(rec, n / m, m)) %= m;
    }
    for (int i = n % m; i > 1; --i)
      if (i % p)
        (ret *= i) %= m;
  }
  return ret;
}

int calc(int N, int K, int M, int p, int k) {
  int a = get_fac(N, p, k), b = get_fac(K, p, k), c = get_fac(N - K, p, k),
      cnt = 0;
  for (int i = N; i; i /= p)
    cnt += i / p;
  for (int i = K; i; i /= p)
    cnt -= i / p;
  for (int i = N - K; i; i /= p)
    cnt -= i / p;
  int m = power(p, k, 1e9),
      ret = a * get_inv(b, m) % m * get_inv(c, m) % m * power(p, cnt, m) % m;
  return ret * (M / m) % M * get_inv(M / m, m) % M;
}

signed main() {
  int N, K, M, ans = 0;
  scanf("%lld%lld%lld", &N, &K, &M);
  for (int i = 2, t = M; t > 1; ++i)
    if (!(t % i)) {
      int k = 0;
      for (; !(t % i); ++k, t /= i)
        ;
      (ans += calc(N, K, M, i, k)) %= M;
    }
  printf("%lld\n", ans);
  return 0;
}

Strange Way to Express Integers

题目描述

Elina is reading a book written by Rujia Liu, which introduces a strange way to express non-negative integers. The way is described as following:

  • Choose $k$ different positive integers $a_1, a_2, \ldots, a_k$. For some non-negative $m$, divide it by every $a_i$ to find the remainder $r_i$. If $a_1, a_2, \ldots, a_k$ are properly chosen, $m$ can be determined, then the pairs $(a_i, r_i)$ can be used to express $m$.

“It is easy to calculate the pairs from $m$,” said Elina. “But how can I find $m$ from the pairs?”
Since Elina is new to programming, this problem is too difficult for her. Can you help her?

题意概述

给定同余方程组
$$
\left\{
\begin{array}{c}
x \equiv r_1 \pmod{a_1} \\
x \equiv r_2 \pmod{a_2} \\
\cdots \\
x \equiv r_k \pmod{a_k}
\end{array}
\right.
$$
求$x$的最小非负整数解。
数据范围:所有输入输出均可以表示为64位整型。

算法分析

对于两个方程
$$
\begin{align}
x &\equiv r_1 \pmod{a_1} \\
x &\equiv r_2 \pmod{a_2}
\end{align}
$$
考虑将它们合并。易知
$$
x=r_1+k_1a_1=r_2+k_2a_2
$$
移项得
$$
k_1a_1=r_2-r_1+k_2a_2
$$
由贝祖定理可知,这个方程有整数解的充要条件为$(a_1, a_2) \mid r_2-r_1$。方程两边同时除以$(a_1, a_2)$
$$
\begin{align}
k_1{a_1 \over (a_1, a_2)}&={r_2-r_1 \over (a_1, a_2)}+k_2{a_2 \over (a_1, a_2)} \\
k_1{a_1 \over (a_1, a_2)} &\equiv {r_2-r_1 \over (a_1, a_2)} \pmod{{a_2 \over (a_1, a_2)}} \\
k_1 &\equiv \left( {a_1 \over (a_1, a_2)} \right)^{-1} \cdot {r_2-r_1 \over (a_1, a_2)} \pmod{{a_2 \over (a_1, a_2)}}
\end{align}
$$
将其代回$x=r_1+k_1a_1$,得
$$
x \equiv r_1+k_1a_1 \pmod{{a_1a_2 \over (a_1, a_2)}}
$$
至此已成功将两个方程合并为一个相同形式的方程。

代码

/*
 * Your life would be very empty if you had nothing to regret.
 */

#include <algorithm>
#include <cstdio>
#include <cstring>

#define int long long

int ex_gcd(int a, int b, int &x, int &y) {
  if (!b) {
    x = 1, y = 0;
    return a;
  }
  int ret = ex_gcd(b, a % b, y, x);
  y -= a / b * x;
  return ret;
}

int get_gcd(int a, int b) {
  int x, y;
  return ex_gcd(a, b, x, y);
}

int get_inv(int a, int b) {
  int x, y;
  return ex_gcd(a, b, x, y), (x += b) %= b;
}

signed main() {
  for (int k, a, r; ~scanf("%lld%lld%lld", &k, &a, &r);) {
    int flag = 0;
    for (int ai, ri; --k;) {
      scanf("%lld%lld", &ai, &ri);
      int gcd = get_gcd(a, ai);
      if ((r - ri) % gcd)
        flag = 1;
      r += get_inv(a / gcd, ai / gcd) * ((ri - r) / gcd) % (ai / gcd) * a;
      (a /= gcd) *= ai, ((r %= a) += a) %= a;
    }
    if (flag)
      puts("-1");
    else
      printf("%lld\n", r);
  }
  return 0;
}

The Sum of Unitary Divisors

题目描述

A natural number $d$ is a unitary divisor of $n$ if $d$ is a divisor of $n$ and if $d$ and ${n \over d}$ are coprime.
Let $\sigma^{\ast}(n)$ be the sum of the unitary divisors of $n$. For example, $\sigma^{\ast}(1)=1, \sigma^{\ast}(2)=3$ and $\sigma^{\ast}(6)=12$.
Let $S(n)=\sum_{i=1}^n \sigma^{\ast}(i)$.
Given $n$, find $S(n) \bmod 2^{64}$.

题意概述

定义$d$为$n$的单位因子当且仅当$d \mid n$且$(d, {n \over d})=1$。
令$\sigma^{\ast}(n)$为$n$的所有单位因子之和,$S(n)=\sum_{i=1}^n \sigma^{\ast}(i)$。给定$n$,求$S(n) \bmod 2^{64}$。
数据范围:$1 \le n \le 5 \times 10^{13}$。

算法分析

推式子
$$
\begin{align}
S(n)&=\sum_{i=1}^n \sigma^{\ast}(i) \\
&=\sum_{i=1}^n i \sum_{j=1}^{\lfloor n/i \rfloor} [(i, j)=1]
\end{align}
$$
利用莫比乌斯反演
$$
\begin{align}
S(n)&=\sum_{d=1}^n \mu(d) \sum_{i=1}^{\lfloor n/d \rfloor} di \sum_{j=1}^{\lfloor n/d^2i \rfloor} 1 \\
&=\sum_{d=1}^n d \cdot \mu(d) \sum_{i=1}^{\lfloor n/d \rfloor} i \cdot \left\lfloor {n \over d^2i} \right\rfloor
\end{align}
$$
显然,当$d^2 \gt n$时不产生贡献。那么
$$
\begin{align}
S(n)=\sum_{d=1}^{\lfloor \sqrt{n} \rfloor} d \cdot \mu(d) \sum_{i=1}^{\lfloor n/d^2 \rfloor} i \cdot \left\lfloor {n \over d^2i} \right\rfloor
\end{align}
$$
问题变成了求$F(n)=\sum_{i=1}^n i \cdot \mu(i)$和$G(n)=\sum_{i=1}^n i \cdot \left\lfloor {n \over i} \right\rfloor$。
对于$F$,我们有$(I \cdot \mu) \ast I=\varepsilon$。所以
$$
\begin{align}
F(n)=1-\sum_{i=2}^n i \cdot F\left(\left\lfloor {n \over i} \right\rfloor\right)
\end{align}
$$
这样$F$和$G$都可以用分块的技巧求出。而求$S$只要枚举$d$,将$\left\lfloor {n \over d^2} \right\rfloor$相等的$d$一起计算。

代码

#include <cmath>
#include <cstdio>

#define int long long

static const int N = 10000005;
static const int M = 50005;
int prime[N], mui[N], g[M];
bool vis[N];

void init() {
  int top = 0; mui[1] = 1;
  for (int i = 2; i < N; ++ i) {
    if (! vis[i]) prime[top ++] = i, mui[i] = -1;
    for (int j = 0; j < top; ++ j) {
      int k = i * prime[j]; if (k >= N) break; vis[k] = true;
      if (i % prime[j]) mui[k] = -mui[i]; else { mui[k] = 0; break; }
    }
  }
  for (int i = 2; i < N; ++ i) (mui[i] *= i) += mui[i - 1];
}

int get_sum(int n) {
  return n & 1 ? ((n + 1) >> 1) * n : (n >> 1) * (n + 1);
}

static const int HASH = 23456789;
int numr, h[HASH];
struct Record {
  int n, v, nxt;
} r[HASH];

void add_rec(int n, int v) {
  r[++ numr] = (Record) { n, v, h[n % HASH] }, h[n % HASH] = numr;
}

int get_f(int n) {
  if (n < N) return mui[n];
  for (int i = h[n % HASH]; i; i = r[i].nxt)
    if (r[i].n == n) return r[i].v;
  int ret = 1;
  for (int i = 2, j; i <= n; i = j + 1)
    j = n / (n / i), ret -= (get_sum(j) - get_sum(i - 1)) * get_f(n / i);
  return add_rec(n, ret), ret;
}

int get_g(int n) {
  if (n < M && g[n]) return g[n];
  int ret = 0;
  for (int i = 1, j; i <= n; i = j + 1)
    j = n / (n / i), ret += (get_sum(j) - get_sum(i - 1)) * (n / i);
  return ret;
}

int calc(int n) {
  int ret = 0;
  for (int i = 1, j; i * i <= n; i = j + 1)
    j = sqrt(n / (n / i / i) + 0.01), ret += (get_f(j) - get_f(i - 1)) * get_g(n / i / i);
  return ret;
}

signed main() {
  int T; scanf("%lld", &T), init();
  for (int i = 1; i < M; ++ i) g[i] = get_g(i);
  for (int n; T --;) scanf("%lld", &n), printf("%llu\n", calc(n));
  return 0;
}

GCD Extreme (Hard)

题目描述

Let $G(n)=\sum_{i=1}^n \sum_{j=i+1}^n (i, j)$.
For example, $G(1)=0, G(2)=(1, 2)=1, G(3)=(1, 2)+(1, 3)+(2, 3)=3$.
Given $N$, find $G(N)$ modulo $2^{64}$.

题意概述

令$G(n)=\sum_{i=1}^n \sum_{j=i+1}^n (i, j)$。给定$N$,求$G(N) \bmod 2^{64}$。
数据范围:$1 \le N \le 235711131719$。

算法分析

推式子
$$
\begin{align}
G(n)&=\sum_{i=1}^n \sum_{j=i+1}^n (i, j)=\sum_{i=2}^n \sum_{j=1}^{i-1} (i, j) \\
&=\sum_{d=1}^n d \sum_{i=2}^{\lfloor n/d \rfloor} \sum_{j=1}^{i-1} [(i, j)=1] \\
&=\sum_{d=1}^n d \sum_{i=2}^{\lfloor n/d \rfloor} \varphi(i)
\end{align}
$$
令$S(n)=\sum_{i=1}^n \varphi(i)$。由于$\varphi \ast 1=I$,因此
$$
\begin{align}
\sum_{i=1}^n (\varphi \ast 1)(i)&={n(n+1) \over 2} \\
&=\sum_{i=1}^n \sum_{j \mid i} \varphi(j) \\
&=\sum_{ij \le n} \varphi(j) \\
&=\sum_{i=1}^n \sum_{j=1}^{\lfloor n/i \rfloor} \varphi(j) \\
&=\sum_{i=1}^n S\left(\left\lfloor {n \over i} \right\rfloor\right)
\end{align}
$$
那么
$$
\begin{align}
S(n)&={n(n+1) \over 2}-\sum_{i=2}^n S\left(\left\lfloor {n \over i} \right\rfloor\right) \\
G(n)&=\sum_{d=1}^n d \cdot \left(S\left(\left\lfloor {n \over d} \right\rfloor\right)-1\right)
\end{align}
$$
接下来就可以利用分块的技巧计算。

代码

#include <map>
#include <cstdio>

#define int long long

static const int N = 10000000;
int prime[N], phi[N];
bool vis[N];

void init() {
  int top = 0; phi[1] = 1;
  for (int i = 2; i < N; ++ i) {
    if (! vis[i]) prime[top ++] = i, phi[i] = i - 1;
    for (int j = 0; j < top; ++ j) {
      int k = i * prime[j]; if (k >= N) break; vis[k] = true;
      if (i % prime[j]) phi[k] = phi[i] * (prime[j] - 1); else { phi[k] = phi[i] * prime[j]; break; }
    }
  }
  for (int i = 2; i < N; ++ i) phi[i] += phi[i - 1];
}

int get_sum(int n) { return n & 1 ? ((n + 1) >> 1) * n : (n >> 1) * (n + 1); }

std :: map <int, int> mp;

int get_phi(int n) {
  if (n < N) return phi[n];
  if (mp.count(n)) return mp[n];
  int ret = get_sum(n);
  for (int i = 2, j; i <= n; i = j + 1)
    j = n / (n / i), ret -= (j - i + 1) * get_phi(n / i);
  return mp[n] = ret;
}

int calc(int n) {
  int ret = 0;
  for (int i = 1, j; i <= n; i = j + 1)
    j = n / (n / i), ret += (get_sum(j) - get_sum(i - 1)) * (get_phi(n / i) - 1);
  return ret;
}

signed main() {
  int T; scanf("%lld", &T), init();
  for (int n; T --;) scanf("%lld", &n), printf("%llu\n", calc(n));
  return 0;
}

Jumping Joe

题目描述

Joe is a frog who likes to jump a lot. In fact, that’s all he does: he jumps forwards and backwards on the integer axis (a straight line on which all the integer numbers, both positive and negative are marked). At first, Joe sits next to the point marked with $0$. From here, he can jump in the positive or the negative direction a distance equal to either $x_1$ or $x_2$. From the point where he arrived, he can jump again a distance equal to $x_1$ or $x_2$, in the positive or the negative direction and so on.. Joe wants to arrive next to the point marked with the number $P$, after exactly $K$ jumps. You have to decide whether such a thing is possible.

题意概述

给定$x_1, x_2, P, K$,求一组非负整数$(P_1, N_1, P_2, N_2)$满足:
$$
\left\{
\begin{align}
P_1x_1-N_1x_1+P_2x_2-N_2x_2&=P \\
P_1+N_1+P_2+N_2&=K
\end{align}
\right.
$$
数据范围:$0 \lt x_1, x_2 \lt 40000, \; -40000 \lt P \lt 40000, \; 0 \le K \lt 2 \times 10^9$。

算法分析

设$g=(x_1, x_2)$,如果$g \not \mid P$则无解,否则可以将$x_1, x_2, P$同时除以$g$,对答案没有影响。
设$a=P_1-N_1, \; b=P_2-N_2$,则$ax_1+bx_2=P$,可用扩展GCD求出它的一组解。若$(a, b)$是一组解,则$(a+kx_2, b-kx_1)$也是一组解。只要在$\pm P$的范围内枚举$k$,判断是否存在$(P_1, N_1, P_2, N_2)$满足第二个方程即可。

代码

#include <cmath>
#include <cctype>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>

#define int long long

using namespace std;

struct IOManager {
  template <typename T> inline bool read(T &x) { char c = '\n'; bool flag = false; x = 0; while (~ c && ! isdigit(c = getchar()) && c != '-') ; c == '-' && (flag = true, c = getchar()); if (! ~ c) return false; while (isdigit(c)) x = (x << 3) + (x << 1) + c - '0', c = getchar(); return (flag && (x = -x), true); }
  inline bool read(char &c) { c = '\n'; while (~ c && ! (isprint(c = getchar()) && c != ' ')) ; return ~ c; }
  inline int read(char s[]) { char c = '\n'; int len = 0; while (~ c && ! (isprint(c = getchar()) && c != ' ')) ; if (! ~ c) return 0; while (isprint(c) && c != ' ') s[len ++] = c, c = getchar(); return (s[len] = '\0', len); }
  template <typename T> inline IOManager operator > (T &x) { read(x); return *this; }
  template <typename T> inline void write(T x) { if (x < 0) putchar('-'), x = -x; if (x > 9) write(x / 10); putchar(x % 10 + '0'); }
  inline void write(char c) { putchar(c); }
  inline void write(char s[]) { int pos = 0; while (s[pos] != '\0') putchar(s[pos ++]); }
  template <typename T> inline IOManager operator < (T x) { write(x); return *this; }
} io;

struct Solver {
private:
  int x1, x2, p, k;
  void input() { io > x1 > x2 > p > k; }
  int ex_gcd(int a, int b, int &x, int &y) {
    if (! b) { x = 1, y = 0; return a; }
    int ret = ex_gcd(b, a % b, y, x); y -= a / b * x; return ret;
  }
  bool check(int x, int y) {
    if (abs(x) + abs(y) > k) return false;
    if (k - abs(x) - abs(y) & 1) return false;
    io < (char *) "YES\n";
    int t = k - abs(x) - abs(y) >> 1;
    if (x > 0) io < x + t < ' ' < t < ' ';
    else io < t < ' ' < -x + t < ' ';
    if (y > 0) io < y < ' ' < 0 < '\n';
    else io < 0 < ' ' < -y < '\n';
    return true;
  }
  void process() {
    int x, y, gcd = ex_gcd(x1, x2, x, y);
    if (p % gcd) io < (char *) "NO\n";
    else {
      x1 /= gcd, x2 /= gcd, p /= gcd, x *= p, y *= p;
      for (int i = -abs(p) - 1; i <= abs(p) + 1; ++ i)
        if (check(x + x2 * i, y - x1 * i)) goto label;
      io < (char *) "NO\n";
      label: ;
    }
  }

public:
  void solve() { input(), process(); }
} solver;

signed main() {
  solver.solve();
  return 0;
}

Integer Sequences

题目描述

A sequence $A$ is called an integer sequence of length $N$ if all its elements $A_1, A_2, \ldots, A_N$ are non-negative integers less than $2000000000$. Consider two integer sequences of length $N, A$ and $X$. The result of their multiplication $(A \cdot X)$ is an integer number $R=\sum_{i=1}^N A_iX_i$. Your task is to solve the equation $A \cdot X \equiv B \pmod P$, given the integer sequence $A$ and the integer numbers $B$ and $P$.

题意概述

给定一个长度为$N$的向量$A$、两个整数$B, P$,求一个长度为$N$的向量$X$使得$A \cdot X \equiv B \pmod P$。向量中的数均为非负整数。
数据范围:$1 \le N \le 100, \; 0 \le B \lt P \le 10000$。

算法分析

设$g_i$表示$(A_1, A_2, \ldots, A_i)$,则有如下方程:
$$
\left\{
\begin{align}
A_1x_1+A_2y_1&=g_2 \\
g_2x_2+A_3y_2&=g_3 \\
\cdots& \\
g_{N-1}x_{N-1}+A_Ny_{N-1}&=g_N \\
g_Nx_N+Py_N&=(g_N, P)
\end{align}
\right.
$$
若$B \nmid (g_N, P)$,则无解;否则,解出所有$x_i, y_i$,即可还原出所有$X_i={B \over (g_N, P)} \cdot y_{i-1} \cdot \prod_{j=i}^N x_j$。

代码

#include <cctype>
#include <cstdio>

using namespace std;

struct IOManager {
  template <typename T> inline bool read(T &x) { char c = '\n'; bool flag = false; x = 0; while (~ c && ! isdigit(c = getchar()) && c != '-') ; c == '-' && (flag = true, c = getchar()); if (! ~ c) return false; while (isdigit(c)) x = (x << 3) + (x << 1) + c - '0', c = getchar(); return (flag && (x = -x), true); }
  inline bool read(char &c) { c = '\n'; while (~ c && ! (isprint(c = getchar()) && c != ' ')) ; return ~ c; }
  inline int read(char s[]) { char c = '\n'; int len = 0; while (~ c && ! (isprint(c = getchar()) && c != ' ')) ; if (! ~ c) return 0; while (isprint(c) && c != ' ') s[len ++] = c, c = getchar(); return (s[len] = '\0', len); }
  template <typename T> inline IOManager operator > (T &x) { read(x); return *this; }
  template <typename T> inline void write(T x) { if (x < 0) putchar('-'), x = -x; if (x > 9) write(x / 10); putchar(x % 10 + '0'); }
  inline void write(char c) { putchar(c); }
  inline void write(char s[]) { int pos = 0; while (s[pos] != '\0') putchar(s[pos ++]); }
  template <typename T> inline IOManager operator < (T x) { write(x); return *this; }
} io;

struct Solver {
private:
  static const int N = 101;
  int n, p, m, a[N], x[N], y[N];
  int ex_gcd(int a, int b, int &x, int &y) {
    if (! b) { x = 1, y = 0; return a; }
    int ret = ex_gcd(b, a % b, y, x); y -= a / b * x; return ret;
  }
  void input() {
    io > n > p > m;
    for (int i = 1; i <= n; ++ i) io > a[i], a[i] %= p;
  }
  void process() {
    int g = a[1]; x[0] = y[0] = 1;
    for (int i = 2; i <= n; ++ i) g = ex_gcd(g, a[i], x[i - 1], y[i - 1]);
    g = ex_gcd(g, p, x[n], y[n]);
    if (m % g) io < (char *) "NO\n";
    else {
      io < (char *) "YES\n", g = m / g;
      for (int i = 1; i <= n; ++ i) {
        int t = g * (y[i - 1] + p);
        for (int j = i; j <= n; ++ j) t = 1ll * t * (x[j] + p) % p;
        io < t < ' ';
      }
      io < '\n';
    }
  }

public:
  void solve() { input(), process(); }
} solver;

int main() {
  solver.solve();
  return 0;
}

Boxes

题目描述

There are two boxes. There are $A$ balls in the first box, and $B$ balls in the second box. It is possible to move balls from one box to another. From one box into another one should move as many balls as the other box already contains. You have to determine, whether it is possible to move all balls into one box.

题意概述

给定两个数$A$和$B$,每次操作可以将两个数中较大的数减去较小的数,再将较小的数翻倍。问至少经过多少次操作后两个数中有一个为$0$,若不可能则输出$-1$。
数据范围:$1 \le A+B \le 2147483647$。

算法分析

可以发现在任何时候,将两个数同时除以它们的最大公约数并不会影响答案;若两个数的和为奇数,则永远不可能结束。基于以上两点,答案不会超过$\log(A+B)$,因此直接模拟即可。

The Equation

题目描述

There is an equation $ax+by+c=0$. Given $a,b,c,x_1,x_2,y_1,y_2$ you must determine, how many integer roots of this equation are satisfy to the following conditions: $x_1 \le x \le x_2, \; y_1 \le y \le y_2$. Integer root of this equation is a pair of integer numbers $(x,y)$.

题意概述

求二元一次方程$ax+by+c=0$满足$x_1 \le x \le x_2 \land y_1 \le y \le y_2$的整数解的个数。
数据范围:$-10^8 \le a, b, c, x_1, x_2, y_1, y_2 \le 10^8$。

算法分析

先判断是否有解。若$a=b=0 \land c \neq 0 \lor (a, b) \nmid c$则无解。
之后我们计算出$g=(a, b)$,令$a’={a \over g}, \; b’={b \over g}, \; c’={c \over g}$,方程变为$a’x+b’y=-c’$,用扩展GCD求出它的一组解$(x, y)$,那么该方程的通解为$(x+kb’, y-ka’)$。
分别计算使$x, y$满足要求的$k$的范围,其交集的大小就是答案。在计算范围时可以二分查找。

代码

#include <cmath>
#include <cctype>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>

#define int long long

using namespace std;

namespace std {
  template <typename T>
  void maxify(T &a, T b) { b > a && (a = b); }
  template <typename T>
  void minify(T &a, T b) { b < a && (a = b); }
}

struct IOManager {
  template <typename T>
  inline bool read(T &x) {
    char c; bool flag = false; x = 0;
    while (~ c && ! isdigit(c = getchar()) && c != '-') ;
    c == '-' && (flag = true, c = getchar());
    if (! ~ c) return false;
    while (isdigit(c)) x = (x << 3) + (x << 1) + c - '0', c = getchar();
    return (flag && (x = -x), true);
  }
  inline bool read(char &c) {
    c = '\n';
    while (~ c && ! (isprint(c = getchar()) && c != ' ')) ;
    return ~ c;
  }
  inline int read(char s[]) {
    char c; int len = 0;
    while (~ c && ! (isprint(c = getchar()) && c != ' ')) ;
    if (! ~ c) return 0;
    while (isprint(c) && c != ' ') s[len ++] = c, c = getchar();
    return (s[len] = '\0', len);
  }
  template <typename T>
  inline IOManager operator > (T &x) {
    read(x); return *this;
  }
  template <typename T>
  inline void write(T x) {
    x < 0 && (putchar('-'), x = -x);
    x > 9 && (write(x / 10), true);
    putchar(x % 10 + '0');
  }
  inline void write(char c) {
    putchar(c);
  }
  inline void write(char s[]) {
    int pos = 0;
    while (s[pos] != '\0') putchar(s[pos ++]);
  }
  template <typename T>
  inline IOManager operator < (T x) {
    write(x); return *this;
  }
} io;

struct Solver {
private:
  int a, b, c, x1, x2, y1, y2, x, y, gcd;
  void input() { io > a > b > c > x1 > x2 > y1 > y2; }
  int ex_gcd(int a, int b, int &x, int &y) {
    if (! b) { x = 1, y = 0; return a; }
    int ret = ex_gcd(b, a % b, y, x); y -= a / b * x; return ret;
  }
  void process() {
    if (! b) swap(a, b), swap(x1, y1), swap(x2, y2);
    if (! b) {
      if (c) io < 0 < '\n';
      else io < (x2 - x1 + 1) * (y2 - y1 + 1) < '\n';
      return;
    }
    gcd = ex_gcd(a, b, x, y);
    if (c % gcd) { io < 0 < '\n'; return; }
    a /= gcd, b /= gcd, c /= gcd, x *= -c, y *= -c;
    if (b < 0) b = -b, a = -a;
    x += 1000000000 * b, y -= 1000000000 * a;
    y += ((x - x1) / b + 1) * a, x -= ((x - x1) / b + 1) * b;
    x += b, y -= a;
    if (x > x2) { io < 0 < '\n'; return; }
    int ans = (x2 - x) / b + 1;
    if (! a) {
      if (y < y1 || y > y2) ans = 0;
    } else {
      int l = -1000000000, r = 1000000000, cnt = 0;
      while (l + 1 < r) {
        int mid = l + r >> 1;
        if (y - a * mid < y1) if (a > 0) r = mid; else l = mid;
        else if (a < 0) r = mid; else l = mid;
      }
      if (a < 0 && y - a * l == y1) -- l, -- r;
      if (a > 0 && y - a * r == y1) ++ r, ++ l;
      int l1 = l, r1 = r; l = -1000000000, r = 1000000000;
      while (l + 1 < r) {
        int mid = l + r >> 1;
        if (y - a * mid < y2) if (a > 0) r = mid; else l = mid;
        else if (a < 0) r = mid; else l = mid;
      }
      if (a < 0 && y - a * r == y2) ++ r, ++ l;
      if (a > 0 && y - a * l == y2) -- l, -- r;
      if (l > l1) swap(l, l1), swap(r, r1);
      maxify(r, 0ll), minify(l1, ans - 1);
      ans = l1 - r + 1;
    }
    io < ans < '\n';
  }

public:
  void solve() { input(), process(); }
} solver;

signed main() {
  solver.solve();
  return 0;
}