2021牛客多校2－F

F－Girlfriend

题目大意

给出四个点A,B,C,DA, B, C, D
另有两点P1,P2P_1, P_2分别满足∣AP1∣≥k1∣BP1∣|AP_1| \geq k_1|BP_1|，∣CP2∣≥k2∣DP2∣|CP_2| \geq k_2|DP_2|
求P1,P2P_1, P_2相交范围的体积
Time : 1000 ms
Memory: 262144 kB

解题思路及分析

阿波罗尼斯圆，在三维情况下成为一个球，根据这个可以直接写出球的球心坐标和半径
不知道公式的也可以现场去推公式（因为是轮换对称的，只推一个坐标就可以写出另外两个）
我写的时候因为忘了结论就是现场推的但是推错好多次…… 这里给出结论：

若，A(x1,y1,z1),B(x2,y2,z2)，且P点满足∣AP∣=k∣BP∣　　则P点的轨迹方程为一个球　　　球心坐标：O(k2x2−x1k2−1,k2y2−y1k2−1,k2z2−z1k2−1)　半径：r=k∣AB∣∣k2−1∣若，A(x_1, y_1, z_1), B(x_2,y_2,z_2) ，且P点满足|AP| = k|BP| 　\\　\\ 则P点的轨迹方程为一个球　\\　\\ 　 球心坐标：O(\frac{k^2x_2-x_1}{k^2-1},\frac{k^2y_2-y_1}{k^2-1},\frac{k^2z_2-z_1}{k^2-1} )\\　\\ 半径：r = \frac{k|AB|}{|k^2-1|}

找到两个球的球心和半径后，可根据球心距与半径的关系判断两球的关系，进而求出两球相交的体积。
两球相交部分的体积可以利用积分去求，涉及到球冠等概念，可以当作结论记住，具体证明这里不再赘述，给出结论：（相交但不内含时满足）

对于大球，根据余弦定理，可得cos⁡α=R2+d2−r22Rd记：h1=R(1−cos⁡α)积分可得球冠部分的体积：πRh12−πh133　　同理对于小球cos⁡β=r2+d2−R22rd记：h2=r(1−cos⁡β)积分可得球冠部分的体积：πrh22−πh233对于大球，根据余弦定理，可得 \\ \cos \alpha = \frac{R^2+d^2-r^2}{2Rd} \\ 记：h_1 = R(1-\cos \alpha ) \\ 积分可得球冠部分的体积：\pi Rh_1^2- \frac{\pi h_1^3}{3}\\ 　\\ \\　\\ 同理对于小球\\ \cos \beta = \frac{r^2+d^2-R^2}{2rd} \\ 记：h_2 = r(1-\cos \beta ) \\ 积分可得球冠部分的体积：\pi rh_2^2- \frac{\pi h_2^3}{3}

将两部分求和即为相交部分的体积。

AC代码

#include <bits/stdc++.h>
typedef long long llong;
const int N = 500000 + 5;
const double PI = 3.1415926535;
using namespace std;
 
int x[10], y[10], z[10];
int k11, k22;
double x1, y_1, z1, r1;	// cmath库里有y1，牛客上会CE于是改成了y_1
double x2, y2, z2, r2;
 
 
double dis(double x1, double y_1, double z1, double x2, double y2, double z2) {
    return (x1 - x2)*(x1 - x2) + (y_1 - y2)*(y_1 - y2) + (z1 - z2)*(z1 - z2);
}
 
double getVolume(double r) {
    return 4.0 * PI * r * r * r / 3.0;
}
 
int main() {
    #ifdef LOCAL
        freopen("E:/code/ACM/in.in", "r", stdin);
        freopen("E:/code/ACM/out.out", "w", stdout);
    #endif
     
    int T;  cin >> T;
    while (T--) {
        for (int i = 0; i < 4; i++) {
            scanf("%d%d%d", &x[i], &y[i], &z[i]);
        }
        scanf("%d%d", &k11, &k22);
        double k1 = k11;
        double k2 = k22;
        x1 = (k1 * k1 * x[1] - x[0]) / (k1 * k1 - 1);
        y_1 = (k1 * k1 * y[1] - y[0]) / (k1 * k1 - 1);
        z1 = (k1 * k1 * z[1] - z[0]) / (k1 * k1 - 1);
        r1 = k1 * sqrt((x[0] - x[1]) * (x[0] - x[1]) + (y[0] - y[1]) * (y[0] - y[1]) + (z[0] - z[1]) * (z[0] - z[1])) / (k1 * k1 - 1);
 
        x2 = (k2 * k2 * x[3] - x[2]) / (k2 * k2 - 1);
        y2 = (k2 * k2 * y[3] - y[2]) / (k2 * k2 - 1);
        z2 = (k2 * k2 * z[3] - z[2]) / (k2 * k2 - 1);
        r2 = k2 * sqrt((x[3] - x[2]) * (x[3] - x[2]) + (y[3] - y[2]) * (y[3] - y[2]) + (z[3] - z[2]) * (z[3] - z[2])) / (k2 * k2 - 1);
 
        double dd = dis(x1, y_1, z1, x2, y2, z2);
        double d = sqrt(dd);
 
        double r = min(r1, r2);
        double R = max(r1, r2);
 
        if (d >= r + R) {
            printf("0.000\n");
        } else if (d + r <= R) {
            printf("%.3f\n", getVolume(r));
        } else {
            double ca = (R * R + dd - r * r) / (2.0 * R * d);
            double h1 = R * (1 - ca);
            double cb = (r * r + dd - R * R) / (2.0 * r * d);
            double h2 = r * (1 - cb);
            double volume = ((3.0 * R - h1) * h1 * h1 * PI) / 3.0 + ((3.0 * r - h2) * h2 * h2 * PI) / 3.0;
            printf("%.3f\n", volume);
        }
         
    }
    return 0;
}