类欧几里得学习

类欧几里得算法学习

模型

类欧几里得是用来快速的求以下式子：

• $f(a,b,c,n)={\sum }_{i=0}^n[\frac{ai+b}{c}]$
• $g(a,b,c,n)={\sum }_{i=0}^{n}i\times [\frac{ai+b}{c}]$
• $h(a,b,c,n)={\sum }_{i=0}^n[\frac{ai+b}{c}{]}^2$

因其求解的方式与欧几里得算法类似，故被称为类欧几里得算法，同样时间复杂度为O(logn)。

计算 f(a,b,c,n)

当a>=c||b>=c

• ${\sum }_{i=0}^n[\frac{ai+b}{c}]={\sum }_{i=0}^n[\frac{a\%c\times i+b\%c}{c}]+[\frac{a}{c}]\sum i+[\frac{b}{c}]\sum$
• $f(a,b,c,n)=f(a\%c,b\%c,c,n)+[\frac{a}{c}](n\times (n+1)/2)+[\frac{b}{c}]\times (n+1)$

当a<c&&b<c

• 对于${\sum }_{i=0}^n[\frac{ai+b}{c}]$，其还有一个几何意思是：求一条直线下方的整点数量。
• 设$m=\lfloor \frac{an+b}{c}\rfloor -1$
• $f(a,b,c,n)$
  $={\sum }_{i=0}^n[\frac{ai+b}{c}]$
  $={\sum }_{i=0}^n{\sum }_{j=1}^{m+1}[j<=\lfloor \frac{ai+b}{c}\rfloor ]$
  $={\sum }_{j=1}^{m+1}{\sum }_{i=0}^n[j<=\lfloor \frac{ai+b}{c}\rfloor ]$
  $={\sum }_{j=0}^m{\sum }_{i=0}^n[j+1<=\lfloor \frac{ai+b}{c}\rfloor ]$
  $={\sum }_{j=0}^m{\sum }_{i=0}^n[j+1<=\frac{ai+b}{c}]$
  $={\sum }_{j=0}^m{\sum }_{i=0}^n[jc+c<=ai+b]$
  $={\sum }_{j=0}^m{\sum }_{i=0}^n[ai>jc+c-b-1]$
  $={\sum }_{j=0}^m{\sum }_{i=0}^n[i>\frac{jc+c-b-1}{a}]$
  $={\sum }_{j=0}^m(n-\lfloor \frac{jc+c-b-1}{a}\rfloor )$
  $={\sum }_{i=0}^{m}n-f(c,c-b-1,a,m)$

当n=0：$f(a,b,c,n)=[\frac{b}{c}]$
当a=0：$f(a,b,c,n)=(n+1)[\frac{b}{c}]$

代码如下：

ll f(ll a,ll b,ll c,ll n){
	if(n==0)return b/c;
	if(a==0)return b/c*(n+1);
	if(a>=c&&b>=c)return f(a%c,b%c,c,n)+a/c*(n*n+1)/2+b/c*(n+1);
    if(a<c||b<c)return (a*n+b)/c*n+f(c,c-b-1,a,(an+b)/c-1);
}

计算h(a,b,c,n)

当a>=c||b>=c

• $h(a,b,c,n)={\sum }_{i=0}^n[\frac{ai+b}{c}{]}^2$
• $={\sum }_{i=0}^n([\frac{a}{c}]\times i+[\frac{b}{c}]+[\frac{a\%c\times i+b\%c}{c}]{)}^2$
• 根据三个的完全平方数展开得到
• 原式$={\sum }_{i=0}^n([\frac{a}{c}]\times i{)}^2+{\sum }_{i=0}^n([\frac{b}{c}]{)}^2+{\sum }_{i=0}^n([\frac{a\%c\times i+b\%c}{c}]{)}^2+2[\frac{a}{c}][\frac{b}{c}]{\sum }_{i=0}^{n}i+2[\frac{a}{c}]{\sum }_{i=0}^{n}i[\frac{a\%c\times i+b\%c}{c}]+2[\frac{b}{c}]{\sum }_{i=0}^n[\frac{a\%c\times i+b\%c}{c}]$
• 化简得到$=[\frac{a}{c}{]}^2\times \frac{(n(n+1)(2n+1))}{6}+[\frac{b}{c}{]}^2\times (n+1)+h(a\%c,b\%c,c,n)+2[\frac{a}{c}][\frac{b}{c}]\times \frac{n(n+1)}{2}+2[\frac{a}{c}]g(a\%c,b\%c,c,n)+2[\frac{b}{c}]f(a\%c,b\%c,c,n)$

