莫比乌斯反演入门 – [SDOI2015]约数个数和

道理我都懂,考试时这东西能推?

本人菜鸡,有问题请指出

鉴于不同博客对于整除符号$|$的定义不同, 特此表明本博客的整除定义

$a | b$ 表明 $b$ 是 $a$ 的倍数

0x01 什么是反演

对于一个数列${f_n}​$,如果有另外一个数列${g_n}​$满足如下条件
$$
g_n = \sum_{i = 1}^n a_if_i
$$
反演的过程是用$ g_n$来表示 $f_n$
$$
f_n = \sum_{i = 0}^n b_ig_i
$$

0x02 莫比乌斯函数

设$n = p_1^{a_1} \times p_2^{a_2} \times p_3^{a_3} \cdots \times p_k^{a_k}$ ,其中 $p$ 为质数,则定义如下
$$
\mu(n) =
\begin{cases} 1 & n = 1 \\
(-1) ^ m & \prod\limits_{i = 1} ^ {m} a_i = 1 \\
0 & \textrm{otherwise}(a_i \gt 1)
\end{cases}
$$

性质 1

$\mu​$ 函数是积性函数

没啥可以说的,线性筛就完事了

void get_pri(){
    vis[1] = 1; mu[1] = 1;
    for(int i = 2; i <= N; i++){
        if( vis[i] == false ){
            pri[ ++pcnt ] = i;
            mu[i] = -1;
        }
        for(int j = 1; j <= pcnt; j++){
            if( pri[j] * i > N ) 
                break;
            vis[ i * pri[j] ] = true;
            if( i % pri[j] == 0 ){
                mu[ pri[j] * i ] = 0;
                break;
            }
            else 
                mu[ pri[j] * i ] = -mu[i];    
        }
    }
}

性质 2

对于任意的非负整数
$$
\sum_{d|n}\mu(d)=
\begin{cases}
1 & n = 1\\
0 & n > 1
\end{cases}
$$

0x03 莫比乌斯反演

有两个数论函数$f(n)$和$g(n)$ ,且满足条件
$$
g(n) = \sum_{d|n}f(d)
$$
则一定存在
$$
f(n) = \sum_{d|n}\mu(d)g(\lfloor \frac{n}{d}\rfloor)
$$

证明

我们这里用性质证明QAQ
$$
\begin{align}
f(n) & = \sum_{d|n}\mu(d)g(\lfloor \frac{n}{d} \rfloor)\\
& = \sum_{d|n}\mu(d) \sum_{d’|\lfloor \frac{n}{d} \rfloor} f(d’)\\
& = \sum_{d|n} f(d) \sum_{d’|\lfloor \frac{n}{d} \rfloor} \mu(d’)\\
&= f(n)
\end{align}
$$

0x10 P3327 [SDOI2015]约数个数和

0x11 题目

题目链接: https://www.luogu.org/problemnew/show/P3327

题目需要让我们求 $\sum_{i=1}^n \sum_{j=1}^m d(i j)​$

其中$ d ​$为约数个数

显然$d(ij) = \sum_{x|i} \sum_{y|j} [gcd(x,y) = 1]​$

给一个比较感性的解释

显然,这个式子在枚举$i,j​$ 的所有约数的 gcd

如果 两个约数 gcd 为 1,则代表在 $i \cdot j$ 里可以有一个新的约数$x \cdot y$,否则表明这个已经被之前枚举过了

然后开始我们的魔法
$$
\sum_{i=1}^n \sum_{j=1}^m d(i j) \\
\begin{align}
& = \sum_{i=1}^n \sum_{j=1}^m \sum_{x|i} \sum_{y|j} [gcd(x,y) = 1]\\
& = \sum_{i=1}^n \sum_{j=1}^m \lfloor \frac{n}{i} \rfloor \lfloor \frac{m}{j} \rfloor [gcd(i,j) = 1] \\
\end{align}
$$
接下来的魔法需要莫比乌斯反演


$$
f(x) = \sum_{i=1}^n \sum_{j=1}^m \lfloor \frac{n}{i} \rfloor \lfloor \frac{m}{j} \rfloor [gcd(i,j) = x]\\
g(x) = \sum_{d|x} f(d)
$$
显然
$$
g(x) = \sum_{i=1}^n \sum_{j=1}^m \lfloor \frac{n}{i} \rfloor \lfloor \frac{m}{j} \rfloor [x|gcd(i,j)]
$$
然后提一下$x​$
$$
g(x) = \sum_{i=1}^{\lfloor \frac{n}{x} \rfloor} \sum_{j=1}^{\lfloor \frac{m}{x} \rfloor} \lfloor \frac{n}{xi} \rfloor \lfloor \frac{m}{xj} \rfloor
$$

然后因为我们求的是$f(1)$,则
$$
\begin{align}
f(1)& = \sum_{1|d} \mu(\frac{d}{1}) g(d)\\
& = \sum_{d=1}^{n} \mu(d) g(d)\\
&= \sum_{d=1}^{n} \mu(d) \sum_{i=1}^{\lfloor \frac{n}{x} \rfloor} \sum_{j=1}^{\lfloor \frac{m}{x} \rfloor} \lfloor \frac{n}{xi} \rfloor \lfloor \frac{m}{xj} \rfloor
\end{align}
$$
后面的显然可以数论分块,先预处理出$\sum_{i=1}^{n}\lfloor \frac{n}{i} \rfloor​$即可

显然,$\sum_{i=1}^{n}\lfloor \frac{n}{i} \rfloor​$ 和 $\sum_{i=1}^n d(i)​$ 是一样的

然而$d(i)$是积性函数,一起放给线性筛处理即可

0x02 代码

#include <cstdio>

inline int Min(int a, int b){return a < b? a: b;}

const int N = 51000;

int _case;
int n, m, pcnt;
long long ans;
bool vis[N + 20];
int pri[N + 20], mu[N + 20], smu[N + 20], d[N + 20], sum_d[N + 20], Min_pri[N + 20];

