更向天涯

用户名：20110618li

访问量：42309

签　名：

向天涯进发！

等　　级
排　　名	12
经　　验	25711
参赛次数	28
文章发表	502
年　　龄	45
在职情况	在职
学　　校	柏苑中心小学
专　　业	编程

　　自我简介：

无情未必真豪杰，怜子如何不丈夫。

TA的其他文章

1319: 没有上司的晚会
浏览：63

1322: 沙子合并
浏览：132

1323: 能量项链
浏览：82

你可能喜欢

剪刀石头布（C语言代码）不知道怎么直接在scanf中用枚举变量

WU-C语言程序设计教程（第三版）课后习题11.11 （C++代码）（想学链表的可以看看）

1082: Tetrahedrons and Spheres

　　作者: 更向天涯　　　　发表时间：2023-01-31 11:35:12 浏览：136　|　评论：1　

原题链接：Tetrahedrons and Spheres

```cpp
#include<iostream>
#include<iomanip>
#include<cmath>
#include<algorithm>
#include<vector>
#define sqr(x)(x)*(x)
using namespace std;
const double eps=1e-6,inf=8;
inline bool cmpDouble(const double &a,const double &b)
{
	return fabs(a-b)<eps;
}
struct Tpoint
{
    double x,y;
	Tpoint(){}
    Tpoint(double a,double b)
    {
    	x=a,y=b;
	}
	inline double norm()
	{
		return sqrt(sqr(x)+sqr(y));
	}
};
inline Tpoint operator+(const Tpoint &a,const Tpoint &b)
{
	return Tpoint(a.x+b.x,a.y+b.y);
}
inline Tpoint operator-(const Tpoint &a,const Tpoint &b)
{
	return Tpoint(a.x-b.x,a.y-b.y);
}
inline Tpoint operator*(const double &a,const Tpoint &b)
{
	return Tpoint(a*b.x,a*b.y);
}
inline Tpoint operator*(const Tpoint &b,const double &a)
{
	return Tpoint(a*b.x,a*b.y);
}
inline Tpoint operator/(const Tpoint &a,const double &b)
{
	return Tpoint(a.x/b,a.y/b);
}
inline bool operator<(const Tpoint &a,const Tpoint &b)
{
	return a.x+eps<b.x||fabs(a.x-b.x)<eps&&a.y+eps<b.y;
}
inline bool operator==(const Tpoint &a,const Tpoint &b)
{
	return fabs(a.x-b.x)<eps&&fabs(a.y-b.y)<eps;
}
inline double det(const Tpoint &a,const Tpoint &b)
{
	return a.x*b.y-a.y*b.x;
}
inline double dot(const Tpoint &a,const Tpoint &b)
{
	return a.x*b.x+a.y*b.y;
}
struct P3
{
    double x,y,z;
    P3(){}
    P3(double a,double b,double c)
	{
		x=a;
		y=b;
		z=c;
    }
	inline void read()
	{
		cin>>x>>y>>z;
		return;
	}
};
inline P3 operator+(const P3 &a,const P3 &b)
{
	return P3(a.x+b.x,a.y+b.y,a.z+b.z);
}
inline P3 operator-(const P3 &a,const P3 &b)
{
	return P3(a.x-b.x,a.y-b.y,a.z-b.z);
}
inline P3 operator*(const double &a,const P3 &b)
{
	return P3(a*b.x,a*b.y,a*b.z);
}
inline P3 operator*(const P3 &b,const double &a)
{
	return P3(a*b.x,a*b.y,a*b.z);
}
inline P3 operator/(const P3 &a,const double &b)
{
	return P3(a.x/b,a.y/b,a.z/b);
}
struct Tcir
{
	double r;
	Tpoint o;
	Tcir(){}
	Tcir(Tpoint x,double y)
	{
		o=x,r=y;
	}
};
vector<Tcir>Circle;
typedef vector<Tpoint>TP;
vector<TP>Poly;
struct Tinter
{
	double x,y,Area,mid;
	int delta;
	Tinter(){}
	Tinter(double xx,double yy,double mm,int dd,double aa)
	{
		x=xx;
		y=yy;
		mid=mm;
		delta=dd;
		Area=aa;
	}
};
inline bool operator<(const Tinter &a,const Tinter &b)
{
	return a.mid>b.mid+eps;
}
inline bool operator==(const Tinter &a,const Tinter &b)
{
	return fabs(a.mid-b.mid)<eps;
}
vector<Tinter>inter;
vector<double>bak;
inline double dist(const Tpoint &a,const Tpoint &b)
{
	return sqr(a.x-b.x)+sqr(a.y-b.y);
}
inline void Add(double x)
{
	bak.push_back(x);
	return;
}
inline void CircleIntersectCircle(const Tcir &a,const Tcir &b)
{
	double l=dist(a.o,b.o),s=((a.r-b.r)*(a.r+b.r)/l+1)/2,t=sqrt(-(l-sqr(a.r+b.r))*(l-sqr(a.r-b.r))/(l*l*4)),ux=b.o.x-a.o.x,uy=b.o.y-a.o.y,ix=a.o.x+s*ux+t*uy,iy=a.o.y+s*uy-t*ux,jx=a.o.x+s*ux-t*uy,jy=a.o.y+s*uy+t*ux;
	Add(ix);
	Add(jx);
	return;
}
inline bool InLine(const Tpoint &a,const Tpoint &b,const Tpoint &c)
{
	return fabs(det(b-a,c-a))<eps&&dot(a-c,b-c)<eps;
}
inline void LineToLine(const Tpoint &a,const Tpoint &b,const Tpoint &c,const Tpoint &d)
{
	double s1=det(c-a,b-a),s2=det(b-a,d-a);
	if(s1*s2<-eps)
		return;
	Tpoint e=c+(d-c)*s1/(s1+s2);
	if(InLine(a,b,e)&&InLine(c,d,e))
		Add(e.x);
	return;
}
inline void LineToCircle(const Tpoint &a,const Tpoint &b,const Tcir &c)
{
	double h=fabs(det(c.o-a,b-a))/(b-a).norm();
	if(h>c.r+eps)
		return;
	double lamda=dot(c.o-a,b-a);
	lamda/=dist(a,b);
	Tpoint x=a+(b-a)*lamda;
	double d=sqrt(sqr(c.r)-sqr(h));
	d/=(b-a).norm();
	Tpoint e=x+(b-a)*d,f=x-(b-a)*d;
	if(InLine(a,b,e))
		Add(e.x);
	if(InLine(a,b,f))
		Add(f.x);
	return;
}
inline void CircleToCircle()
{
	for(int i=0;i<Circle.size();i++)
	{
		Add(Circle[i].o.x-Circle[i].r);
		Add(Circle[i].o.x+Circle[i].r);
		Add(Circle[i].o.x);
		for(int j=i+1;j<Circle.size();j++)
			if(dist(Circle[i].o,Circle[j].o)<=sqr(Circle[i].r+Circle[j].r)&&dist(Circle[i].o,Circle[j].o)>=sqr(Circle[i].r-Circle[j].r))
				CircleIntersectCircle(Circle[i],Circle[j]);
	}
	return;
}
inline void CircleToPoly()
{
	for(int i=0;i<Circle.size();i++)
		for(int j=0;j<Poly.size();j++)
			for(int v=0;v<Poly[j].size();v++)
				LineToCircle(Poly[j][v],Poly[j][(v+1)%Poly[j].size()],Circle[i]);
	return;
}
inline void PolyToPoly()
{
	for(int i=0;i<Poly.size();i++)
		for(int u=0;u<Poly[i].size();u++)
			Add(Poly[i][u].x);
	for(int i=0;i<Poly.size();i++)
		for(int j=i+1;j<Poly.size();j++)
			for(int u=0;u<Poly[i].size();u++)
				for(int v=0;v<Poly[j].size();v++)
					LineToLine(Poly[i][u],Poly[i][(u+1)%Poly[i].size()],Poly[j][v],Poly[j][(v+1)%Poly[j].size()]);
	return;
}
inline void Get(const Tcir &c,double x,double &l,double &r)
{
    double y=fabs(c.o.x-x),d=sqrt(fabs(sqr(c.r)-sqr(y)));
	l=c.o.y+d;
	r=c.o.y-d;
	return;
}
inline double arcArea(const Tcir &a,double l,double x,double r,double y)
{
	double len=sqrt(sqr(l-r)+sqr(x-y)),d=sqrt(sqr(a.r)-sqr(len)/4.0),angle=atan(len/2.0/d);
	return fabs(angle*sqr(a.r)-d*len/2.0);
}
inline void Get_Interval(const Tcir &a,double l,double r)
{
	double L1,L2,R1,R2,M1,M2;
	Get(a,l,L1,L2);
	Get(a,r,R1,R2);
	Get(a,(l+r)/2.0,M1,M2);
	int D1=1,D2=-1;
	double A1=arcArea(a,l,L1,r,R1),A2=arcArea(a,l,L2,r,R2);
	inter.push_back(Tinter(L1,R1,M1,D1,A1));
	inter.push_back(Tinter(L2,R2,M2,D2,A2));
	return;
}
inline double calcSlice(double xl,double xr)
{
	inter.clear();
	double lmost=-inf,rmost=inf;
	for(int i=0;i<Poly.size();i++)
	{
		int cc=0;
		Tinter I[5];
		for(int u=0;u<Poly[i].size();u++)
		{
			Tpoint x=Poly[i][u],y=Poly[i][(u+1)%Poly[i].size()];
			double l=min(x.x,y.x),r=max(x.x,y.x);
			if(l<=xl+eps&&xr<=r+eps)
			{
				if(fabs(l-r)<eps)
					continue;
				Tpoint d=y-x,Left=x+d/d.x*(xl-x.x),Right=x+d/d.x*(xr-x.x),Mid=(Left+Right)/2;
				I[cc++]=Tinter(Left.y,Right.y,Mid.y,1,0);
			}
		}
		sort(I,I+cc);
		if(cc==2)
		{
			I[1].delta=-1;
			inter.push_back(I[0]);
			inter.push_back(I[1]);
			lmost=max(lmost,I[1].mid);
			rmost=min(rmost,I[0].mid);
		}
	}
	for(int i=0;i<Circle.size();i++)
		if(fabs(Circle[i].o.x-xl)<Circle[i].r+eps&&fabs(Circle[i].o.x-xr)<Circle[i].r+eps)
			Get_Interval(Circle[i],xl,xr);
	if(!inter.size())
		return 0;
	double ans=0;
	sort(inter.begin(),inter.end());
	int cnt=0;
	for(int i=0;i<inter.size();i++)
	{
		if(cnt>0)
		{
			ans+=(fabs(inter[i-1].x-inter[i].x)+fabs(inter[i-1].y-inter[i].y))*(xr-xl)/2.0;
			ans+=inter[i-1].delta*inter[i-1].Area;
			ans-=inter[i].delta*inter[i].Area;
		}
		cnt+=inter[i].delta;
	}
	return ans;
}
int n,m;
struct Poly4
{
    P3 p[4];
}a[105];
struct Sphere
{
    P3 o;
    double r;
    inline void read()
	{
        o.read();
		cin>>r;
		return;
	}
}b[105];
inline bool equal(double a,double b)
{
	return fabs(a-b)<eps;
}
inline bool InterSect(const Tpoint &a,const Tpoint &b,const Tpoint &c,const Tpoint&d)
{
    double s1=det(b-a,c-a),s2=det(b-a,d-a);
    if(s1*s2>-eps)
		return false;
    s1=det(d-c,a-c),s2=det(d-c,b-c);
    if(s1*s2>-eps)
		return false;
    return true;
}
inline void ToHull(vector<Tpoint>&a)
{
    sort(a.begin(),a.end());
    int hull[10],len,limit=1;
    hull[len=1]=0;
    for(int i=1;i<4;i++)
	{
        while(len>limit&&det(a[hull[len]]-a[hull[len-1]],a[i]-a[hull[len]])>=0)
			--len;
        hull[++len]=i;
    }
    limit=len;
    for(int i=2;i>=0;i--)
	{
        while(len>limit&&det(a[hull[len]]-a[hull[len-1]],a[i]-a[hull[len]])>=0)
			--len;
        hull[++len]=i;
    }
    vector<Tpoint>b=a;
    a.resize(len-1);
    for(int i=0;i<len-1;i++)
        a[i]=b[hull[i+1]];
    return;
}
inline double calcArea(double z)
{
    Poly.clear();
    Circle.clear();
    bak.clear();
    for(int i=0;i<n;i++)
	{
        vector<Tpoint>cross;
        for(int j=0;j<4;j++)
            for(int k=j+1;k<4;k++)
                if(!equal(a[i].p[j].z,a[i].p[k].z))
				{
                    double l=min(a[i].p[j].z,a[i].p[k].z),r=max(a[i].p[j].z,a[i].p[k].z);
                    if(l<=z+eps&&z<=r+eps)
					{
                        P3 d=a[i].p[k]-a[i].p[j];
                        d=d/d.z;
                        d=d*(z-a[i].p[j].z);
                        d=d+a[i].p[j];
                        Tpoint x(d.x,d.y);
                        cross.push_back(x);
                    }
                }
        sort(cross.begin(),cross.end());
        cross.erase(unique(cross.begin(),cross.end()),cross.end());
        if(cross.size()>2)
		{
			if(cross.size()>4)
			{
                puts("Too Many Points!!!");
                while(true);
            }
            if(cross.size()==4)
                ToHull(cross);
            Poly.push_back(cross);
        }
    }
    for(int i=0;i<m;i++)
        if(fabs(z-b[i].o.z)+eps<b[i].r)
        {
            Tpoint o(b[i].o.x,b[i].o.y);
            double r=sqrt(sqr(b[i].r)-sqr(z-b[i].o.z));
            Circle.push_back(Tcir(o,r));
        }
    CircleToCircle();
    CircleToPoly();
    PolyToPoly();
    sort(bak.begin(),bak.end());
    double res=0;
    for(int i=0;i+1<bak.size();i++)
        if(fabs(bak[i+1]-bak[i])>eps)
            res+=calcSlice(bak[i],bak[i+1]);
	return res;
}
int main()
{
	while(cin>>n>>m&&n+m)
	{
		for(int i=0;i<n;i++)
			for(int j=0;j<4;j++)
				a[i].p[j].read();
		for(int i=0;i<m;i++)
			b[i].read();
		const int Block=4000;
		double last=0,Ans=calcArea(-inf)+calcArea(inf),h=(inf+inf)/(double)Block;
		for(int i=1;i<Block;i+=2)
			Ans+=4*calcArea(-inf+i*h);
		for(int i=2;i<Block;i+=2)
			Ans+=2*calcArea(-inf+i*h);
		Ans*=h/3.0;
		cout<<fixed<<setprecision(3)<<Ans<<endl;
	}
	return 0;
}
```

