2014-2015 ACM-ICPC Northeastern European Regional Contest (NEERC 14) 题解

D. Damage Assessment
I. Improvements

没注意到文件IO让罚时上天，我是 * *
gym100553

回到顶部

D. Damage Assessment

牛客也有这道题

题意：两个球缺和一个圆柱组成了一个倾斜摆放的几何体，求其所盛水的体积。

思路

将原问题分为圆柱体和球缺两个独立部分，分别对水的截面面积做定积分得到体积。

求定积分的数值解用到 自适应Simpson算法

先考虑一个圆上面积问题：

（1）子问题

如图，已知 $R$ 和 $D(0\leq D \leq R)$ ，求左边弓形面积。

显然面积是 $F_1(R,D)=r^2\arccos(1-\frac{D}{R})-\sqrt{2RD-D^2}(R-D)$ ， $R<D\leq 2R$ 的情况同理。

（2）圆柱内水的体积

设倾斜角 $\omega =\arcsin(\frac{t}{l})$ ，先定义圆柱侧面长方形的左上端点的高度 $h_1=\frac{d}{\cos \omega}$ ，

则答案为 $\int_0^{\min(h,t+h_1)}F_2(x)dx$ ，其中 $F_2(x)$ 是图中蓝色水平面截圆柱体的面积， $x$ 是水平面的高度。

具体来说，考虑将一个 $F_2(x)$ 截出来的椭圆的长轴乘以一个 $\sin \omega$ ，得到以 $d$ 为直径的圆（其面积也乘了个 $\sin \omega$ ），对于特定的 $x$ 计算出损失的部分的“长度”（例如图中红线），即根据（1）算出左右两边损失了弓形后剩下的面积，最后再除以 $\sin \omega$ 仿射回去

（3）球缺内水的体积

以左边的球缺为例，不妨设水界面 $h<h_1$ ，换一个与（2）不同的坐标系：

紫色是水界面 $h$ ，红色箭头的指向就是积分方向（ $x$ 正半轴方向），红箭头的根部就是原点，蓝色是其中一个截面，对于固定的一个截面（知道这条蓝线的 $x$ 坐标），可以根据坐标系的知识（比如算一下直线和1/4圆交点）算出这条蓝线的高度，进而转化为 （1） 的情形。

（4）一些细节

本题的数据是带量纲的，注意读入和输出的单位换算。（可以直接把所有输入数除以100）

注意这两个边界：

当 $t=0$ ， $\sin \omega=0$ ，计算（2） $F_2(x)$ 的时候直接返回一个长方形面积就行
当 $t=l$ ， $\cos \omega=0$ ，此时左边的球缺是完整的，右边的球缺（ $h>t$ 时）的截面是个圆，需特判

时间复杂度 $O(能过)$ （与积分算法的实现有关），更多实现细节可以看代码：

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

const double PI = 3.1415926535897932384626433;

struct Integration{ // 自适应Simpson算法求定积分
    typedef double Func(double);
    const Func*f;Integration(const Func&g):f(&g){}
    typedef pair<double,double> pdd;
    pdd add(pdd a,pdd b)const{
        double mid=(a.first+b.first)/2;
        return (pdd){mid,f(mid)};
    }
    #define simpson(p1,p2,p3) (((p3).first-(p1).first)*((p1).second+4*(p2).second+(p3).second)/6)
    double find(pdd p1,pdd p3,pdd p5,double EPS,double res,int dep)const{
        pdd p2=add(p1,p3),p4=add(p3,p5);
        double fl=simpson(p1,p2,p3),fr=simpson(p3,p4,p5),d=(fl+fr-res)/15;
        if(abs(d)<=EPS&&dep<0)return fl+fr+d; 
        return find(p1,p2,p3,EPS/2,fl,dep-1)+find(p3,p4,p5,EPS/2,fr,dep-1);  
    }
    double operator()(double l,double r,double EPS=1e-7/*精度*/,int dep=16/*最小递归深度*/)const{
        pdd p1(l,f(l)),p3(r,f(r)),p2=add(p1,p3);
        return find(p1,p2,p3,EPS,simpson(p1,p2,p3),dep);
    }
    #undef simpson
};


double d, l, r, t, h; 


//(1) 弓形面积
double F1(double r, double d) {
    auto F0 = [](double r, double d) -> double { // D<=R 的情形
        return r * r * acos(1 - d / r) - sqrt(d * (2 * r - d)) * (r - d);
    };
    if (d == r)
        return PI * r * r / 2;
    else if (d < r)
        return F0(r, d);
    else
        return PI * r * r - F0(r, 2 * r - d);
}


//(2) 柱体截面积
double h1, cs, b, sc; 
double F2(double x) {
    if (t < 1e-6) return 2 * sqrt(x * (d - x)) * l; // 特判 直接返回一个长方形
    double res = sc;
    if (x < h1) res -= F1(b, (h1 - x) / cs);
    if (x > t) res -= F1(b, (x - t) / cs);
    return res / t * l;
}