当a<c&&b<c

• 同理根据证明直接写出最后结论，设$m=[\frac{an+b}{c}]-1$
• $h(a,b,c,n)=n(m+1{)}^2-2g(c,c-b-1,a,m)-f(c,c-b-1,a,m)$

当n=0：$h(a,b,c,n)=[\frac{b}{c}{]}^2$
当a=0：$h(a,b,c,n)=(n+1)\times [\frac{b}{c}{]}^2$

计算g(a,b,c,n)

当a>=c||b>=c

• $g(a,b,c,n)={\sum }_{i=0}^{n}i[\frac{ai+b}{c}]$
• $=[\frac{a}{c}]\frac{n(n+1)(2n+1)}{6}+[\frac{b}{c}]\frac{n(n+1)}{2}+g(a\%c,b\%c,c,n)$

当a<c&&b<c

• 设$m=[\frac{an+b}{c}]-1$
• $g(a,b,c,n)=\frac{n(n+1)(m+1)-h(c,c-b-1,a,m)-f(c,c-b-1,a.m)}{2}$

当n=0：$f(a,b,c,n)=0$
当a=0：$f(a,b,c,n)=[\frac{b}{c}](\frac{n(n+1)}{2})$

计算三者

由于后两者均需要其他形式，故递归的时候可以一起算比较快
模板如下：

#include <bits/stdc++.h>
#define ll long long
using namespace std;
const ll P = 998244353;
ll i2 = 499122177, i6 = 166374059;
struct data {
	ll f, g, h;
}; 
data calc(ll n, ll a, ll b, ll c) {
	ll ac = a / c, bc = b / c, m = (a * n + b) / c, n1 = n + 1, n21 = n * 2 + 1;
	data ans;
	ans.f=ans.g=ans.h=0;
	if (a == 0) {  // 迭代到最底层
		ans.f = bc * n1 % P;
		ans.g = bc * n % P * n1 % P * i2 % P;
		ans.h = bc * bc % P * n1 % P;
		return ans;
	}
	if (a >= c || b >= c) {  // 取模
		ans.f = n * n1 % P * i2 % P * ac % P + bc * n1 % P;
		ans.g = ac * n % P * n1 % P * n21 % P * i6 % P + bc * n % P * n1 % P * i2 % P;
		ans.h = ac * ac % P * n % P * n1 % P * n21 % P * i6 % P +
		      bc * bc % P * n1 % P + ac * bc % P * n % P * n1 % P;
		ans.f %= P, ans.g %= P, ans.h %= P;

		data e = calc(n, a % c, b % c, c);  // 迭代

		ans.h += e.h + 2 * bc % P * e.f % P + 2 * ac % P * e.g % P;
		ans.g += e.g, ans.f += e.f;
		ans.f %= P, ans.g %= P, ans.h %= P;
		return ans;
	} else {
		data e = calc(m - 1, c, c - b - 1, a);
		ans.f = n * m % P - e.f, ans.f = (ans.f % P + P) % P;
		ans.g = m * n % P * n1 % P - e.h - e.f, ans.g = (ans.g * i2 % P + P) % P;
		ans.h = n * m % P * (m + 1) % P - 2 * e.g - 2 * e.f - ans.f;
		ans.h = (ans.h % P + P) % P;
		return ans;
	}
}
ll T, n, a, b, c;
int main() {
	scanf("%lld", &T);
	while (T--) {
		scanf("%lld%lld%lld%lld", &n, &a, &b, &c);
		data ans = calc(n, a, b, c);
		printf("%lld %lld %lld\n", ans.f, ans.h, ans.g);
	}
	return 0;
}

例题

题目1：

• 题目描述：求${\sum }_{i=0}^{HM-1}|i/HM-imodM/|<=A$