CodeChef2015July挑战赛 NTHCIR

https://www.codechef.com/problems/NTHCIR

会发现所有的圆都跟c1和c2相切

那么以c1,c2的切点(0,c1.r)为反演点，反演半径为2*c1.r(这个随便选)，反演出两条直线y=ky和y=ry，c3-cn都跟ci-1和c1,c2相切，那么就是半径为(uy-ly)/2的圆从左到右排列过去了。

由于c1.r>=c2.r+c3.r，那么可以用余弦公式求出角度的cos和sin值(杨程好习惯，尽量不用反三角函数)，然后把c3.o的坐标算出来，反演过去，就知道c3的反演圆的坐标在哪了，那么cn的反演圆圆心坐标就横坐标加一下就好了

 

#include<bits/stdc++.h>
using namespace std;

const double eps=1e-10;

inline int sgn(int x)
{
	if(x>-eps || x<eps) return 0;
	if(x>0) return 1;
	else	return -1;
}

struct point
{
	double x,y;
	point(double a=0,double b=0)
	{
		x=a;y=b;
	}
	point operator + (const point &b)const
	{
		return point(x+b.x,y+b.y);
	}
	point operator - (const point &b)const
	{
		return point(x-b.x,y-b.y);
	}
	point operator * (const double &t)const
	{
		return point(x*t,y*t);
	}
};
struct circle
{
	point o;
	double r;
};
inline double sqr(double x){return x*x;}
inline double det(const point &a,const point &b)
{
	return a.x*b.y-b.x*a.y;
}
inline double dist(const point &a,const point &b)
{
	return sqrt(sqr(a.x-b.x)+sqr(a.y-b.y));
}
inline double mysqrt(double x){return sqrt(max(0.0,x));} 
inline circle inversion(circle a,point p,double R)
{
    circle res;
    double d1=dist(p,a.o);
    res.r=a.r*R*R/(sqr(d1)-sqr(a.r));
    double d2=R*R/(d1-a.r)-res.r;
    res.o=p+(a.o-p)*(d2/d1);
    return res;
}

int t,n,P,m,b;
double ans;
circle a,c[5];

inline void prework()
{
	scanf("%d",&t);
	scanf("%d%d%d%d",&n,&P,&m,&b);
	for(int i=1;i<=4;i++)
		scanf("%lf",&c[i].r);
}

inline void mainwork()
{
	c[1].o=point(0,0);c[2].o=point(0,c[1].r-c[2].r);
	double len1=c[1].r-c[3].r,len2=c[1].r-c[2].r;
	double len3=c[2].r+c[3].r;
	double cosa=(sqr(len1)+sqr(len2)-sqr(len3))/(2*len1*len2);
	cosa=min(max(-1.0,cosa),1.0);
	double sina=mysqrt(1.0-sqr(cosa));
	sina=min(max(-1.0,sina),1.0);
	c[3].o=point(-len1*sina,len1*cosa);
	a.o=point(0,c[1].r);a.r=2*c[1].r;
	double uy=-c[1].r,ly=c[1].r-sqr(a.r)/(2*c[2].r);
	circle c3=inversion(c[3],a.o,a.r);
	circle now;
	ans=0;
	for(int i=1;i<=t;i++)
	{
		n=(1ll*P*n)%m+b;
		if(n<=3)
			ans+=c[n].r;
		else
		{
			now.r=(uy-ly)/2;
			now.o=point(c3.o.x+(n-3)*(uy-ly),c3.o.y);
			now=inversion(now,a.o,a.r);
			ans+=now.r;
		}
	}
}

inline void print()
{
	printf("%.6f",ans);
}

int main()
{
	prework();
	mainwork();
	print();
	return 0;
}