DYNHUL - Editorial

PROBLEM LINK:

Practice
Div-2 Contest
Div-1 Contest

Author: Alei Reyes
Tester: Encho Mishinev

DIFFICULTY:

TIE-BREAK

PREREQUISITES:

Convex Hull

PROBLEM:

We are given N points in a plane, the i-th with coordinates (X_i,Y_i). After a point is removed the penalty increases by twice the area of the convex hull of the remaining points. Find the optimal order to remove all points while the penalty is minimized

EXPLANATION:

Given a set of points S, the greedy algorithm always removes the point P \in S that minimizes the area of convex\_hull(S-\{ P \}).

To find such P we could try all points in S, but that takes O(N^3 log(N)) and exceeds the time limit. However the point P must belong to the perimeter of the convex hull, because if we remove a point that is in the interior the area of the new convex hull will not change. Note that since the points are generated randomly, the expected number of points in the perimeter of the convex hull is small.

Let P_1, P_2, ..., P_M be the points in the perimeter of the convex hull in clockwise order, sometimes is better to remove the point P_x that maximizes the area of the triangle T with vertices P_{x-1}, P_x, P_{x+1}. Let’s call this the triangle heuristic, the idea can be improved further if we also take into account the number of points inside T.

We can use random search and keep flipping a coint to decide if we should use the greedy or the triangle heuristic when calculating the next point to remove.

When the number of points becomes small we can use the following dynamic programming: f(S) = area(convex\_hull(S) ) + min_{P \in S}(f(S-\{P\})).

Ryodan heuristics

  • Maximum collinear points

The convex hull of a set of collinear points have zero area, thea idea is to find the line L with maximum number of points, apply the greedy (or the triangle heuristic) on the set S-L, and then remove the points that belongs to L.

  • Look backwards

If we read the final permutation in reverse order, the process becomes starting from an empty set keep adding the points in some order. We can start from an arbitrary point and keep greedily adding the point that minimizes the new convex hull area.

SOLUTIONS:

Setter's Solution
#include<bits/stdc++.h>
using namespace std;
typedef long long int uli;
vector<vector<pair<int,int> >>all;


int rint(char nxt){
  char ch=getchar();
  int v=0;
  int sgn=1;
  if(ch=='-')sgn=-1;  
  else{
    assert('0'<=ch&&ch<='9');
    v=ch-'0';
  }
  while(true){
    ch=getchar();
    if('0'<=ch && ch<='9')v=v*10+ch-'0';
    else{
      assert(ch==nxt);
      break;
    }
  }
  return v*sgn;
}

#define REMOVE_REDUNDANT

typedef double T;
const T EPS = 1e-7;
struct PT { 
  T x, y; 
  int idx;
  PT() {} 
  PT(T x, T y) : x(x), y(y) {}
  bool operator<(const PT &rhs) const { return make_pair(y,x) < make_pair(rhs.y,rhs.x); }
  bool operator==(const PT &rhs) const { return make_pair(y,x) == make_pair(rhs.y,rhs.x); }
  PT operator -(const PT p){
    PT ans={x-p.x,y-p.y};
    return ans;
  }
};

T cross(PT p, PT q) { return p.x*q.y-p.y*q.x; }
T area2(PT a, PT b, PT c) { return cross(a,b) + cross(b,c) + cross(c,a); }

#ifdef REMOVE_REDUNDANT
bool between(const PT &a, const PT &b, const PT &c) {
  return (fabs(area2(a,b,c)) < EPS && (a.x-b.x)*(c.x-b.x) <= 0 && (a.y-b.y)*(c.y-b.y) <= 0);
}
#endif

