WC2010 能量场

问题简述

能量场中有N个粒子,每个粒子都有一个质量m和结合系数c,两个粒子a,b合并时会产生mamb(ca - cb)的能量。(1)找出两个粒子相结合,使得产生的能量最大。(2)从中找出任意k个粒子排列成一个环,相邻两个粒子分别合并,使得总能量最大,产生负能量的粒子必须是环上连续的一段。

问题分析

乍一看这个问题的第一问好像很简单,枚举每对粒子即可,但是时间复杂度是O(N2)的,而且难以想到如何优化。第二问则更加困难,搜索、动态规划是肯定不行的,贪心、图论也难以找到建模方式。分析发现这个问题可以归约到一个0/1规划问题,如果没有特殊性将无法解决,而特殊性无非在于能量产生的公式,因此将目光聚焦到这个公式上,对公式进行一些变形:

  • mamb(ca - cb)
  • = mambca - mambcb
  • = macamb - mbcbma
  • 设xi=mici,yi=mi则原式
  • = xayb - xbya

经过变形,我们可以明显地看出公式变形为了两个向量叉积的公式,这给我们以启发:把每个粒子看做平面上的一个点,两个粒子合并产生的能量就是原点指向两个点的向量的叉积。因此问题的第一问就转化为了:找到两个向量,使它们的叉积最大。而第二问找到一个环合并的公式恰好对应了多边形的面积公式,再加上“产生负能量的粒子必须是环上连续的一段”这个限制,这个多边形必须是简单多边形[1]

要使第二问要求的环对应的多边形面积尽量大,应是平面上这些点能组成的最大的简单多边形,那么就应该是这些点的凸包。

相较之下第一问反而更难求解,不过有了对应的几何意义,就容易下手了。两个向量的叉积对应了两个向量所夹的平行四边形的有向面积,要使之最大首先应该是正值,即让第一个向量在第二个向量顺时针方向。当我们确定了第一个向量[2],即确定了平行四边形的一个底边,要使面积最大,应最大化平行四边形的高。于是我们可以做一条与该向量平行的直线,不断向上平移,直到遇到距离最远的点为止,这样的高最大。第一个向量与这个最远的点对应的向量做叉积就是对应的最大面积,很显然这个最远的点一定在凸包上,反过来考虑第一个向量的终点也一定在凸包上,因此查找这对向量时只需考虑凸包上的点。

image

有这个性质以后,如果直接枚举这对顶点,可能会快不少,但时间复杂度依然是O(N2)的。这时如果我们以某种特定的顺序枚举第一个向量,可以减少不少枚举量。具体方法是将从原点到凸包上所有的点的向量按照逆时针方向排序,按顺序枚举,这时候我们枚举的向量就是逆时针方向移动的,对应的第二个向量的终点在凸包上也是逆时针顺序移动(从最上点到最左点),因此枚举就是均摊线性的时间复杂度了。

算法描述

  1. 将所有粒子抽象为平面第一象限内的一个点(mici,mi)。

  2. 求平面上点的凸包。

  3. 把凸包上的点按照向量极角的顺序排序依次枚举作为第一个向量i。

  4. 找到对应第二个向量的终点j,应满足向量<j,j+1>在向量i逆时针方向。

复杂度分析

求凸包的时间复杂度为O(NlogN),枚举最优向量对的时间复杂度为O(N),因此总体时间复杂度为O(NlogN)。

参考程序

/* 
 * Problem: NOI Winter Camp 2010 efield
 * Author: Guo Jiabao
 * Time: 2010.3.15 12:21
 * Label: Solved
 * Memo: Computing Geometry + Convex Hull + Graham Scan
*/
#include <cstdio>
#include <cstdlib>
#include <cmath>
#include <cstring>
#include <iostream>
#include <sstream>
#include <vector>
#include <list>
#include <deque>
#include <string>
#include <queue>

#define var(a,b) typeof(b) a(b)
#define foreach(a,b) for (var(a,b.begin());a!=b.end();++a)
using std::sort;

const int MAXN = 50002;

struct point
{
    double x,y;
    int id;
};

int N,ipole;
point P[MAXN],pole,convex[MAXN],A[MAXN];

point operator -(const point &a,const point &b)
{
    point t={a.x-b.x,a.y-b.y,0};
    return t;
}

double operator *(const point &a,const point &b)
{
    return a.x*b.y - b.x*a.y;
}

void init()
{
    freopen("efield.in","r",stdin);
    freopen("efield.out","w",stdout);
    scanf("%d",&N);
    ipole=0;
    for (int i=1;i<=N;i++)
    {
        double m,c;
        scanf("%lf%lf",&m,&c);
        P[i].x = m*c;
        P[i].y = m;
        P[i].id = i;
        if (P[i].y > P[ipole].y)
            ipole = i;
    }
    pole = P[ipole];
    P[ipole] = P[N--];
}

bool cmp(const point &a,const point &b)
{
    point p = a - pole;
    point q = b - pole;
    double fc = p * q;
    if (fc > 0)
        return true;
    if (fc < 0)
        return false;
    return (p.x * p.x + p.y * p.y > q.x * q.x + q.y * q.y);
}

void graham()
{
    int top;
    sort(&P[1],&P[N+1],cmp);
    convex[1] = pole;
    convex[top=2] = P[1];
    for (int i=2;i<=N;i++)
    {
        point p = convex[top] - convex[top-1];
        point q = P[i] - convex[top];
        if (p*q < 0)
        {
            top--;
            i--;
        }
        else
            convex[++top] = P[i];
    }
    N = top;
    convex[N+1] = convex[1];
}

void getSumArea()
{
    int i;
    printf("%d\n",N);
    for (i=1;i<N;i++)
        printf("%d ",convex[i].id);
    printf("%d\n",convex[i].id);
/*    double SumArea = 0;
    for (i=1;i<=N;i++)
        SumArea += convex[i] * convex[i+1];
    SumArea *= 0.5;*/
}

bool cmp2(const point &a,const point &b)
{
    return a * b > 0;
}

void getPair()
{
    int M = 1,i,j,ai,aj;
    double PairArea = 0;
    for (i=1;i<=N;i++)
    {
        A[i] = convex[i];
        if (convex[i].x < convex[i-1].x)
            M = i;
    }
    sort(&A[1],&A[N+1],cmp2);
    for (i=1,j=1;i<=N;i++)
    {
        while (j<M && A[i]*(convex[j]-convex[j+1]) <= 0)
            j++;
        if (j == M)
            break;
        double cur = A[i] * convex[j];
        if (cur > PairArea)
        {
            PairArea = cur;
            ai = A[i].id;
            aj = convex[j].id;
        }
    }
    printf("%d %d\n",ai,aj);
}

void solve()
{
    graham();
    getPair();
    getSumArea();
}

int main()
{
    init();
    solve();
    return 0;
}

</>[1] 如果不是简单多边形,则面积公式一定会有多段连续的负值。

[2] 即顺时针方向下方的那个向量,下同。

附录

本题答案不唯一,因此需要一个Special Judge,我写了一个Cena的。使用方法就是编译后放进数据文件夹,添加题目的时候设置自定义校验器,如果实在不会看帮助吧。 下载地址:efield-check-cena.zip image

相关日志