Top

球与球之间的相交体积与面积


1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
const double PI = acos(-1.0);
struct point {
double x,y,z;
point() {}
void Init(double a, double b,double c){
x = a;
y = b;
z = c;
}
point(double a, double b,double c){
x = a;
y = b;
z = c;
}
point operator -(const point &b)const { //返回减去后的新点
return point(x - b.x, y - b.y,z-b.z);
}
point operator +(const point &b)const { //返回加上后的新点
return point(x + b.x, y + b.y,z+b.z);
}
//数乘计算
point operator *(const double &k)const { //返回相乘后的新点
return point(x * k, y * k,z*k);
}
point operator /(const double &k)const { //返回相除后的新点
return point(x / k, y / k,z/k);
}
double operator *(const point &b)const { //点乘
return x*b.x + y*b.y+z*b.z;
}
};
double dist(point p1, point p2) { //返回平面上两点距离
return sqrt((p1 - p2)*(p1 - p2));
}
struct sphere {//球
double r;
point centre;
};
void SphereInterVS(sphere a, sphere b,double &v) {
double d = dist(a.centre, b.centre);//球心距
if(d>=a.r+b.r)return;
if(a.r-d>=b.r){
v+=PI*4.0/3.0*b.r*b.r*b.r;
return;
}
if(a.r<=b.r-d){
v+=PI*4.0/3.0*a.r*a.r*a.r;
return;
}
double t = (d*d + a.r*a.r - b.r*b.r) / (2.0 * d);//
double h = sqrt((a.r*a.r) - (t*t)) * 2;//h1=h2,球冠的高
double angle_a = 2 * acos((a.r*a.r + d*d - b.r*b.r) / (2.0 * a.r*d)); //余弦公式计算r1对应圆心角,弧度
double angle_b = 2 * acos((b.r*b.r + d*d - a.r*a.r) / (2.0 * b.r*d)); //余弦公式计算r2对应圆心角,弧度
double l1 = ((a.r*a.r - b.r*b.r) / d + d) / 2;
double l2 = d - l1;
double x1 = a.r - l1, x2 = b.r - l2;//分别为两个球缺的高度
double v1 = PI*x1*x1*(a.r - x1 / 3);//相交部分r1圆所对应的球缺部分体积
double v2 = PI*x2*x2*(b.r - x2 / 3);//相交部分r2圆所对应的球缺部分体积
v += v1 + v2;//相交部分体积
}
文章版权为Anoyer博客所有,转载请以链接形式标明本文地址