void ConvexHull(vector<PT> &pts) {
  sort(pts.begin(), pts.end());
  pts.erase(unique(pts.begin(), pts.end()), pts.end());
  vector<PT> up, dn;
  for (int i = 0; i < (int)pts.size(); i++) {
    while (up.size() > 1 && area2(up[up.size()-2], up.back(), pts[i]) >= 0) up.pop_back();
    while (dn.size() > 1 && area2(dn[dn.size()-2], dn.back(), pts[i]) <= 0) dn.pop_back();
    up.push_back(pts[i]);
    dn.push_back(pts[i]);
  }
  pts = dn;
  for (int i = (int) up.size() - 2; i >= 1; i--) pts.push_back(up[i]);

#ifdef REMOVE_REDUNDANT
  if (pts.size() <= 2) return;
  dn.clear();
  dn.push_back(pts[0]);
  dn.push_back(pts[1]);
  for (int i = 2; i < (int)pts.size(); i++) {
    if (between(dn[dn.size()-2], dn[dn.size()-1], pts[i])) dn.pop_back();
    dn.push_back(pts[i]);
  }
  if (dn.size() >= 3 && between(dn.back(), dn[0], dn[1])) {
    dn[0] = dn.back();
    dn.pop_back();
  }
  pts = dn;
#endif
}
bool lcmp(PT a, PT b){
  if(a.x!=b.x)return a.x<b.x;
  return a.y<b.y;
}
bool rem[2123];


vector<PT>filter(vector<pair<int,int> >p){
  vector<PT>all;
  for(int j=0;j<int(p.size());j++)if(!rem[j]){
    PT a;
    a.x=p[j].first;
    a.y=p[j].second;
    a.idx=j;
    all.push_back(a);
  }
  return all;
}
vector<PT> getHull(vector<pair<int,int> >p){
  auto all=filter(p);
  ConvexHull(all);
  return all;
}

