水果忍者

题目描述

    你知道水果忍者 吗 ？ 这 是一个智能手机上面很流行的游 戏 。

    这 个出色的游 戏 迎合了人类最喜 爱 的运 动轨 迹：抛物 线 。就是各种各 样 物体被扔出去之后都会形成的曲 线 运 动 。 现 在 这 里的 问题 也是关于抛物 线 的。

    现 在有 N 个水果，假 设 它 们 就是一个个在平面上的的 圆 ，有着同 样 的半径。每个水果有着各自的位置和移 动 速度。由于重力的影响（ 这 里假 设 重力常数 g=10) ，所有水果都从一开始就沿着抛物 线 运 动 。你的目 标 就是在 给 定的 时间 限制之内，一次切到最多的水果，即是找到在某一 时 刻的一条直 线 穿 过 最多的 圆 。

输入格式

输 入包含多 给测试 数据。

输 入的第一行 给 出了 测试 数据个数 T 。

每个 测试 数据首先 给 出三个整数，水果数目 N(<=100) ，水果的半径R(<=1000)， 还 有 时间 限制 E(<=1000) 。

接下来有N行，每行 给 出 4 个 实 数，分别 对应圆 心坐 标 (xi,yi) 和速度(vxi,vyi)。

所有数都是 采用标 准 单 位，所以你可以放心使用以下著名的物理公式来 计 算某一 时 刻水果所 处 的位置：

x' = x + vx * t

y' = y + vy * t - 0.5 * 10 * t * t

数据保 证 所有 的浮点数的 绝对值 小于 100000 。

输 出格式

对应 每个 测试 数据， 输 出一个整数，即在 时 限内可以同 时 切到的水果数目的最大 值 。

样例输入

2

4 1 2

2.0 4.0 1.0 0.0

5.0 4.0 2.0 0.0

2.0 12.0 1.5 -4.0

-2.0 -5.0 10.0 10.0

4 1 1

2.0 4.0 1.0 0.0

5.0 4.0 2.0 0.0

2.0 12.0 1.5 -4.0

-2.0 -5.0 10.0 10.0

样例输出

4

3

提示

请 注意一条直 线 “ 穿 过 ” 一个 圆 意味 着， 圆 心与直 线 的距离小于 圆 的半径加上精度 误 差 (1e-5) 。

#include <cstdio>

#include <cstring>

#include <cmath>

#include <cassert>

#include <algorithm>

#include <vector>

using namespace std;

const int maxn = 101;

const double INF = 1e8;

const double EPS = 1e-5;

const double MAXLOOP = 100;

struct seg {

    seg() { l = 1; r = 0; }

    seg(double il, double ir) { l = il; r = ir; }

    bool empty() { return r < l - EPS; }

    double l, r;

};

struct point {

    int mrk;

    double pos;

}points[maxn * 2];

struct vect {

   

    double d[20];

    int n;

    int size() {

        return n;

    }

    vect() {

        n= 0;

    }

    void push_back(double x) {

        d[n++] = x;

    }

    double operator [] (int index) {

        return d[index];

    }

};

int pn;

double x[maxn], y[maxn], vx[maxn],vy[maxn];

int n, r, e; // number, radius, time limit

int best; // maximum number of cutting

void add_seg(seg &s) {

    if (abs(s.l - s.r)<EPS) s.r = s.l;

    points[pn].mrk = 1;

    points[pn].pos = s.l;

    ++pn;

    points[pn].mrk = 0;

    points[pn].pos = s.r;

    ++pn;

}

vectdifferentiate(vect co) {

    vect res;

    for (int i = 1; i < co.size(); ++i) {

        res.push_back(co * i);

    }

    return res;

}

double substitute(vect co, double x) {

    double res = 0, tmp = 1;

    for (int i = 0; i < co.size(); ++i) {

        res += tmp * co;

        tmp *= x;

    }

    return res;

}

vectremove_dup(vect vec) {

    vect res;

    if (vec.size() == 0) return res;

    //sort(vec.begin(), vec.end());

    sort(vec.d, vec.d+vec.n);

    res.push_back(vec[0]);

    for (int i=1;i<vec.size();++i) {

        if (vec - vec[i-1] > EPS * 10)res.push_back(vec);

    }

    return res;

}

bool newton(vect co, vect dco, double x, double &res_x) {

    double last = -INF;

    int count = 0;

    while (count++ < MAXLOOP) {

        double f = substitute(co, x);

        double fp = substitute(dco, x);

        if ( abs(fp) < EPS ) {

            // no idea whether this works

            x += EPS;

            continue;

        }

        x -= f / fp;

        if ( abs(x - last) < EPS ) {

            res_x = x;

            return true;

        }

        last = x;

    }

    return false;

}