0.0分

3 人评分

看不懂代码解释一下代码? 或者生成一段代码？试试AI编程助手吧

　　评论区

更向未来(ikun)

妙啊

2023-02-24 18:20:53

更向天涯

TA的其他文章

你可能喜欢

C语言程序设计教程（第三版）课后习题5.7 （C语言代码）

C语言程序设计教程（第三版）课后习题1.5 （C语言代码）

分糖果（C++代码）

剪刀石头布（C语言代码）不知道怎么直接在scanf中用枚举变量

WU-C语言程序设计教程（第三版）课后习题11.11 （C++代码）（想学链表的可以看看）

1082: Tetrahedrons and Spheres

评论区

精彩推荐

字符串输入输出函数（C语言代码）

蓝桥杯历届试题-翻硬币（C++代码）

A+B for Input-Output Practice (III) （C语言代码）

C语言程序设计教程（第三版）课后习题11.8 （C语言代码）

C二级辅导-温度转换（C语言代码）

众数问题（C语言代码）

字符串的修改（C语言代码）

1025题初学者，求帮忙看下,不知道哪错了

C语言程序设计教程（第三版）课后习题10.7 指针（C语言代码）

1227题解

更向天涯

TA的其他文章

你可能喜欢

C语言程序设计教程（第三版）课后习题5.7 （C语言代码）