void pre(){
    mu[1] = 1; d[1] = 1; vis[1] = true;
    for(int i = 2; i <= N; i++){
        if( vis[i] == false ){
            pri[ ++pcnt ] = i;
            d[i] = 2; Min_pri[i] = 1;
            mu[i] = -1;
        }
        for(int j = 1; j <= pcnt; j++){
            if( pri[j] * i > N )
                break;
            vis[ pri[j] * i ] = true;
            if( i % pri[j] == 0 ){
                mu[ pri[j] * i ] = 0;    
                d[ pri[j] * i ] = d[i] / (Min_pri[i] + 1) * (Min_pri[i] + 2);
                Min_pri[ pri[j] * i ] = Min_pri[i] + 1;
            }
            else {
                mu[ pri[j] * i ] = -mu[i];
                d[ pri[j] * i ] = d[i] << 1;
                Min_pri[ pri[j] * i ] = 1;
            }
        }
    }
    for(int i = 1; i <= N; i++)
        smu[i] = smu[i - 1] + mu[i];
    for(int i = 1; i <= N; i++)
        sum_d[i] = sum_d[i - 1] + d[i];
}

int main(){
    pre();
    scanf("%d", &_case);
    while( _case -- ){
        scanf("%d%d", &n, &m);
        ans = 0;
        if(n > m) {int tmp = m; m = n; n = tmp;}
        for(int left = 1, rig; left <= n; left = rig + 1){
            rig = Min(n / (n / left), m / (m / left));    
            ans += 1LL * (smu[rig] - smu[left - 1]) * sum_d[n / left] * sum_d[m / left];
        }
        printf("%lld\n", ans);
    }
}

Luogu P1403 [AHOI2005]约数研究

题目

实际上一开始看到题目我是打算写筛法的,后来发先是个除法分块秒切题目

显然

$$
\sum_{i = 1}^{n} f(i) = \sum_{i = 1}^{n} \lfloor \frac{i}{n} \rfloor
$$

除法分块即可

代码

#include <cstdio>

int n;
long long ans;

int main(){
    scanf("%d", &n);
    for(int left = 1, rig; left <= n; left = rig + 1){
        rig  = n / (n / left);
        ans += 1LL * (n / left) * (rig - left + 1);
    }
    printf("%lld", ans);
}


[CQOI2007]余数求和

题目

一开始老是把思路和做​ 和​ 弄混……

然后陷入思维死角,感谢题解把我拉了出来​

$$
\sum_{i = 1}^{n} k \mod i \
= \sum_{i=1}^{n} k – i \times \lfloor \frac{k}{i} \rfloor \
= n \times k – \sum_{i=1}^{n} i \times \lfloor \frac{k}{i} \rfloor​
$$

然后就是除法分块,这个转化是真的要命QAQ

代码

#include <cstdio>

long long n,k,ans=0;

inline long long Min(long long a,long long b){return a<b?a:b;}

int main(){
   scanf("%lld%lld",&n,&k);
   ans=n*k;
   for(long long l=1,r;l<=n;l=r+1){
       if(k/l!=0) r=Min(k/(k/l),n);
       else r=n;
       ans-=(k/l)*(r-l+1)*(l+r)/2;
  }
   printf("%lld",ans);
}

Luogu 1248 加工生产调度

题目

题目链接: https://www.luogu.org/problemnew/show/P1248

这道题目暴力显然不可取$ 1000!$ , DP 的状态不太容易定义,只能贪心大法

分析贪心

我们可以先假设只有两个物品$ A $ 和$ B $,则存在以下两种情况($ x $代表物品在 A 工厂加工所需时间,$ y $ 代表 B 工厂,以下同理)

  1. 先加工$ A $后加工$ B $,则时间为$ A.x + B.y + \max(A.y, B.x) $
  2. 先加工$ B $后加工$ A $,则时间为$ A.y + B.x + \max(A.x, B.y) $

则若先加工 A 比先加工 B 优,应满足
$$
A.x + B.y + \max(A.y, B.x) < A.y + B.x + \max(A.x, B.y) \\
\begin{align}
&= \max(A.y, B.x) – A.y – B.x < \max(A.x, B.y) – A.y – B.x \\
&= -\min(A.y, B.x) < -\min(A.x, B.y) \\
&= \min(A.x, B.y) < \min(A.y, B.x)
\end{align}
$$
然后放到cmp里,A 了?

等一下,这个方法……是错误的

为什么错误?

首先,我们需要知道在 C++ 标准库中,sort 要求排序运算符必须保证严格排序

什么是严格排序

若一个比较运算符,则应当满足

非自反性,非对称性,传递性,不可比性的传递性 这四点

好消息是,我们的函数并不具有不可比性的传递性

即我们的排序函数可能汇出形如下面的状况:

$ x < y $ 且 $y < z$ 但$ x \geq z$

这样一来,我们排出来的序列就不满足之前的要求了

排序的改进

所以我们的目的很简单,使我们的排序符号变成保证严格排序

我们可以把相等的情况单独提出来讨论

  1. $A.x < B.x$ 且 $A.y < B.y$ 按 $A.x < A.y$ 排序
  2. $A.x = B.x $ 且 $A.y = B.y $ 不做处理
  3. $ A.x > B.x $ 且 $A.y > B.y$ 按 $B.y < B.x$ 排序

多关键字排序即可

代码

#include <cstdio>
#include <algorithm>

inline int Max(int a, int b){return a > b? a: b;}

const int N = 1100;

int n;
struct node{
    int x, y, d, id;
}a[N];

bool cmp(node x, node y){
    if(x.d == y.d){
        if(x.d <= 0) 
            return x.x < y.x;
        else 
            return x.y > y.y;
    }
    else return x.d < y.d;
}

int main(){
#ifdef woshiluo
    freopen("luogu.1248.in", "r", stdin);
    freopen("luogu.1248.out", "w", stdout);
#endif
    scanf("%d", &n);
    for(int i = 1; i <= n; i++){ scanf("%d", &a[i].x); }
    for(int i = 1; i <= n; i++){ 
        scanf("%d", &a[i].y); 
        a[i].id = i;
        a[i].d = (a[i].x == a[i].y ? 0: (a[i].x > a[i].y ? 1: - 1));
    }

    std::sort(a + 1, a + n + 1, cmp);

    int time_a = a[1].x, time_b = a[1].x + a[1].y;
    for(int i = 2; i <= n; i++){
        time_a += a[i].x;
        time_b = Max(time_a, time_b) + a[i].y;
    }
    printf("%d\n", time_b);
    for(int i = 1; i <= n; i++){ printf("%d ", a[i].id); }
}

