• 【XSY2903】B 莫比乌斯反演


    题目描述

      有一个(n imes n)的网格,除了左下角的格子外每个格子的中心里都有一个圆,每个圆的半径为(R),问你在左下角的格子的中心能看到多少个圆。

      (nleq {10}^9,R_0leq {10}^6,R=frac{R_0}{{10}^6})

    题解

      先枚举格子,判断是否能看到。

      显然能看到就一定能看到圆的中心。

      条件是:

    [egin{align} frac{lvert Ax+By vert}{sqrt{A^2+B^2}}&leq R\ lvert Ax+By vert&leq Rsqrt{A^2+B^2}\ {(Ax+By)}^2&leq R^2(A^2+B^2)\ frac{{(Ax+By)}^2}{A^2+B^2}&leq R^2\ end{align} ]

      如果(gcd(a,b) eq 1),那么一定会被挡住。

      否则一定存在(x,y,s.t.ax+by=1)

    [egin{align} frac{1}{A^2+B^2}&leq R^2\ A^2+B^2&geq frac{1}{R^2} end{align} ]

      那么现在我们要求的就是

    [egin{align} &sum_{i=1}^nsum_{j=1}^n[gcd(i,j)=1][i^2+j^2<R^2]\ =&sum_{d=1}^nmu(d)sum_{1leq i,jleq frac{n}{d}}[i^2+j^2<frac{R^2}{D^2}] end{align} ]

      时间复杂度:(O(frac 1Rlog frac 1R))

    题解

    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #include<cstdlib>
    #include<ctime>
    #include<utility>
    #include<cmath>
    #include<functional>
    using namespace std;
    typedef long long ll;
    typedef unsigned long long ull;
    typedef pair<int,int> pii;
    typedef pair<ll,ll> pll;
    void sort(int &a,int &b)
    {
    	if(a>b)
    		swap(a,b);
    }
    void open(const char *s)
    {
    #ifndef ONLINE_JUDGE
    	char str[100];
    	sprintf(str,"%s.in",s);
    	freopen(str,"r",stdin);
    	sprintf(str,"%s.out",s);
    	freopen(str,"w",stdout);
    #endif
    }
    int rd()
    {
    	int s=0,c;
    	while((c=getchar())<'0'||c>'9');
    	do
    	{
    		s=s*10+c-'0';
    	}
    	while((c=getchar())>='0'&&c<='9');
    	return s;
    }
    void put(int x)
    {
    	if(!x)
    	{
    		putchar('0');
    		return;
    	}
    	static int c[20];
    	int t=0;
    	while(x)
    	{
    		c[++t]=x%10;
    		x/=10;
    	}
    	while(t)
    		putchar(c[t--]+'0');
    }
    int upmin(int &a,int b)
    {
    	if(b<a)
    	{
    		a=b;
    		return 1;
    	}
    	return 0;
    }
    int upmax(int &a,int b)
    {
    	if(b>a)
    	{
    		a=b;
    		return 1;
    	}
    	return 0;
    }
    ll ans;
    ll calc(ll n,ll r)
    {
    	ll j=0;
    	while(j<n&&(j+1)*(j+1)<=r)
    		j++;
    	ll s=0;
    	for(int i=1;i<=n&&(ll)i*i<=r;i++)
    	{
    		while((ll)i*i+j*j>r)
    			j--;
    		s+=j;
    	}
    	return s;
    }
    int miu[1000010];
    int b[1000010];
    int pri[1000010];
    int cnt;
    int main()
    {
    	open("b");
    	ll n,r;
    	scanf("%lld%lld",&n,&r);
    	r=(ll)1e12/r/r;
    	miu[1]=1;
    	n--;
    	for(int i=2;i<=1000000;i++)
    	{
    		if(!b[i])
    		{
    			pri[++cnt]=i;
    			miu[i]=-1;
    		}
    		for(int j=1;j<=cnt&&i*pri[j]<=1000000;j++)
    		{
    			b[i*pri[j]]=1;
    			if(i%pri[j]==0)
    				break;
    			miu[i*pri[j]]=-miu[i];
    		}
    	}
    	for(int i=1;(ll)i*i<=r;i++)
    		if(miu[i])
    			ans+=miu[i]*calc(n/i,r/i/i);
    	ans+=2;
    	printf("%lld
    ",ans);
    	return 0;
    }
    
  • 相关阅读:
    C++获取时间函数
    平滑算法:三次样条插值(Cubic Spline Interpolation)
    为什么想要交谈?
    c++日常小问题
    看板娘
    世界碰撞算法原理和总结(sat gjk)
    转载c++默认初始化文章--google翻译
    从4行代码看右值引用(转载 《程序员》2015年1月刊)
    c++模板特例化 函数模板(非法使用显式模板参数 )
    InverseTransformPoint 函数问题
  • 原文地址:https://www.cnblogs.com/ywwyww/p/9045627.html
Copyright © 2020-2023  润新知