int removePoint(int idx,vector<pair<int,int> >p){
  rem[idx]=true;
  auto nhull=getHull(p);
  uli area=0;
  for(int k=2;k<int(nhull.size());k++){
    area+=cross(nhull[k]-nhull[0], nhull[k-1]-nhull[0]);
  }
  if(area<0)area=-area;
  rem[idx]=false;
  return area;
}
uli score;
uli bestScore=-1;
vector<int>greedy(vector<pair<int,int> >p,double prob1, double prob2){
  vector<int>ans;
  int n=p.size();
  for(int i=0;i<n;i++)rem[i]=false;  
  score=0;
  for(int i=0;i<n;i++){
    vector<PT>h=getHull(p);
    sort(h.begin(),h.end(),lcmp);
    int r=rand()%100;
    int idx=-1;
    uli best=0;
    if(r<prob1*100.0){
      idx=rand()%h.size();      
    }
    else if(r<=prob2*100.0){
      for(int j=0;j<int(h.size());j++){
        uli area=removePoint(h[j].idx,p);
        if(idx==-1 || area<best){
          idx=j;
          best=area;
        }
      }
    }
    else{
      int sz=h.size();
      int maxi=0;
      for(int j=0;j<sz;j++){
        int bef=(j-1+sz)%sz;
        int nxt=(j+1)%sz;
        int tri=fabs(cross(h[nxt]-h[i],h[i]-h[bef]));
        if(idx==-1 || tri>maxi){
          maxi=tri;
          idx=j;
        }
      }      
    }
    best=removePoint(h[idx].idx,p);
    ans.push_back(h[idx].idx);
    rem[h[idx].idx]=true;
    score+=best;
  }
  return ans;
}
clock_t tStart;
double getTime(){
  return (double)(clock() - tStart)/CLOCKS_PER_SEC;
}
int main(){
  tStart = clock();  
  int n=rint('\n');
  // assert(n==1024);
  vector<pair<int,int> >all(n);
  for(int i=0;i<n;i++){
    all[i].first=rint(' ');
    all[i].second=rint('\n');
  }
  uli g=rint('\n');
  vector<int>ans;
  double prob1=0.0;
  double prob2=1.0;
  while(getTime()<4){
    vector<int>bet=greedy(all,prob1,prob2);
    if(bestScore==-1 || score<bestScore){
      bestScore=score;
      ans=bet;
    }
    prob1=0.3;
    prob2=0.7;
  }
  cerr<<"best="<<g-bestScore<<endl;
  for(int x:ans){
    cout<<x+1<<" ";
  }
  cout<<endl;
  return 0;
}
Tester's Solution
  #include <iostream>
  #include <stdio.h>
  #include <vector>
  #include <assert.h>
  #include <algorithm>
  using namespace std;

  struct pt
  {
      int x, y;
      int id;

      pt(int X, int Y): x(X), y(Y) {}
      pt(): x(0), y(0) {}

      pair<int,int> asPair()
      {
          return make_pair(x, y);
      }

      pt operator-(const pt &other) const
      {
          return pt(x - other.x, y - other.y);
      }
  };

  int n;
  pt pts[2011];

  int sq(int a)
  {
      return a*a;
  }

  int sdist(pt A, pt B)
  {
      return sq(A.x-B.x) + sq(A.y-B.y);
  }

  bool LL(pt A, pt B)
  {
      return A.asPair() < B.asPair();
  }

  bool DD(pt A, pt B)
  {
      return make_pair(A.y, A.x) < make_pair(B.y, B.x);
  }

  int CP(pt A, pt B)
  {
      return A.x * B.y - A.y * B.x;
  }

  int CP(pt A, pt B, pt C)
  {
      return CP(B - A, C - A);
  }

  pt pivot;
  bool SP(pt A, pt B)
  {
      if (A.id == pivot.id)
          return true;
      else if (B.id == pivot.id)
          return false;
      else
      {
          int cp = CP(pivot, A, B);

          if (cp != 0)
              return cp > 0;
          else
              return sdist(pivot, A) < sdist(pivot, B);
      }
  }

  vector<pt> hull;
  vector<pt> active;

  bool alive[2011];

  void showHull()
  {
      int i,j;

      fprintf(stderr, "[");
      for (i=1;i<hull.size();i++)
      {
          fprintf(stderr, "[(%d, %d), (%d, %d)],", hull[i-1].x, hull[i-1].y, hull[i].x, hull[i].y);
      }
      fprintf(stderr, "]\n");
  }

  void assertHull()
  {
      int i,j;

      for (i=1;i<=n;i++)
      {
          if (!alive[ pts[i].id ])
              continue;

          for (j=1;j<hull.size();j++)
          {
              if (CP(hull[j-1], hull[j], pts[i]) < 0)
              {
                  printf("(%d, %d) when looking at line (%d, %d)-(%d, %d)\n", pts[i].x, pts[i].y, hull[j-1].x, hull[j-1].y, hull[j].x, hull[j].y);
                  printf("Gives %d\n", CP(hull[j-1], hull[j], pts[i]));
                  printf("Id is %d\n", pts[i].id);
              }
              assert(CP(hull[j-1], hull[j], pts[i]) >= 0);
          }
          assert(CP(hull.back(), hull[0], pts[i]) >= 0);
      }
  }

  int getHull()
  {
      int i,j;

      hull.clear();
      active.clear();
      for (i=1;i<=n;i++)
      {
          if (alive[ pts[i].id ])
          {
              active.push_back(pts[i]);
          }
      }

      if (active.empty())
          return 0;

      sort(active.begin(), active.end(), DD);
      pivot = active[0];
      sort(active.begin(), active.end(), SP);

      assert(active[0].id == pivot.id);

      //printf("Try\n");
      for (i=0;i<active.size();i++)
      {
          while(hull.size() > 1 && CP(hull[ (int)hull.size() - 2 ], hull[ (int)hull.size() - 1 ], active[i]) <= 0)
              hull.pop_back();

          hull.push_back(active[i]);
      }
      while(hull.size() > 1 && CP(hull[ (int)hull.size() - 2 ], hull[ (int)hull.size() - 1 ], pivot) <= 0)
          hull.pop_back();

      int area = 0;
      for (i=1;i<hull.size();i++)
      {
          area += CP(hull[i-1], hull[i]);
      }
      area += CP(hull.back(), hull[0]);

      /*assert(CP(hull[ (int)hull.size()-2 ], hull.back(), hull[0]) >= 0);
      assert(CP(hull.back(), hull[0], hull[1]));
      assert(hull[0].id == pivot.id);
      for (i=2;i<hull.size();i++)
      {
          assert(CP(hull[i-2], hull[i-1], hull[i]) >= 0);
      }
      */
      //assertHull();

      return abs(area);
  }

  int main()
  {
      //freopen("1.in.txt", "r", stdin);
      //freopen("res.txt", "w", stdout);

      int i,j;

      scanf("%d", &n);

      for (i=1;i<=n;i++)
      {
          alive[i] = true;

          scanf("%d %d", &pts[i].x, &pts[i].y);
          pts[i].id = i;
      }

      vector<int> seq;
      int ans = 0;
      int iter = 0;
      while(getHull() > 0)
      {
          iter++;

          int best = -1, bestScore = -1;
          pt bp;

          vector<pt> cands;
          for (i=0;i<hull.size();i++)
          {
              cands.push_back(hull[i]);
          }
          sort(cands.begin(), cands.end(), LL);

          for (i=0;i<cands.size();i++)
          {
              if (!alive[ cands[i].id ])
                  continue;

              alive[ cands[i].id ] = false;
              int score = getHull();
              alive[ cands[i].id ] = true;

              if (best == -1 || score < bestScore)
              {
                  best = cands[i].id;
                  bestScore = score;
                  bp = cands[i];
              }
          }

          alive[best] = false;
          seq.push_back(best);
          ans += bestScore;

          //fprintf(stderr, "Choosing %d with %d\n", best, bestScore);

          /*if (best == 236)
          {
              getHull();
              showHull();

              for (i=0;i<hull.size();i++)
              {
                  if (hull[i].asPair() == make_pair(88, 67))
                  {
                      fprintf(stderr, "ITS IN THE HULL\n");
                  }
              }

              for (i=1;i<=n;i++)
              {
                  if (pts[i].asPair() == make_pair(88, 67))
                  {
                      fprintf(stderr, "EXISTS THOUGH\n");
                  }
              }
              exit(1);
          }*/
      }

      for (i=1;i<=n;i++)
      {
          if (alive[i])
              seq.push_back(i);
      }

      for (i=0;i<seq.size();i++)
      {
          if (i != 0)
              printf(" ");

          printf("%d", seq[i]);
      }
      printf("\n");

      int realAns;
      scanf("%d", &realAns);
      fprintf(stderr, "Score: %d vs %d\n", ans, realAns);

      return 0;
  }

Ryodan solution

3 Likes

I tried the approach mentioned in this paper but could not make it due to WA (may be because of not checking collinear points correctly).
http://morpheo.inrialpes.fr/~wuhrer/data/uploads/publications/outliers-cccg.pdf

@alei @enchom Do you think we can apply the similar mentioned in paper especially for case of removing 1 outlier?

maybe using that paper to find the best c points to remove N/c times, for some small c (c<=4), but I think heuristics or beam-search like used by other contestants works better.

for c=1, you just have to try the points in the convex hull.

1 Like

what a clean code from Ryodan!

Yeah I tried with c=1 and in each step I calculated triangle area considering a point on the convex hull and adjacent vertices in the same hull. Then calculated exposed vertices from second hull (excluding first hull points) which lie inside the triangle and make another convex hull from these points including earlier two adjacent vertices. Subtract this area from the triangle we calculated earlier. Then we can subtract this remaining area from first convex hull total area to check for minimum remaining convex hull area.