后记

对与邻项排序的各种各样问题 ouuan 大佬的博客里有一篇讲的十分清楚 https://ouuan.github.io/%E6%B5%85%E8%B0%88%E9%82%BB%E9%A1%B9%E4%BA%A4%E6%8D%A2%E6%8E%92%E5%BA%8F%E7%9A%84%E5%BA%94%E7%94%A8%E4%BB%A5%E5%8F%8A%E9%9C%80%E8%A6%81%E6%B3%A8%E6%84%8F%E7%9A%84%E9%97%AE%E9%A2%98/

博主也是根据这篇博文学的,%%% ouuan

Luogu P1471 方差

我都快忘记还有一种公式叫完全平方公式了…

$$ (a + b) ^ 2 = a ^ 2 + b ^ 2 + 2ab $$

对,就是它,这道题目的核心

0x01 分析题目

题目链接: [https://www.luogu.org/problemnew/show/P1471]

让我们看看什么是方差

对于一个长度为 $ n $ 的数列 $ A $ ,其方差为
$$ \frac{1}{n} \sum^{n}_{i=1}(A_i – k)^2 $$
其中,$ k $ 为 $ A $ 数列的平均数

emmm, 看起来好复杂啊,等一下,这里有个平方

$$
\begin{align}
\frac{1}{n} \sum^{n}_{i=1}(A_i – k)^2 &= \frac{1}{n} \sum^{n}_{i=1}{A_i ^ 2 + k ^ 2 + 2A_ik} \\
&= \frac{1}{n} [\sum^{n}_{i=1}{A_i ^ 2} + nk^2 + 2k(\sum^{n}_{i=1}A_i)]
\end{align}
$$

恩,舒服多了

看一下操作,区间加,区间平均数(区间求和),区间方差?

区间方差怎么求?

先回到刚才的式子

然后,$ nk^2 $可以直接线段树求和处理出来
$2k(\sum^{n}_{i=1}A_i)$ 同上

$ \sum^{n}_{i=1}{A_i ^ 2} $怎么办?

回到完全平方公式,我们求平方和的时候,给数字$A_i$加$x$相当于把这个原$A_i^2$加上 $ (A_i + x) ^ 2 – A_i ^ 2 = (A_i ^ 2 + x ^ 2 + 2A_ix) – A_i ^ 2 = x ^ 2 + 2A_ix $

这不归根结底还是求和线段树吗……直接双标记线段树即可

0x02 代码

#include <cstdio>
const int N = 110000;
int op, n, m, x, y;
double k, z, tmp_double;
struct node{
    double sum, square;
    void operator+=(const node &b){
        this -> sum += b.sum;
        this -> square += b.square;
    }
}tmp;
node tree[N << 2];
double lazy[N << 2];
inline void pushup(int now){
    tree[now].sum = tree[now << 1].sum + tree[now << 1 | 1].sum;
    tree[now].square = tree[now << 1].square + tree[now << 1 | 1].square;
}
inline void pushdown(int now, int lson, int rson){
    if(lazy[now]){
        tree[now << 1].square  = tree[now << 1].square + (lazy[now] * tree[now << 1].sum * 2.0) + (lazy[now] * lazy[now]) * lson;
        tree[now << 1].sum += lazy[now] * lson;
        tree[now << 1 | 1].square  = tree[now << 1 | 1].square + (lazy[now] * tree[now << 1 | 1].sum * 2.0) + (lazy[now] * lazy[now]) * rson;
        tree[now << 1 | 1].sum += lazy[now] * rson;
        lazy[now << 1] += lazy[now];
        lazy[now << 1 | 1] += lazy[now];
        lazy[now] = 0;
    }
}
inline void query_add(int now, int left, int rig, int from, int to, double val){
    if(from <= left && rig <= to){
        tree[now].square = tree[now].square + ((val * tree[now].sum ) * 2) + (val * val) * (rig - left + 1);
        tree[now].sum = tree[now].sum + val * (rig - left + 1);
        lazy[now] += val;
        return ;
    }
    int mid = (left + rig) >> 1;
    pushdown(now, mid - left + 1, rig - mid);
    if(from <= mid) query_add(now << 1, left, mid, from, to, val);
    if(to > mid) query_add(now << 1 | 1, mid + 1, rig, from, to, val);
    pushup(now);
}
inline node query_sum(int now, int left, int rig, int from, int to){
    if(from <= left && rig <= to){
        return tree[now];
    }
    int mid = (left + rig) >> 1; node res = (node){0, 0};
    pushdown(now, mid - left + 1, rig - mid);
    if(from <= mid) res += query_sum(now << 1, left, mid, from, to);
    if(to > mid) res += query_sum(now << 1 | 1, mid + 1, rig, from, to);
    return res;
}
int main(){
#ifdef woshiluo
    freopen("luogu.1471.in", "r", stdin);
    freopen("luogu.1471.out", "w", stdout);
#endif
    scanf("%d%d", &n, &m);
    for(int i = 1; i <= n; i++){
        scanf("%lf", &z);
        query_add(1, 1, n, i, i, z);
    }
    while(m--){
        scanf("%d", &op);
        if(op == 1){
            scanf("%d%d%lf", &x, &y, &z);
            k += (y - x + 1) * z;
            query_add(1, 1, n, x, y, z);
        }
        else if(op == 2){
            scanf("%d%d", &x, &y);
            tmp = query_sum(1, 1, n, x, y);
            printf("%.4lf\n", tmp.sum / (y - x + 1));
        }
        else if(op == 3){
            scanf("%d%d", &x, &y);
            tmp = query_sum(1, 1, n, x, y);
            printf("%.4lf\n", ( (tmp.square + (tmp.sum / (y - x + 1)) * tmp.sum - ( 2 * tmp.sum * (tmp.sum / (y - x + 1)) ) ) / (y - x + 1) ) );
        }
    }
}
]]>

Educational Codeforces Round 60 (Rated for Div. 2) 解题报告

A. Best Subsegment

给定一串数列, 求区间最大平均值

显然,最大的区间平均值可以就是数列中最大的那个,然后基于最大的点想两边枚举区间长度

