【数学】[Baltic2016][BZOJ5184] Spiral

分析:

恶心至极。。。
首先,把图拆成一个个的环,
第一层是:1
第二层是:2,3,4,5,6,7,8,9
……
然后,就可以发现,单独的一个环的贡献,可以在O(1)复杂度内算出来。(拆成四条边)

然而这要T

问题就在于如何快速求多个环的贡献。

首先,我们可以把原矩阵也拆分了,拆分成数个部分环的组合:
在这里插入图片描述
这样一来,每个部分环中,被剖分的矩形只有两种可能:
在这里插入图片描述
在这里插入图片描述
可以通过推倒证明,在这两种情况下,每一个环的贡献一定满足:
S i = A ∗ i 3 + B ∗ i 2 + C ∗ i 1 + D ∗ i 0 Si=A*i^3+B*i^2+C*i^1+D*i^0 Si=Ai3+Bi2+Ci1+Di0
(其中i表示是从上一层开始的第i个环)

所以,要求K个环的贡献和,就只需要算:
A ∗ ( 0 3 + 1 3 + 2 3 + … … + ( k − 1 ) 3 ) + B ∗ ( 0 2 + 1 2 + … … + ( j − 1 ) 2 ) + C ∗ ( 0 1 + 1 1 + … … + ( k − 1 ) 1 ) + D ∗ ( 0 0 + 1 0 + 2 0 + … … + ( k − 1 ) 0 ) A*(0^3+1^3+2^3+……+(k-1)^3)+B*(0^2+1^2+……+(j-1)^2)+C*(0^1+1^1+……+(k-1)^1)+D*(0^0+1^0+2^0+……+(k-1)^0) A(03+13+23++k1)3)+B(02+12++(j1)2)+C(01+11++k1)1)+D(00+10+20++(k1)0)
这几个都是有公式的,自己去套就好了。

然后,从证明过程中可以看出,推出A,B,C,D是几乎不可能的。但有一个更好的方法:我们可以暴力算出前4个环的贡献,然后将A,B,C,D反解出来。

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<vector>
#define SF scanf
#define PF printf
#define MAXN 100010
#define MOD 1000000007
#define y0 warspitewarspitewarspitewarspitewarspitewarspite
using namespace std;
typedef long long ll;
vector<ll> st;
ll ans;
ll inv6=(MOD+1ll)/6ll,inv2=(MOD+1ll)/2ll,v[5];
ll pres(ll x){
	return (1ll+x)*x/2ll%MOD;
}
ll calc(ll l1,ll r1,ll l2,ll r2,ll val){
	ll l=max(l1,l2);
	ll r=min(r1,r2);
	if(l<=r)
		return (pres(r)-pres(l-1)+(r-l+1)*(val%MOD)%MOD)%MOD;
	return 0;
}
ll count(ll L,ll l,ll r,ll d,ll u){
	if(L==0)
		return 	l<=0&&r>=0&&u>=0&&d<=0;
	ll val=(2ll*L+1ll)*(2ll*L+1ll);
	ll res=0;
	val-=L;
	if(u>=-L&&d<=-L)
		res=res+calc(l,r,-L,L,val);
	val-=2ll*L;
	if(l<=-L&&r>=-L)
		res=res+calc(d,u,-L+1,L-1,val);
	val-=2ll*L;
	if(u>=L&&d<=L)
		res=res+calc(l,r,-L,L,val);
	val-=2ll*L;
	if(l<=L&&r>=L)
		res=res+calc(d,u,-L+1,L-1,val);
	res=(res%MOD+MOD)%MOD;
	return res;
}	
int main(){
//	freopen("t.in.4","r",stdin);
//	freopen("wa.out","w",stdout);
	int n,q;
	SF("%d%d",&n,&q);
	for(int i=1;i<=q;i++){
		ans=0;
		ll x0,x1,y0,y1;
		SF("%lld%lld%lld%lld",&x0,&y0,&x1,&y1);
		st.clear();
		st.push_back(0);
		st.push_back(1);
		st.push_back(2);
		st.push_back(abs(x0)-1);
		st.push_back(abs(x0));
		st.push_back(abs(x0)+1);
		st.push_back(abs(x1)-1);
		st.push_back(abs(x1));
		st.push_back(abs(x1)+1);
		st.push_back(abs(y0)-1);
		st.push_back(abs(y0));
		st.push_back(abs(y0)+1);
		st.push_back(abs(y1)-1);
		st.push_back(abs(y1));
		st.push_back(abs(y1)+1);
		sort(st.begin(),st.end());
		for(int i=1;i<int(st.size());i++){
			int las=st[i-1];
			int now=st[i];
			if(las<0||las==now)
				continue;
			int k=now-las;
			if(k<4){
				for(int i=las;i<now;i++)
					(ans+=count(i,x0,x1,y0,y1))%=MOD;
			}
			else{
				for(int i=0;i<4;i++)
					v[i]=count(i+las,x0,x1,y0,y1);
				ll A=(((v[3]-2ll*v[2]+v[1])-(v[2]-2ll*v[1]+v[0]))%MOD+MOD)%MOD;
				A=A*inv6%MOD;
				v[1]=(v[1]-A+MOD)%MOD;
				v[2]=(v[2]-8ll*A%MOD+MOD)%MOD;
				v[3]=(v[3]-27ll*A%MOD+MOD)%MOD;
				ll B=((v[2]-2ll*v[1]+v[0])%MOD+MOD)%MOD;
				B=B*inv2%MOD;
				v[1]=(v[1]-B+MOD)%MOD;
				v[2]=(v[2]-4ll*B%MOD+MOD)%MOD;
				v[3]=(v[3]-9ll*B%MOD+MOD)%MOD;
				ll C=v[1]-v[0];
				ll D=v[0];
				ll K=k;
				ll SUM=K*(K-1ll)%MOD*inv2%MOD;
//				PF("{%lld %lld %lld %lld(%lld,%lld,%lld,%lld) [%d %d]}\n",A,B,C,D,v[0],v[1],v[2],v[3],las,now);
				ans=(ans+K*D%MOD)%MOD;
				ans=(ans+SUM*C%MOD)%MOD;
				ans=(ans+K*(K-1ll)%MOD*(2ll*k-1ll)%MOD*inv6%MOD*B%MOD)%MOD;
				ans=(ans+SUM*SUM%MOD*A%MOD)%MOD;
			}
		}
		ans=(ans+MOD)%MOD;
		PF("%lld\n",ans);
	}
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值