C语言程序设计教程（第三版）课后习题1.5 （C语言代码）

分糖果 （C++代码）

剪刀石头布 （C语言代码）不知道怎么直接在scanf中用枚举变量

WU-C语言程序设计教程（第三版）课后习题11.11 （C++代码）（想学链表的可以看看）

1082: Tetrahedrons and Spheres

评论区

精彩推荐

字符串输入输出函数 （C语言代码）

蓝桥杯历届试题-翻硬币 （C++代码）

A+B for Input-Output Practice (III) （C语言代码）

C语言程序设计教程（第三版）课后习题11.8 （C语言代码）

C二级辅导-温度转换 （C语言代码）

众数问题 （C语言代码）

字符串的修改 （C语言代码）

1025题 初学者，求帮忙看下,不知道哪错了

C语言程序设计教程（第三版）课后习题10.7 指针（C语言代码）

1227题解

分糖果（C++代码）

剪刀石头布（C语言代码）不知道怎么直接在scanf中用枚举变量

　　评论区

字符串输入输出函数（C语言代码）

蓝桥杯历届试题-翻硬币（C++代码）

C二级辅导-温度转换（C语言代码）

众数问题（C语言代码）

字符串的修改（C语言代码）

1025题初学者，求帮忙看下,不知道哪错了