#include <cstdio>
#include <algorithm>
inline int Max(int a, int b){return a > b? a: b;}
const int N = 1e5 + 1e4;
int n, max, max_id, len, rig, max_len;
int a[N];
int main(){
    scanf("%d", &n);
    for(int i = 1; i <= n; i++){
        scanf("%d", &a[i]);
        if(a[i] > max)
            max = a[i];
    }
    for(int i = 1; i <= n; i++){
        if(a[i] == max && ( a[i - 1] != max || i == 1 ) ){
            rig = i; len = 0;
            while( a[rig] == a[i] && rig <= n ){
                rig++; len++;
            }
            max_len = Max(len, max_len);
        }
    }
    printf("%d", max_len);
}

B. Emotes

总共 $n$ 个表情, 你需要输出一个长度为 $m$ 的长度序列,但是同一个表情不能连续出现 $k$ 次

排个序,每个周期按,最大的出现 $k$ 次,然后出现一次第二大的,枚举周期即可

#include <cstdio>
#include <algorithm>
const int N = 2e5 + 1e4;
int n, m, k, tmp;
long long ans;
int a[N];
bool cmp(int a, int b){
    return a > b;
}
int main(){
    scanf("%d%d%d", &n, &m, &k);
    for(int i = 1; i <= n; i++){
        scanf("%d", &a[i]);
    }
    std::sort(a + 1, a + n + 1, cmp);
    if( m <= k ){
        printf("%lld", 1LL * a[1] * m);
        return 0;
    }
    else {
        ans = 1LL * a[1] * k + (long long) a[2];
        ans *= (m / (k + 1) );
        m %= (k + 1);
        ans += 1LL * a[1] * m;
        printf("%lld", ans);
    }
}

C. Magic Ship

本来是到直接两点坐标差的绝对值之和就可以的

然后有洋流

我们设每一轮洋流为一个周期,然后二分最少需要多少周期

然后 $O(n)$ 判断不在整数周期内部的就可以

#include <cstdio>
inline long long Aabs(long long a){return a < 0? (0 - a): a;}
const long long N = 1e5 + 1e4;
long long x1, y1, x2, y2, dx, dy, n, tmp, ans, mid;
char str[N];
int main(){
    scanf("%lld%lld", &x1, &y1);
    scanf("%lld%lld", &x2, &y2);
    scanf("%lld", &n);
    scanf("%s", str + 1);
    for(long long i = 1; i <= n; i++){
        if(str[i] == 'U') dy++;
        if(str[i] == 'D') dy--;
        if(str[i] == 'L') dx--;
        if(str[i] == 'R') dx++;
    }
    long long left = 0, rig = (int)2e9;
    while(left <= rig){
        mid = (left + rig) >> 1;
        if( Aabs(x1 + dx * mid - x2) + Aabs(y1 + dy * mid - y2) <= n * mid ) rig = mid - 1;
        else left = mid + 1;
    }
    x1 += dx * rig;
    y1 += dy * rig;
    ans = n * rig;
    if( Aabs(x1 - x2) + Aabs(y1 - y2) <= ans){
        printf("%lld", ans);
        return 0;
    }
    for(long long i = 1; i <= n; i++){
        if(str[i] == 'U') y1++;
        if(str[i] == 'D') y1--;
        if(str[i] == 'L') x1--;
        if(str[i] == 'R') x1++;
        ans++;
        if( Aabs(x1 - x2) + Aabs(y1 - y2) <= ans ){
            printf("%lld", ans);
            return 0;
        }
    }
    printf("-1");
}


]]>

Codeforces Round #538 (Div. 2) 解题报告

A Got Any Grapes?

顺序判断即

比赛的时候忘记写 else 既然 PP 了,然后就 FST

#include <cstdio>
int an,dm,mi;
int gr,pu,bl;
int main(){
    scanf("%d%d%d", &an, &dm, &mi);
    scanf("%d%d%d", &gr, &pu, &bl);
    if(an > gr) {
        printf("NO\n");
        return 0;
    }
    else gr -= an;
    if(pu + gr < dm){
        printf("NO\n");
        return 0;
    }
    else {
        if(pu <= dm) {dm -= pu; pu = 0;}
        else {pu -= dm; dm = 0;}
        if(gr < dm) {
            printf("NO\n");
            return 0;
        }
        else gr -= dm;
    }
    if(gr + pu + bl < mi) printf("NO\n");
    else printf("YES\n");
}

B Yet Another Array Partitioning Task

反正最后还是要前$ m \times k $ 个数字

直接离散化,然后见当前区间有$ m $ 个就收

#include <cstdio>
#include <algorithm>
const long long N = 2e5 + 1e4;
struct node{
    long long now, id;
}b[N];
long long n, m, k, cnt, time;
long long a[N];
bool vis[N];
bool cmp(node a, node b){
    return a.now > b.now;
}
int main(){
    scanf("%lld%lld%lld", &n, &m, &k);
    for(long long i = 1; i <= n; i++){
        scanf("%lld", &a[i]);
        b[i].now = a[i], b[i].id = i;
    }
    std::sort(b + 1, b + n + 1, cmp);
    long long tmp = m * k;
    for(long long i = 1; i <= tmp; i++) vis[ b[i].id ] = true, cnt += b[i].now;
    printf("%lld\n",cnt);
    cnt = 0;
    for(long long i = 1; i <= n; i++){
        cnt += vis[i];
        if(cnt == m){
            time++;
            time < k && printf("%lld ",i);
            cnt = 0;
        }
    }
}

C. Trailing Loves (or L’oeufs?)

求 $ b $ 进制下 $ n! $ 的末尾的 0 的个数

显然我们要求中间乘出来有多少个$b$

接下来就是分解质因数,然后枚举求最小值即可

