Interval GCD（线段树区间增加、区间查询）

题目描述
给定一个长度为N的数列A，以及M条指令 (N≤5*10^5, M<=10^5)，每条指令可能是以下两种之一：
“C l r d”，表示把 A[l],A[l+1],…,A[r] 都加上 d。
“Q l r”，表示询问 A[l],A[l+1],…,A[r] 的最大公约数(GCD)。
输入
第一行两个整数N,M，第二行N个整数Ai，接下来M行每条指令的格式如题目描述所示。
输出
对于每个询问，输出一个整数表示答案。
样例输入
5 5
1 3 5 7 9
Q 1 5
C 1 5 1
Q 1 5
C 3 3 6
Q 2 4
样例输出
1
2
4
算法提示：
这个题的关键是，如果执行“C l r d”后，[l，r]区间的每个数都都加上d后，然后再求最大公约数。有的学生总问我，直接一个for语句修改后，再用线段数组至下而上维护最大公约数不就可以了，这个问题说过好多次了，如果这样，时间复杂度是O（N）。
还有的同学这样想，既然用线段数组维护最大公约数，那么就在线段树上[l，r]区间每一个A[ i]都增加一个d，每找到a[i]就的用log（N），那么把这个区间都修改完就用（r-l+1） * log（N），时间复杂度N * log(N)。
如果多次询问 M * N * log（N），代价太大。
所以我们用树状数组求前缀和的功能维护每个数组元素的变化。不会的同学参考A Tiny Problem with Integers(树状数组区间增加、单点查询)这篇博文。
这样只需要用线段树维护，在l处增加d，在r+1处减去d，就可以了。
问题又产生了，因为区间修改，只是在l处增加d，在r+1处减去d，如果直接在线段树上查询最大公约数，当然是错误答案了。
这个时候gcd(x,y,z)=gcd(x,y-x,z-y)登场了。

我们知道gcd(x,y)=gcd(x,y-x)。它可以进一步扩展到三个数的情况：gcd(x,y,z)=gcd(x,y-x,z-y)。实际上，该性质对任意多个整数都成立。
再回顾上面“C l r d”这条命令的影响：
A[l]+d,A[l+1]+d,A[l+2]+d,A[l+3]+d,…,A[r-1]+d,A[r]+d。
求A[l，r]最大公约数，根据上面的公式，我们也可以求
gcd(A[l],A[L+1]-A[l],A[l+2]-A[l+1],A[l+3]-A[l+2],…,A[r]-A[r-1]。
这样一来，A[l，r]区间增加d后，两两相减，除了第一项A[l],其它向相减后增加的d，正好抵消掉。例如(A[r]+d)-（A[r-1]+d）=A[r]-A[r-1]。

这样求最大公约数就变为gcd（A[l]+变化量，A[l+1],A[l+2],…,A[r]）,除了第一项A[l]+变化量，其它都可以直接在线段树上查询ask(1,l+1,r)，变化量用树状数组前缀和可以查询到，然后求gcd（A[l]+变化量，ask（1，l+1,r））即可。

综上所述，我们可以构造一个长度为N的新数列B，其中
B[i]=a[i]-a[i-1],B[1]可为任意值。数列B称为A的差分序列。
用线段树维护序列B的区间最大公约数。

这样一来，询问“Q l r”，就等于求出gcd（A[l],ask(1,l+1,r)）。
在指令“C l r d”下，只有B[r+1]被减掉了d，所以在维护B的线段树上只进行两次单点修改即可。另外，询问时需要数列A中的值，可额外用一个支持“区间增加，单点查询”的树状数组对数列A进行维护。

完整代码：

#include<iostream>
#include<cstdio>
#include<algorithm> 
using namespace std;
const int SIZE = 500010;
struct SegmentTree {
	int l, r;
	long long dat;
} t[SIZE * 4];
long long a[SIZE], b[SIZE], c[SIZE];
int n, m;

long long gcd(long long a, long long b) {
	return b ? gcd(b, a%b) : a;
}

// 维护区间gcd的线段树 
void build(int p, int l, int r) {
	t[p].l = l, t[p].r = r;
	if (l==r) { t[p].dat = b[l]; return; }
	int mid = (l + r) / 2;
	build(p*2, l, mid);
	build(p*2+1, mid+1, r);
	t[p].dat = gcd(t[p*2].dat, t[p*2+1].dat);
}

void change(int p, int x, long long v) {
	if (t[p].l == t[p].r) { t[p].dat += v; return; }
	int mid = (t[p].l + t[p].r) / 2;
	if (x <= mid) change(p*2, x, v);
	else change(p*2+1, x, v);
	t[p].dat = gcd(t[p*2].dat, t[p*2+1].dat);
}

long long ask(int p, int l, int r) {
	if (l <= t[p].l && r >= t[p].r) return abs(t[p].dat);
	int mid = (t[p].l + t[p].r) / 2;
	long long val = 0;
	if (l<=mid) val = gcd(val, ask(p*2, l, r));
	if (r>mid) val = gcd(val, ask(p*2+1, l, r));
	return abs(val);
}

long long sum(int x) {
	long long y = 0;
	for (; x; x -= x & -x) y += c[x];
	return y;
}

void add(int x, long long y) {
	for (; x <= n; x += x & -x) c[x] += y;
}

int main() {
	cin >> n >> m;
	for (int i = 1; i <= n; i++) {
		scanf("%lld", &a[i]);
		b[i] = a[i] - a[i-1];
	}
	build(1, 1, n);
	for (int i = 1; i <= m; i++) {
		char str[2]; scanf("%s", str);
		int l, r; scanf("%d%d", &l, &r);
		if (str[0] == 'Q') {
			long long al = a[l] + sum(l);
			long long val = l < r ? ask(1, l + 1, r) : 0;
			printf("%lld\n", gcd(al, val));
		} else {
			long long delta; scanf("%lld", &delta);
			change(1, l, delta);
			if (r < n) change(1, r + 1, -delta);
			add(l, delta);
			add(r + 1, -delta);
		}
	}
}