vectsolveeqs(vect co) {

    vect root;

   

    if (co.size() <= 2) {

        if (co.size() == 2 && co[1] != 0) {

            root.push_back(-co[0] / co[1]);

        }

        return root;

    }

    vect dco = differentiate(co); // differentail coefficient

    vect droot = solveeqs(dco);

    droot.push_back(-INF);

    droot.push_back(INF);

    sort(droot.d, droot.d+droot.n);

    //sort(droot.begin(), droot.end());

   

    for (int i=1;i<droot.size(); ++i) {

        double initx = (droot[i-1] + droot)* 0.5, res;

        if (newton(co, dco, initx, res)) {

            root.push_back(res);

        }

    }

    return remove_dup(root);

}

//sigma(ax^n) <= 0

vector<seg> solvele(vect co) {

    vect root = solveeqs(co);

    vector<seg> res;

    bool neg = substitute(co, -INF) < 0 ;

    double last = -INF;

    for (int i=0;i<root.size();++i) {

        if (neg) {

            res.push_back(seg(last, root));

        }

        last = root;

        neg ^=1;

    }

    if (neg) res.push_back(seg(last, INF));

    return res;

}

vectnegpoly(vect co) {

    for (int i = 0; i < co.size(); ++i) {

        //co *= -1;

        co.d *= -1;

    }

    return co;

}

vector<seg> solvege(vect co) {

    return solvele(negpoly(co));

}

vectmultipoly(vect co1, vect co2) {

    vect res;

    for (int i = 0; i<co1.size()+co2.size()-1;++i) res.push_back(0);

    for (int i = 0; i < co1.size(); ++i) {

        for (int j = 0; j < co2.size(); ++j) {

            //res[i+j] += co1 * co2[j];

            res.d[i+j] += co1 * co2[j];

        }

    }

    return res;

}

vectsumpoly(vect co1, vect co2) {

    vect res;

    for (int i = 0; i < co1.size() || i < co2.size(); ++i) {

        double ta, tb;

        if (i < co1.size()) ta = co1; else ta = 0;

        if (i < co2.size()) tb = co2; else tb = 0;

        res.push_back(ta + tb);

    }

    return res;

}

// ax +b

vect lco(double a, double b) {

    vect res;

    res.push_back(b);

    res.push_back(a);

    return res;

}

segintersect(seg a, seg b) {

    seg res;

    if (a.l > b.l) res.l = a.l; else res.l = b.l;

    if (a.r < b.r) res.r = a.r; else res.r = b.r;

    return res;

}

void fix_time_limit(seg &ret) {

    if (ret.l < 0) ret.l = 0;

    if (ret.r > e) ret.r = e;

}

void calc_valid_seg( int a, int b, int c) {

    double x1 = x[b] - x[a], x2 = x[c] - x[a];

    double y1 = y[b] - y[a], y2 = y[c] - y[a];

    double vx1 = vx[b] - vx[a], vx2 = vx[c] - vx[a];

    double vy1 = vy[b] - vy[a], vy2 = vy[c] - vy[a];

    vect tmp;   

    vect co = sumpoly(multipoly(lco(vy1, y1), lco(vx2, x2)),negpoly(multipoly(lco(vy2, y2), lco(vx1, x1))));

    //  int start2 =clock();

    vector<seg> seg1 = solvege(co);

    tmp.push_back(-4*r*r);

    vect co2 = sumpoly(multipoly(co, co),multipoly(sumpoly(multipoly(lco(vx1,x1), lco(vx1,x1)), multipoly(lco(vx1,x1),lco(vx1,x1))), tmp));

    vector<seg> seg2 = solvele(co2);

    for (int i = 0; i < seg1.size(); ++i) {

        fix_time_limit(seg1);

        if (seg1.empty()) continue;

        for (int j = 0; j < seg2.size(); ++j){

            fix_time_limit(seg2[j]);

            if (seg2[j].empty()) continue;

            seg tmpseg = intersect(seg1,seg2[j]);

            if (!tmpseg.empty())add_seg(tmpseg);

        }

    }

}

int cmp(const point & l, const point & r) {

    if (l.pos != r.pos) return l.pos < r.pos;

    return l.mrk < r.mrk;

}

void sort_seg() {

    sort(points, points + pn, cmp);

}

void calc_max() {

    int i, tot = 0, maxx = 0;

    for (i = 0; i < pn; ++i) {

        if (points.mrk) ++tot; else --tot;

        if (tot > maxx) maxx = tot;

    }

    if (maxx + 2 > best) best = maxx + 2;

}

int main() {

    //freopen("D:\\Contest\\4+2\\2011\\ninjapro\\input1.txt","r", stdin);

    int t, i, j, k;

    scanf("%d", &t);

    while (t--) {

        scanf("%d%d%d", &n, &r, &e);

        for (i=0;i<n;++i) {

            scanf("%lf%lf%lf%lf", &x, &y,&vx, &vy);

        }

        best = 2;

        // enumerate the fixed line

        for (i=0;i<n;++i) {

            for (j=0;j<n;++j) {

                pn = 0;

                if (i == j) continue;

                for (k = 0; k < n; ++k) {

                    if (i == k || j == k) continue;

                    calc_valid_seg(i, j, k);

                }

                sort_seg();

                calc_max();

                //printf("%d %d\n", i,j);

            }

        }

        printf("%d\n", best);

    }

    return 0;

}