#include <cstdio>
#include <cmath>
const long long INF = 1e18 + 1e17;
inline long long Min(long long a, long long b){return a < b? a: b;}
struct node{
    long long now, cnt;
}pri[(int)(1e6)];
int pcnt;
long long n, b, tmp, cnt, ans = INF;
inline void get_pri(long long b){
    tmp = std::sqrt(b);
    for(long long i = 2; i <= tmp; i++){
        if(b % i == 0){
            pri[ ++pcnt ] = (node){i, 0};
            while(b % i == 0){
                pri[pcnt].cnt++;
                b/=i;
            }
        }
    }
    if(b != 1) pri[ ++pcnt ] = (node){b, 1};
}
int main(){
    scanf("%lld%lld", &n, &b);
    get_pri(b);
    for(int i = 1; i <= pcnt; i++){
        tmp = 1; cnt = 0;
        while(tmp * pri[i].now <= n){
            tmp *= pri[i].now;
            if(tmp < 0 || tmp % pri[i].now != 0) break;
            cnt += n / tmp;
        }
        ans = Min(ans, cnt/pri[i].cnt);
    }
    printf("%lld\n", ans);
}

D. Flood Fil

有一个非常显然的地方,就是我们每次有两个状态,当前联通部分向左对其或者想右对其

然后我们先把数列中的联通部分预处理出来,然后区间 dp 即可

#include <cstdio>
#include <cstring>
inline int Min(int a, int b){return a < b? a: b;}
const int N = 5100;
int n;
int a[N], l[N], r[N], f[N][N];
int main(){
    scanf("%d", &n);
    for(int i = 1; i <= n; i++) scanf("%d", &a[i]);
    l[1] = 1;
    for(int i = 2; i <= n; i++){
        if(a[i] == a[i - 1]) l[i] = l[i - 1];
        else l[i] = i;
    }
    r[n] = n;
    for(int i = n - 1; i >= 1; i--){
        if(a[i] == a[i + 1]) r[i] = r[i + 1];
        else r[i] = i;
    }
    memset(f, 0x3f, sizeof(f));
    for(int i = 1; i<= n;i ++) f[ l[i] ][ r[i] ] = 0;
    for(int len = 0; len < n; len++){
        for(int i = 1,j; i + len <= n; i++){
            j = i + len;
            if(i > 1 && j < n && a[i - 1] == a[j + 1])
                f[ l[i - 1] ][ r[j + 1] ] = Min(f[ l[i - 1] ][ r[j + 1] ],f[i][j] + 1);
            if(i > 1)
                f[ l[i - 1] ][j] = Min(f[ l[i - 1] ][j], f[i][j] + 1);
            if(j < n)
                f[i][ r[j + 1] ] = Min(f[i][ r[j + 1] ], f[i][j] + 1);
        }
    }
    printf("%d", f[1][n]);
}

E. Arithmetic Progression

交互题目,60次询问内知道当前乱序等差数列的首项和公差

先二分找最大的,然后random_shuffle随机化询问,求与最大项差的 gcd 即为公差

#include <cstdio>
#include <algorithm>
int gcd(int a, int b){return b? gcd(b ,a%b): a;}
inline int Min(int a, int b){return a < b? a: b;}
inline int Aabs(int a){return a < 0? (0 - a): a;}
const int N = 1e6 + 1e5;
int n, global_tmp, d, max, chance_cnt = 60;
int a[N];
inline bool has_num(int now){
    printf("> %d\n",now);
    fflush(stdout);
    scanf("%d", &global_tmp);
    chance_cnt--;
    return global_tmp;
}
inline int binary_search_max(int max_limit){
    int left = 0, rig = max_limit, mid, res;
    while(left <= rig){
        mid = (left + rig) >> 1;
        if( has_num(mid) ) left = mid + 1;
        else rig = mid - 1, res = mid;
    }
    return res;
}
int main(){
    scanf("%d", &n);
    max = binary_search_max(1e9);
    for(int i = 1; i <= n; i++) a[i] = i;
    for(int i = 1; i <= 5; i++) std::random_shuffle(a + 1, a + n + 1);
    for(int i = 1; i <= Min(n ,chance_cnt); i++){
        printf("? %d\n", a[i]);
        fflush(stdout);
        scanf("%d", &global_tmp);
        if(global_tmp == max) continue;
        if(d == 0) d = Aabs(max-global_tmp);
        d = gcd(d , Aabs(max-global_tmp));
    }
    fflush(stdout);
    printf("! %d %d\n", max - (n - 1) * d, d);
}

F. Please, another Queries on Array?

这个题目首先得知道$\varphi(n)$的求法

然后就是乘积线段树和或线段树维护一下

就没有然后了