//(3) 的左球缺
double k, y3, c;
double F3(double x) {
    double ym = sqrt(x * (2 * r - x));
    if (t + 1e-6 > l) return PI * ym * ym; // 特判 直接返回一个圆
    return F1(ym, ym + max(min(ym, y3 - k * x), -ym));
}


//(3) 的右球缺
double y4;
double F4(double x) {
    double ym = sqrt(x * (2 * r - x));
    return F1(ym, ym + max(min(ym, y4 + k * x), -ym));
}


signed main() {
    freopen("damage.in", "r", stdin);freopen("damage.out", "w", stdout); // 牛客上提交似乎要去掉这句？

    scanf("%lf%lf%lf%lf%lf", &d, &l, &r, &t, &h);
    d /= 100; l /= 100; r /= 100; t /= 100; h /= 100;// 所有输入数除以100以和输出的单位同步

    // 求解全部辅助参数
    cs = sqrt(1 - (t / l) * (t / l));  // cos w
    h1 = d * cs;                       // 圆柱侧面长方形的左上端点的高度
    b = d / 2;                         // 圆柱侧面半径
    sc = b * b * PI;                   // 圆柱侧面(直径为d的圆)积       
    c = r - sqrt(r * r - b * b);       // 球缺的厚度           
    k = t / sqrt(l * l - t * t);               
    y3 = b - (h1 - h) / cs + k * c;             
    y4 = (h - t) / cs - b - k * c;           


    // 根据几何情况统计答案
    double ans = (Integration(F2))(0, min(h, t + h1)) + (Integration(F3))(0, c);
    if (t + 1e-6 > l) {
        if (h > t) ans += (Integration(F3))(0, h - t);
    } else {
        ans += (Integration(F4))(0, c);
    }
    
    printf("%.9lf\n", ans);
    return 0;
}

回到顶部

I. Improvements

想一下合法的调整结果应该长啥样。

想象一下一个人在X0..n上一步一步地跳跃，将轨迹连成一些弧线，可以发现这些弧线不能相交，形象的说，只能是:新的弧线不超过一个旧的弧线（被一个旧的弧线包含了），或者新的弧线是旧的弧线的基础上再往“外部”跳
倒过来想这个跳跃过程，从Xn一步步跳回X0，设每一步跳跃经过的所有点都被覆盖，维护一个“当前被覆盖过的点/线段集合”，我们发现这个集合起初是一个点，跳了一次后就是一个线段，此后始终是一个每一步都在不缩小或者扩展的线段（每一步都不能跳进这个线段的严格内部，要么往外部跳，要么往端点跳），（联想到最长上升/下降子序列）我们称其为“不缩序列”

答案就是序列Xn，Xn-1，Xn-2，....，X0形成的以X0结尾的最长不缩子序列
证明：一方面，原地不动的那些保留点本身必须构成一个不缩子序列；另一方面，选择一个最长不缩子序列作为保留后，我们一定有办法将其余的点调整到不影响这个不缩子序列的位置上（比如，将调整点向最近的保留点靠拢，这就保证调整完整个序列是一个不缩序列）。

注意到序列X其实是0..n的一个排列（Xi两两不同），那么其最长不缩子序列就等价于一个最长上升子序列和一个最长下降子序列的直接叠加。

所以：考虑X0，X1，...，Xn这个序列，选择一个Xi作为结尾，使得其最长上升子序列长度+最长下降子序列长度最大，这个最大值就是答案。

时间复杂度 $O(n\log n)$

#include <bits/stdc++.h>
using namespace std;
const int N = 1000009;

int n, a[N];

int f[N], g[N];  // f[i]是以Xi结尾的最长上升子序列长度 g[i]是最长下降子序列长度


int s[N];  // 树状数组求 LIS
void add(int x, int d) {
    for (int i = x; i <= n; i += (i & (-i))) s[i] = max(s[i], d);
}
int ask(int x) {
    int res = 0;
    for (int i = x; i > 0; i -= (i & (-i))) res = max(res, s[i]);
    return res;
}

signed main() {
    freopen("improvements.in", "r", stdin); 
    freopen("improvements.out", "w", stdout);

    // 将 0..n 的排列变换为 1..(n+1)的排列
    scanf("%d", &n);
    for (int i = 1; i <= n; ++i) scanf("%d", &a[i + 1]);
    ++n;
    for (int i = 1; i <= n; ++i) ++a[i];

    // 最长上升子序列
    for (int i = 1; i <= n; ++i) {
        f[i] = ask(a[i] - 1) + 1;
        add(a[i], f[i]);
    }

    // 最长下降子序列
    for (int i = 1; i <= n; ++i) s[i] = 0;
    for (int i = 1; i <= n; ++i) {
        int u = n - a[i] + 1;
        g[i] = ask(u - 1) + 1;
        add(u, g[i]);
    }

    int ans = f[1] + g[1] - 2;
    for (int i = 2; i <= n; ++i) ans = max(ans, f[i] + g[i] - 2);
    printf("%d\n", ans);
    return 0;
}