#include <cstdio>
const int N = 4e5 + 1e4;
const long long mod = 1e9 + 7;
int n, q, x, y, z, pcnt;
int a[N];
long long p[310], pri[N], inv[310];
char op[10];
inline long long ksm(long long a, long long p){
    long long res = 1;
    while(p){
        if(p & 1) res = (res * a) % mod;
        a = (a * a) % mod;
        p >>= 1;
    }
    return res;
}
inline void prime(){
    for(int i = 2; i <= 300; i++){
        for(int j = 1; j <= pcnt; j++)
            if(i % p[j] == 0) pri[i] |= pri[ p[j] ];
        if(pri[i] == 0) {p[++pcnt] = i; pri[i] = (1LL << (pcnt - 1LL));}
    }
}
inline void get_inv(){
    for(int i = 1; i <= pcnt; i++)
        inv[i] = ksm(p[i], mod - 2);
}
// Segment Tree Start
struct node{
    long long val,pri;
    void operator+= (const node &b){
        this -> val = (this -> val * b.val) % mod;
        pri |= b.pri;
    }
}tree[N << 2], lazy[N << 2];
inline void pushup(int now){
    tree[now].val = (tree[now << 1].val * tree[now << 1 | 1].val) % mod;
    tree[now].pri = tree[now << 1].pri | tree[now << 1 | 1].pri;
}
inline void pushdown(int now, int lson, int rson){
    if(lazy[now].pri != 0){
        tree[now << 1].val  =  (1LL * tree[now << 1].val * ksm(lazy[now].val, lson)) % mod;
        tree[now << 1 | 1].val  =  (1LL * tree[now << 1 | 1].val * ksm(lazy[now].val, rson)) % mod;
        tree[now << 1].pri |= lazy[now].pri;
        tree[now << 1 | 1].pri |= lazy[now].pri;
        lazy[now << 1].val = (1LL * lazy[now << 1].val * lazy[now].val) % mod;
        lazy[now << 1 | 1].val = (1LL * lazy[now << 1 | 1].val * lazy[now].val) % mod;
        lazy[now << 1].pri |= lazy[now].pri;
        lazy[now << 1 | 1].pri |= lazy[now].pri;
        lazy[now].val = 1LL; lazy[now].pri = 0;
    }
}
inline void query_mut(int now, int left, int rig, int from, int to, int val){
    if(from <= left && rig <= to){
        tree[now].val = (1LL * tree[now].val * ksm(val, (rig - left + 1LL))) % mod;
        tree[now].pri |= pri[val];
        lazy[now].val = (lazy[now].val * val) % mod;
        lazy[now].pri |= pri[val];
        return ;
    }
    int mid = (left + rig) >> 1;
    pushdown(now, mid - left + 1LL, rig - mid);
    if(from <= mid) query_mut(now << 1, left, mid, from, to, val);
    if(to > mid) query_mut(now << 1 | 1, mid + 1, rig, from, to, val);
    pushup(now);
}
inline node query_sum(int now, int left, int rig, int from, int to){
    if(from <= left && rig <= to) return tree[now];
    int mid = (left + rig) >> 1;
    node res = (node){1, 0};
    pushdown(now, mid - left + 1LL, rig - mid);
    if(from <= mid) res += query_sum(now << 1, left, mid, from, to);
    if(to > mid) res += query_sum(now << 1 | 1, mid + 1, rig, from, to);
    pushup(now);
    return res;
}
inline void build_tree(int now, int left, int rig){
    lazy[now].val = 1; lazy[now].pri = 0;
    if(left == rig){
        scanf("%lld", &tree[now].val);
        tree[now].pri = pri[ tree[now].val ];
        return ;
    }
    int mid = (left + rig) >> 1;
    build_tree(now << 1, left, mid);
    build_tree(now << 1 | 1, mid + 1, rig);
    pushup(now);
}
// Segment Tree End
int main(){
    prime();
    get_inv();
    scanf("%d%d", &n, &q);
    build_tree(1, 1, n);
    for(int i = 1; i <= n; i++) scanf("%d", &a[i]);
    while(q--){
        scanf("%s", op);
        if(op[0] == 'M'){
            scanf("%d%d%d", &x, &y, &z);
            query_mut(1, 1, n, x, y, z);
        }
        else {
            scanf("%d%d", &x, &y);
            node tmp = query_sum(1, 1, n, x, y);
            for(int i = 1; i <= pcnt; i++){
                if((tmp.pri >> (i - 1LL)) & 1LL)
                    tmp.val = (((tmp.val * (p[i] - 1) ) % mod) * inv[i]) % mod;
            }
            printf("%lld\n", tmp.val);
        }
    }
}


]]>

CCF WC 2019 游记

Day -1

上午起来收拾了一下就前往乌鲁木齐市机场了

在机场里面互相定位是一件困难的事情,我们最终通过奇迹淫巧和瞎挥手聚在了一起

然后是漫长的安检和候机….

在经历各种各样奇怪的娱乐过后,飞机终于落地

飞机上拍的云朵

下去坐地铁,站了一个多小时后有疯狂转圈终于吃上了人生中第一顿麦当劳并到达了宾馆

发现电视有 HDML 口,接之

一顿瞎嗨,用电视看了看鬼畜和老番

Day 0

早上六点起来去萝岗吃早茶

买了一张地铁日卡,从 苏元 苟到坐到 萝岗

c0per 大爷在买日卡的时候刷人家的红包刷出来了 10 块钱…

内地人消费水平过于强大瑟瑟发抖…

精致却贵

平均每人30+却连五分饱都吃不上…

然后坐回 苏元 在宾馆休息一下就再度出发,站了20站到 公园前

动漫星城 get

都想买喵 买不起喵

咱是真的没想到有一整栋楼都是2.5次元类型的产品,这里是天堂吗

然后我们的钱包削减的速度也让咱心疼…

然后步行去animate

我错了,animate 让我感到钱包爆炸….

接着大家在麦当劳吃午饭(虽然我只吃了了冰激凌)

经过一阵时间的的休整,我们又站回了苏元

经过漫长的爬坡,终于到达了广州市第二中学

于宿舍拍摄

晚上就是开幕式,没有 dzd 讲话差评

广州二中真的是舒服啊qwq

Day 1 – 4 听课

在第一课堂颓废(伦敦真是太可爱了)

第一课堂完全不会,第二课堂没有意思

真实冬眠

实际上还是把生成函数的入门了解的差不多并补习许多数学知识

即得易见平凡,仿照上例显然,留作习题答案略,读者自证不难
反之亦然同理,推论自然成立,略去过程Q, E, D,由上可知证毕

对析合树有了一个初步的认识

早饭不是很热乎使人不太欢乐

但是整体而言伙食已经很优秀了

还托第一课堂的福知道一个及其虐人的gal — 交响乐之雨

Day 5 考试 && 文艺汇演

考试

T1 嗯,purfey序列?然后?数学?滚!type = 0 大发好

T2 嗯,这题面怎么跟个文开发档一样?subtask1 两行?subtask2 打表?subtask3 直接bfs?然后?不会

T3 博弈论?数学? …有点困…诶!怎么突然就过去 1.5h 了?我WC怎么睡着了???喵喵喵

考完感觉自己只写了暴力?Fe 预定?

文艺汇演

曾经我认为,学 Oi 的到后期都会逐渐自闭,但是通过这此文艺汇演,我发现,所有 Oier 的心都是连在一起的,大家都不会孤独的!

这次文艺汇演从无到有到开始只有40分钟左右的时间,可以说是效率非常的高了

上来先是各路dalao唱歌

然后 Oi改编 歌曲超级好评,每当听起这些歌曲的时候,心里都有一种认为自己要更努力的冲动

LCA 的女装相声超级好评,禁赛边缘

Day 5 科技馆 && 闭幕式

科技馆

不知道为什么会有社会活动这种奇怪的东西,不过路上景色十分优秀

珠江

科技馆没什么可以说的,每个城市都是一样的qwq

闭幕式

回来知道T2 Subtask2 重新测试,自己分数貌似是个Ag

但是看泄露名单并没有

然而我写的是正解喵?

并不知道发生了什么,所以就随便吧,反正也没指望拿牌子

dzd 终于讲话了,上来就怼教育部也是非常的优秀了,不过他句句在理,十分钦佩这样子的人

意料之外的拿上了Cu

对于一个初中生来说,这应该是一个不错的成绩吧?

不过要是以初中上要求自己的话,还是别努力了吧(况且好像还有好多初中选手吊打我)

所以下次一定要好好努力,起码,不睡觉!

Day 6 有缘再见

正如 wh 教授所言,美好的事物总是短暂的

我们乘坐早上 11:30 的巴士前往白云机场

临走前去找 ouuan 大佬 PY 到一个士力架

到达机场,过了安检在星巴克买了一杯红茶拿铁

果然到了内地就是要喝星巴克

各位 Oier 们, 有缘再见!

祝各位 RP++

]]>

Luogu P3381 最小费用最大流

费用流

最大流 && Dinic 算法 一文中,我们简述了网络流的概念及最大流的求法

但如果每条边不单单有流的限制,还有使用的费用,我们现在不光只求最大流,还要求费用最小

这就是最小费用最大流问题

最小费用最大流

怎么求?最大流我们可以用ek,然后……

我们要保证费用最小……于是我们就可以

bfs换成spfa !

是的你没有看错就这么暴力

别忘记了,反向边的作用的反悔药,所以我们应当把方向边的费用改成原来的相反数

然后跑就对了

Luogu P3381

#include <cstdio>
#include <cstring>
inline int Min(int a,int b){return a<b?a:b;}
const int N = 5100;
const int M = 51000;
const int INF = 0x7f7f7f7f;
int n, m, s, t, u, v, w, f;
// edge start
struct edgE{
    int to,next,flow,val;
}e[M<<1];
int ehead[N],ecnt;
inline void add_edge(int now, int to, int flow, int val){
    ecnt++;
    e[ ecnt ].to = to;
    e[ ecnt ].val = val;
    e[ ecnt ].flow = flow;
    e[ ecnt ].next = ehead[now];
    ehead[ now ] = ecnt;
}
// edge end
// MCMF start
int Maxflow, Mincost, head, tail, head_dep, tail_dep;
int flow[N], dis[N], cur[N], pre[N], q[N];
bool vis[N];
inline bool spfa(){
    memset(vis, false, sizeof(vis));
    memset(dis, 0x7f, sizeof(dis));
    head = tail = head_dep = tail_dep = 0;
    q[head] = s;vis[s] = true;dis[s] = 0;pre[t] = 0;flow[s] = INF;
    while(head <= tail || head_dep < tail_dep){
        for(int i = ehead[ q[head] ]; i; i = e[i].next){
            if(e[i].flow > 0 && dis[ q[head] ] + e[i].val < dis[ e[i].to ]){
                dis[ e[i].to ] = dis[ q[head] ] + e[i].val;
                cur[ e[i].to ] = i;
                pre[ e[i].to ] = q[head];
                flow[ e[i].to ] = Min(flow[ q[head] ], e[i].flow);
                if(!vis[ e[i].to ]){
                    vis[ e[i].to ] = true;
                    tail++;
                    if(tail > n) tail=1, tail_dep++;
                    q[ tail ] = e[i].to;
                }
            }
        }
        vis[ q[head] ]= false;
        head++;
        if(head > n) head = 1, head_dep++;
    }
    return pre[t] != 0;
}
inline void MCMF(){
    while( spfa() ){
        Maxflow += flow[t];
        Mincost += flow[t] * dis[t];
        int now = t;
        while(now != s){
            e[ cur[now] ].flow -= flow[t];
            e[ cur[now] + ( ( cur[now] & 1 )? 1: -1) ].flow += flow[t];
            now = pre[now];
        }
    }
}
// MCMF end
int main(){
    scanf("%d%d%d%d", &n, &m, &s, &t);
    for(int i = 1; i <= m; i++){
        scanf("%d%d%d%d", &u, &v, &w, &f);
        add_edge(u, v, w, f);
        add_edge(v, u, 0, -f);
    }
    MCMF();
    printf("%d %d", Maxflow, Mincost);
}
]]>

最大流 && Dinic 算法

最大流问题

给定一个有向图,有两个特别的点$S$和$T$,每条边都有一个容量限制 ,现询问,从 $ S $ 到 $ T $ 最大可以流多少过去

第一眼看过去,跟水厂子往水管里灌水一样简单

问题在于这个东西并不好计算

有没有办法计算它呢

如果我们每次都枚举我们怎么走的话,时间复杂度是不可估计的

我们有没有办法让我们走过去的流反悔呢?

有,这里就要引入两个概念: 残余网络 && 反向边

残余网络

对于每次我们跑出来的一条流过过后以每条边现能通过的最大流量所成的网络图叫残余网络

反向边

反向边是在给算法一个反悔的机会

如果原图中有一个从 $ u $ 到 $ v $ 的边,则从 $ v $ 到 $ u $ 的边就是这条边的反向边

反向边如何后悔呢?

EK 算法

我们在每一次流过去大小为 $ t $ 一条流时,除了将原有的边上流量减去$ t $之外,也需要将其反向边的流来加上$ t $

这样以来,我们就保证以下的两点

  • 残余网络中总流量和原图相同
  • 在残余网络中走出来的流一定存在

至于怎么解释,咱的看法是,写两边就明白了,看别人的解释之后越看越懵

我们一直寻找可行的流,直到无法流到$ t $为止,已经的流过去的和就一定是一个最大流

复杂度: $O(EK(n,m))=O(m^2n)=O(TLE)$

Dinic 算法

EK算法使用bfs寻找流,但bfs寻找流的一些特性决定其不容易实现当前弧优化,以及多路增广

而且EK算法还会走回路…

这种情况下我们又要引出一个新概念: 分层图

分层图

我们以$ S $为根,可以往下走的边流大于0,然后像树一样记录每个节点的深度(多次访问取第一次)

这样以来我们可以确定以下两点:

  • 避免了回路的可能
  • 如果存在可行流,则点$ T $深度一定存在

等一下,如果我们搜索的时候都限制只能往当前节点深度+1的点去搜了,我们是不是就可以dfs了?
是的qwq

不过我们先来看一下普通的Dinic吧[Luogu P3376]

#include <cstdio>
#include <cstring>
inline int Min(int a,int b){return a<b?a:b;}
const int N=11000;
const int M=110000;
int n,m,s,t,u,v,w;
// edge start
struct edge{
    int next,to,val,un;
}e[M<<1];
int ehead[N],ecnt;
inline void add_edge(int now,int to,int val,int un){
    ecnt++;
    e[ecnt].to=to;
    e[ecnt].val=val;
    e[ecnt].next=ehead[now];
    e[ecnt].un=ecnt+un;
    ehead[now]=ecnt;
}
// edge end
namespace Dinic{
    int depth[N];
    int q[N],head,tail;
    inline bool bfs(){
        memset(depth,0,sizeof(depth));
        head=tail=0;
        q[head]=s;
        depth[s]=1;
        while(head<=tail){
            for(int i=ehead[q[head]];i;i=e[i].next){
                if(e[i].val>0&&depth[e[i].to]==0){
                    q[++tail]=e[i].to;
                    depth[e[i].to]=depth[q[head]]+1;
                }
            }
            head++;
        }
        if(depth[t]==0) return false;
        else return true;
    }
    int dfs(int now,int dist){
        if(now==t) return dist;
        for(int i=ehead[now];i;i=e[i].next){
            if(depth[e[i].to]==depth[now]+1&&e[i].val!=0){
                int di=dfs(e[i].to,Min(dist,e[i].val));
                if(di>0){
                    e[i].val-=di;
                    e[e[i].un].val+=di;
                    return di;
                }
            }
        }
        return 0;
    }
    inline int dinic(){
        int ans=0;
        while(bfs()){
            while(int d=dfs(s,0x7fffffff)) ans+=d;
        }
        return ans;
    }
}
int main(){
    scanf("%d%d%d%d",&n,&m,&s,&t);
    for(int i=1;i<=m;i++){
        scanf("%d%d%d",&u,&v,&w);
        add_edge(u,v,w,1);
        add_edge(v,u,0,-1);
    }
    printf("%d",Dinic::dinic());
}

好了,我们来讲讲两个优化: 当前弧优化 && 多路增广

当前弧优化

在当前状况下已经有一次从$ u $到$ v $的增广,则不存在第二次从$ u $到$ v $的增广

故此在我们遍历边的时候,我们可以不从ehead[now]开始,而是从上一次遍历过的边的下一个开始

namespace Dinic{
    int depth[N],rnow[N];
    int q[N],head,tail;
    inline bool bfs(){
        memset(depth,0,sizeof(depth));
        head=tail=0;
        q[head]=s;
        depth[s]=1;
        while(head<=tail){
            for(int i=ehead[q[head]];i;i=e[i].next){
                if(e[i].val>0&&depth[e[i].to]==0){
                    q[++tail]=e[i].to;
                    depth[e[i].to]=depth[q[head]]+1;
                }
            }
            head++;
        }
        if(depth[t]==0) return false;
        else return true;
    }
    int dfs(int now,int dist){
        if(now==t) return dist;
        for(int& i=rnow[now];i;i=e[i].next){
            if(depth[e[i].to]==depth[now]+1&&e[i].val!=0){
                int di=dfs(e[i].to,Min(dist,e[i].val));
                if(di>0){
                    e[i].val-=di;
                    e[e[i].un].val+=di;
                    return di;
                }
            }
        }
        return 0;
    }
    inline int dinic(){
        int ans=0;
        while(bfs()){
            for(int i=1;i<=n;i++) rnow[i]=ehead[i];
            while(int d=dfs(s,0x7fffffff)) ans+=d;
        }
        return ans;
    }
}

多路增广

很多时候在dfs的时候,我们的dist还没用完就return了,很可惜了有没有

多路增广就可以避免这些问题,我们把剩下的dist,继续往下找!

#include <cstdio>
#include <cstring>
inline int Min(int a,int b){return a<b?a:b;}
const int N=11000;
const int M=110000;
int n,m,s,t,u,v,w;
// edge start
struct edge{
    int next,to,val,un;
}e[M<<1];
int ehead[N],ecnt;
inline void add_edge(int now,int to,int val,int un){
    ecnt++;
    e[ecnt].to=to;
    e[ecnt].val=val;
    e[ecnt].next=ehead[now];
    e[ecnt].un=ecnt+un;
    ehead[now]=ecnt;
}
// edge end
namespace Dinic{
    int depth[N],rnow[N];
    int q[N],head,tail;
    inline bool bfs(){
        memset(depth,0,sizeof(depth));
        head=tail=0;
        q[head]=s;
        depth[s]=1;
        while(head<=tail){
            for(int i=ehead[q[head]];i;i=e[i].next){
                if(e[i].val>0&&depth[e[i].to]==0){
                    q[++tail]=e[i].to;
                    depth[e[i].to]=depth[q[head]]+1;
                }
            }
            head++;
        }
        if(depth[t]==0) return false;
        else return true;
    }
    int dfs(int now,int dist){
        if(now==t) return dist;
        int fl=0,di;
        for(int& i=rnow[now];i;i=e[i].next){
            if(depth[e[i].to]==depth[now]+1&&e[i].val!=0){
                di=dfs(e[i].to,Min(dist,e[i].val));
                if(di>0){
                    e[i].val-=di;
                    e[e[i].un].val+=di;
                    fl+=di;
                    dist-=di;
                    if(!dist) break;
                }
            }
        }
        return fl;
    }
    inline int dinic(){
        int ans=0;
        while(bfs()){
            for(int i=1;i<=n;i++) rnow[i]=ehead[i];
            while(int d=dfs(s,0x7fffffff)) ans+=d;
        }
        return ans;
    }
}
int main(){
    scanf("%d%d%d%d",&n,&m,&s,&t);
    for(int i=1;i<=m;i++){
        scanf("%d%d%d",&u,&v,&w);
        add_edge(u,v,w,1);
        add_edge(v,u,0,-1);
    }
    printf("%d",Dinic::dinic());
}

就这些了,那么优化过后的复杂度是什么呢
$O(Dinic(n,m))=O(n^2m)=O(AC)$

还有一个ISAP算法,不在今天的讨论之内,毕竟复杂度和Dinic差别不大且常数略大

End

]]>