LCA HL分解

heavy light decompositionを試したのでLCAの問題にsubmit。
参考にした資料
http://math314.hateblo.jp/entry/2014/06/24/220107

問題
http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=GRL_5_C

#include<iostream>
#include<vector>
#include<algorithm>

#define REP(i, n) for (int i=0;i<int(n);++i)
#define EACH(i,c) for(__typeof((c).begin())i=(c).begin();i!=(c).end();++i)
#define ALL(A) A.begin(), A.end()

using namespace std;

struct HeavyLight {
  int pathCount,n;
  vector<int>size,parent,in,out,path,pathRoot;
  vector<vector<int> > tree;

  HeavyLight(vector<vector<int> > t)
    :pathCount(0),n(t.size()),size(n),parent(n),in(n),out(n),
     path(n),pathRoot(n),tree(t){
    int time=0;
    dfs(0,-1,time);
    buildPaths(0,newPath(0));
  }

  void dfs(int u, int p, int &k){
    in[u]=k++, parent[u]=p, size[u]=1;
    EACH(v,tree[u])if(*v!=p)dfs(*v, u, k),size[u]+=size[*v];
    out[u]=k++;
  }

  int newPath(int u){ pathRoot[pathCount]=u; return pathCount++; }

  void buildPaths(int u, int pt){
    path[u]=pt;
    EACH(v,tree[u])if(*v != parent[u])
      buildPaths(*v, 2*size[*v] >= size[u] ? pt : newPath(*v));
  }

  bool isAncestor(int p, int ch){
    return in[p] <= in[ch] && out[ch] <= out[p];
  }

  int lca(int a,int b){
    for(int root; !isAncestor(root=pathRoot[path[a]],b);a=parent[root]);
    for(int root; !isAncestor(root=pathRoot[path[b]],a);b=parent[root]);
    return isAncestor(a,b)?a:b;
  }
};

int main(void){
  ios::sync_with_stdio(false);

  int n;
  cin >> n;
  vector<vector<int> >tree(n);
  
  REP(i,n){
    int k;
    cin >> k;
    REP(j,k){
      int ch;
      cin >> ch;
      tree[i].push_back(ch);
    }
  }

  HeavyLight hl = HeavyLight(tree);

  int q;
  cin >> q;
  while(q--){
    int u,v;
    cin >> u >> v;
    cout << hl.lca(u,v) << "\n";
  }
  
  return 0;
}

POJ 2217: Secretary

問題文
http://poj.org/problem?id=2217

Suffix Treeを実装したので最長共通部分文字列の問題にsubmit。

#include <iostream>
#include <vector>
#include <algorithm>
#include <string>
#include <cctype>

using namespace std;

const int SIZE = 128;

struct Node {
	int b,e,d; //begin, end, distance from root
	Node *parent, *children[128], *suffixLink;

        Node():b(0),e(0),d(0),parent(NULL),suffixLink(NULL){
          fill(children, children + SIZE, (Node *)0);
        }
	Node(int b, int e, int d, Node *p):b(b),e(e),d(d),parent(p),suffixLink(NULL){
          fill(children, children + SIZE, (Node *)0);
	}
};

int a[30000];
Node *buildSuffixTree(string s){
  int n = s.size();
  for(int i = 0; i < n; i++)a[i] = (int)s[i];

  Node *root = new Node();
  Node *node = root;
  for(int i = 0, tail = 0; i < n; i++, tail++){
    Node *last = NULL;
    while(tail >= 0){
      Node *ch = node->children[a[i - tail]];
      while(ch != NULL && tail >= ch->e - ch->b){
        tail -= ch->e - ch->b;
        node = ch;
        ch = ch->children[a[i - tail]];
      }

      if(ch == NULL){
        node->children[a[i]] = new Node(i, n, node->d + node->e - node->b, node);
        if(last != NULL)last->suffixLink = node;
        last = NULL;
      } else {
        int afterTail = a[ch->b + tail];
        if(afterTail == a[i]){
          if(last != NULL)last->suffixLink = node;
          break;
        } else {
          Node *splitNode = new Node(ch->b, ch->b + tail, node->d + node->e - node->b, node);
          splitNode->children[a[i]] = new Node(i, n, ch->d + tail, splitNode);
          splitNode->children[afterTail] = ch;
          ch->b += tail;
          ch->d += tail;
          ch->parent = splitNode;
          node->children[a[i - tail]] = splitNode;
          if(last != NULL)last->suffixLink = splitNode;
          last = splitNode;
        }
      }
      if(node == root)tail--;
      else node = node->suffixLink;
    }
  }
  return root;
}

int lcsLength;
int lcsBeginIndex;
//longest common substring
int lcs_traverse(Node *node, int i1, int i2){
  if(node->b <= i1 && i1 < node->e)return 1;
  if(node->b <= i2 && i2 < node->e)return 2;

  int visited = 0;
  for(int f = 0; f < SIZE; f++)
    if(node->children[f] != NULL)
      visited |= lcs_traverse(node->children[f], i1, i2);

  if(visited == 3){
    int curLength = node->d + node->e - node->b;
    if(lcsLength < curLength){
      lcsLength = curLength;
      lcsBeginIndex = node->b - node->d;
    }
  }
  return visited;
}

string lcs(string s1, string s2){
  string s = s1 + '\1' + s2 + '\2';
  Node *tree = buildSuffixTree(s);

  lcsLength = 0;
  lcsBeginIndex = 0;
  lcs_traverse(tree, s1.size(), s1.size() + s2.size() + 1);
  return s.substr(lcsBeginIndex, lcsLength);
}

int main(void){
  int n;
  cin >> n;
  cin.ignore();
  while(n--){
    string s1,s2;
    getline(cin, s1);
    getline(cin, s2);
    lcs(s1, s2);
    cout << "Nejdelsi spolecny retezec ma delku " << lcsLength << "." << endl;
  }
  return 0;
}

Dice4-Scala

問題文
http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=ITP1_11_D

最近Scalaを学び始めた。詳しいことは後で。

import scala.io.StdIn._

object Main {
  case class Rotation(run: Dice => Dice){
    def apply(d: Dice): Dice = run(d)
    def compose(other: Rotation) = Rotation(run compose other.run)
  }

  object Rotation {
    def id = Rotation(d => d)
    def north = Rotation(d => Dice(List(1,5,2,3,0,4).map(d.pip)))
    def east = Rotation(d => Dice(List(3,1,0,5,4,2).map(d.pip)))
    def ccw = Rotation(d => Dice(List(0,3,1,4,2,5).map(d.pip)))

    final def axial(rot: Rotation): Seq[Rotation] =
      for(i <- 0 to 3) yield List.fill(i)(rot).fold(id)(_ compose _)

    def enumerateTop: Seq[Rotation] =
      for(i <- 1 to 6) yield ((0 to i).map(j => if(j % 2 == 0) north else east)).fold(id)(_ compose _)

    def enumerate: Seq[Rotation] = enumerateTop.flatMap(e => axial(ccw).map(e compose _))
  }

  case class Dice(pip: List[Int]){
    lazy val all: Seq[Dice] = Rotation.enumerate.map(_.run(this))
    def ===(that: Dice): Boolean = pip.corresponds(that.pip)(_ == _)
    def ==(that: Dice): Boolean = all.exists(_ === that)
  }

  def main(args: Array[String]): Unit = {
    val Array(n,_*) = readLine.split(" ").map(_.toInt)
    val dices = for(i <- 1 to n) yield Dice(readLine.split(" ").map(_.toInt).toList)

    val res = for {
      i <- 0 until n
      j <- i+1 until n
    } yield dices(i) == dices(j)

    println(if(res.forall(_ == false)) "Yes" else "No")
  }
}

LCA LinkCutTree

問題文
http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=GRL_5_C&lang=jp

参考スライド
http://www.slideshare.net/iwiwi/2-12188845

Link Cut TreeでLCAをやってみました。
Link Cut Tree部分はほとんど参考スライドのものと同じです。

#include<iostream>
#include<vector>
#include<algorithm>
#include<queue>
 
#define INF (1<<29)
#define max_n 10000
 
using namespace std;
 
struct node_t{
 
  node_t *pp,*lp,*rp;
  int id,val,mini,minId,lazy;
 
  node_t(int id,int v):id(id),val(v){
    pp=lp=rp=NULL; lazy=0;
  }
 
  bool is_root(){
    return !pp || (pp->lp != this && pp->rp != this);
  }
   
  void rotr(){
    node_t *q=pp,*r=q->pp;
    if((q->lp=rp))rp->pp=q;
    rp=q;q->pp=this;
    if((pp=r)){
      if(r->lp==q)r->lp=this;
      if(r->rp==q)r->rp=this;
    }
  }
 
  void rotl(){
    node_t *q=pp,*r=q->pp;
    if((q->rp=lp))lp->pp=q;
    lp=q;q->pp=this;
    if((pp=r)){
      if(r->lp==q)r->lp=this;
      if(r->rp==q)r->rp=this;
    }
  }
 
  void splay(){
    while(!is_root()){
      node_t *q=pp;
      if(q->is_root()){
    if(q->lp==this)rotr();
    else rotl();
      } else {
    node_t *r=q->pp;
    if(r->lp==q){
      if(q->lp==this){q->rotr();rotr();}
      else {rotl();rotr();}
    } else {
      if(q->rp==this){q->rotl();rotl();}
      else {rotr();rotl();}
    }
      }
    }
  }
};
 
node_t *expose(node_t *x){
  node_t *rp=NULL;
  for(node_t *p=x;p;p=p->pp){
    p->splay();
    p->rp=rp;
    rp=p;
  }
  x->splay();
  return rp;
}
 
node_t *find_root(node_t *x){
  expose(x);
  while(x->lp)x=x->lp;
  return x;
}
 
void link(node_t *c,node_t *p){
  expose(c);
  expose(p);
  c->pp=p;
  p->rp=c;
}
 
node_t *lowestCommonAncestor(node_t *x,node_t* y){
  expose(x);
  return expose(y);
}
 
node_t *node[100001];
 
int main(void){
  int n;
  cin >> n;
   
  for(int i=0;i<n;i++)node[i]=new node_t(i,0);
  for(int i=0;i<n;i++){
    int k,t;
    cin >> k;
    for(int j=0;j<k;j++){
      cin >> t;
      link(node[t],node[i]);
    }
  }
 
  int q;
  cin >> q;
  for(int i=0;i<q;i++){
    int u,v;
    cin >> u >> v;
    cout << lowestCommonAncestor(node[u],node[v])->id << endl;
  }
   
  return 0;
}

POJ:3580 SuperMemo

問題文
http://poj.org/problem?id=3580

前Treapで通した問題。今回はRBST。

#include<iostream>
#include<vector>
#include<algorithm>
#include<ctime>
#include<cstdlib>

#define INF (1<<29)

using namespace std;

typedef long long ll;

template<class T>
struct RBST{
public:
  struct node_t{
    T val,mini,lazy;
    int size;
    node_t *lch,*rch;
    bool rev;
    
    node_t(int v):val(v),mini(v){
      size=1;
      lazy=0;
      lch=rch=0;
      rev=false;
    }
  };

  node_t *root;  
  RBST():root(0){srand(time(NULL));}

private:

  int val(node_t *t){return t->val+t->lazy;}
  int size(node_t *x){return !x?0:x->size;}
  int mini(node_t *t){return !t?INF:t->mini+t->lazy;}
  
  node_t *update(node_t *t){
    
    if(!t)return t;
    
    t->size=size(t->lch)+size(t->rch)+1;
    t->mini=min(min(mini(t->lch),mini(t->rch)),val(t));
    
    if(t->lazy){
      t->val+=t->lazy;
      if(t->lch)t->lch->lazy+=t->lazy;
      if(t->rch)t->rch->lazy+=t->lazy;
      t->lazy=0;
    }
    
    if(t->rev){
      swap(t->lch,t->rch);
      if(t->lch)t->lch->rev^=true;
      if(t->rch)t->rch->rev^=true;
      t->rev=false;
    }  
    return t;
  }
  
  node_t* merge(node_t *l,node_t *r){
    l=update(l),r=update(r);
    if(!l || !r)return l?l:r;
    
    int m=l->size,n=r->size;
    if(rand()%(m+n)<m){
      l->rch=merge(l->rch,r);
      return update(l);
    } else {
      r->lch=merge(l,r->lch);
      return update(r);
    }
  }
  
  typedef pair<node_t*,node_t*>pnn;
  
  pnn split(node_t *t,int k){
    if(!update(t))return pnn(0,0);
    if(k<=size(t->lch)){
      pnn tmp=split(t->lch,k);
      t->lch=tmp.second;
      return pnn(tmp.first,update(t));
    } else {
      pnn tmp=split(t->rch,k-size(t->lch)-1);
      t->rch=tmp.first;
      return pnn(update(t),tmp.second);
    }
  }
  
  node_t *insert(node_t *t,int k,int val){
    pair<node_t *,node_t *>s=split(t,k);
    t=merge(s.first, new node_t(val));
    t=merge(t,s.second);
    return update(t);
  }
  
  node_t *erase(node_t *t,int k){
    pair<node_t *,node_t *>s1,s2;
    s2=split(t,k+1);
    s1=split(s2.first,k);
    delete s1.second;
    return update((t=merge(s1.first,s2.second)));
  }
  
  node_t *find(node_t *t, int k){
    int c=size(t->lch);
    if(k<c)return find(t->lch, k);
    if(k==c)return t;
    return find(t->rch, k-c-1);
  }

public:
  
  void insert(int k,T val){ root=insert(root,k,val); }
  void erase(int k){ root=erase(root,k); }
  node_t *find(int k){ return find(root,k); }

  void add(int id,T val){
    node_t *a=find(id);
    T tmp=val(a);
    erase(id);
    insert(id,tmp+val);
  }

  void rangeAdd(int l,int r,T v){
    pair<node_t *,node_t *>s1,s2;
    s2=split(root,r);
    s1=split(s2.first,l);
    s1.second->lazy+=v;
    root=merge(s1.first,s1.second);
    root=merge(root,s2.second);
  }

  void update(int id,T val){
    erase(id);
    insert(id,val);
  }

  T rangeMinimumQuery(int l,int r){
    pair<node_t *,node_t *>s1,s2;
    s2=split(root,r);
    s1=split(s2.first,l);
    T res=mini(s1.second);
    root = merge(s1.first,merge(s1.second,s2.second));
    return res;
  }

  T rangeSumQuery(int l,int r){
    pair<node_t *,node_t *>s1,s2;
    s2=split(root,r);
    s1=split(s2.first,l);
    T res=sum(s1.second);
    root = merge(s1.first,merge(s1.second,s2.second));
    return res;
  }

  void reverse(int l,int r){
    pair<node_t *,node_t *>s1,s2;
    s2=split(root,r);
    s1=split(s2.first,l);
    node_t *a=s1.first,*b=s1.second,*c=s2.second;
    b->rev^=true;
    root=merge(a,b);
    root=merge(root,c);
  }

  void revolve(int l,int r,int t){
    t%=(r-l);
    pair<node_t*, node_t*>a,b,c;
    c=split(root, r);
    b=split(c.first, r-t);
    a=split(b.first, l);
    root=merge(a.first, b.second);
    root=merge(root, a.second);
    root=merge(root, c.second);
  }
};

  int main(void){
  
  int n,m,in;
  RBST<ll> t;
  
  scanf("%d",&n);
  for(int i=0;i<n;i++){
    scanf("%d",&in);
    t.insert(i,in);
  }
  
  cin >> m;
  
  char com[10];
  int x,y,T,D,P;

  for(int i=0;i<m;i++){
    scanf("%s",com);

    if(com[0]=='A'){
      scanf("%d %d %d",&x,&y,&D);
      t.rangeAdd(x-1,y,D);
    }
    else if(com[0]=='R' && com[3]=='E'){
      scanf("%d %d",&x,&y);
      t.reverse(x-1,y);
    }
    else if(com[0]=='R' && com[3]=='O'){
      scanf("%d %d %d",&x,&y,&T);
      t.revolve(x-1,y,T);
    }
    else if(com[0]=='I'){
      scanf("%d %d",&x,&P);
      t.insert(x,P);
    }
    else if(com[0]=='D'){
      scanf("%d",&x);
      t.erase(x-1);
    }
    else if(com[0]=='M'){
      scanf("%d %d",&x,&y);
      printf("%lld\n",t.rangeMinimumQuery(x-1,y));
    }
  }    
  return 0;
}

SPOJ 6044 最小包含円

問題文
http://www.spoj.com/problems/QCJ4/

最小包含円のライブラリ検証問題。
自分の実装ではならしO(n)になっているはず…

#include<cmath>
#include<algorithm>
#include<iostream>
#include<vector>
#include<climits>
#include<cfloat>
#include<cstdio>

using namespace std;

typedef double Real;

Real EPS = 1e-8;
const Real PI = acos(-1);

int sgn(Real a, Real b=0){return a<b-EPS?-1:a>b+EPS?1:0;}
Real sqr(Real a){return sqrt(max(a,(Real)0));}

struct Point{  
  Real add(Real a, Real b){
    if(abs(a+b) < EPS*(abs(a)+abs(b)))return 0;
    return a+b;
  }

  Real x, y;
  Point(){}
  Point(Real x,Real y) : x(x) , y(y){}

  Point operator + (Point p){return Point(add(x,p.x), add(y,p.y));}
  Point operator - (Point p){return Point(add(x,-p.x), add(y,-p.y));}
  Point operator * (Real d){return Point(x*d,y*d);}
  Point operator / (Real d){return Point(x/d,y/d);}
  bool operator == (Point p){return !sgn(dist(p));}
  //昇順
  bool operator < (Point p)const{return (p.x!=x)?x<p.x:y<p.y;}
  Real norm(){return sqr(x*x+y*y);}
  Real dist(Point a){return (*this-a).norm();}
  Real dot(Point a){return x*a.x+y*a.y;}
  Real cross(Point a){return x*a.y-y*a.x;}
};

struct Circle{
  Point p;
  Real r;
  Circle(){}
  Circle(Point p, Real r) : p(p) , r(r){}
  
  bool contain(Point a){
    return sgn((a-p).norm(),r)<=0;
  }
  
  Circle circumCircle(Point a,Point b){
    Point q=(a+b)/2;
    return Circle(q,(a-q).norm());
  }
  
  Circle circumscribedCircle(Point p, Point q, Point r){
    Point a=(q-p)*2,b=(r-p)*2;
    Point c(p.dot(p)-q.dot(q),p.dot(p)-r.dot(r));
    Circle res;
    res.p.x=a.y*c.y-b.y*c.x;
    res.p.y=b.x*c.x-a.x*c.y;
    res.p=res.p/a.cross(b);
    return Circle(res.p, p.dist(res.p));
  }

  Circle minEnclosingCircle(vector<Point>ps){
    if(ps.size()==0)return Circle(Point(0,0),0);
    if(ps.size()==1)return Circle(ps[0],0);

    Circle circle=circumCircle(ps[0],ps[1]);
    for(int i=2;i<ps.size();i++){
      if(!circle.contain(ps[i])){
	circle=circumCircle(ps[0],ps[i]);
	for(int j=1;j<i;j++){
	  if(!circle.contain(ps[j])){
	    circle=circumCircle(ps[j],ps[i]);
	    for(int k=0;k<j;k++){
	      if(!circle.contain(ps[k])){
		circle=circumscribedCircle(ps[i],ps[j],ps[k]);
	      }
	    }
	  }
	}
      }
    }
    return circle;
  }

};

int main(void){

  int n;
  cin >> n;
  vector<Point>ps(n);

  for(int i=0;i<n;i++)cin >> ps[i].x >> ps[i].y;
  Circle c;
  printf("%.2f",c.minEnclosingCircle(ps).r*2);
  cout << endl;

  return 0;
}

CGL_6/A: 線分交差

問題文
http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=CGL_6_A


書籍「プログラミングコンテスト攻略のためのアルゴリズムとデータ構造」で紹介されいるが、セグ木を使って効率よく解けると書いてあったので実装。
x座標を座標圧縮すると平面走査で使っている二分探索木の部分にセグ木を使えてO(n log n)になる……ということだと思う。

#include<iostream>
#include<vector>
#include<algorithm>
#include<map>

#define all(c) (c).begin(),(c).end()
#define x second
#define y first
#define a first
#define b second
#define PUSH 0
#define SUM 1
#define POP 2

using namespace std;

typedef pair<int,int> pii;
typedef pair<pii,int>piii;
typedef pair<pii,pii> segment;

const int MAX_N = 1<<17;
int n,dat[2*MAX_N];

void init(int n_){
  n=1;
  while(n<n_)n*=2;
  fill(dat,dat+2*MAX_N,0);
}

void update(int k,int a){
  k+=n-1;
  dat[k]=a;
  while(k>0){
    k=(k-1)/2;
    dat[k]=dat[k*2+1]+dat[k*2+2];
  }
}

int sum(int a,int b,int k=0,int l=0,int r=n){
  if(r<=a || b<=l)return 0;
  if(a<=l && r<=b)return dat[k];
  
  int m=(l+r)/2;
  int vl=sum(a,b,k*2+1,l,m);
  int vr=sum(a,b,k*2+2,m,r);
  return vl+vr;
}

map<int,int> xmp;

struct Q{
  int y,x1,x2,q;
  Q(){}
  Q(int y,int x1,int x2,int q):y(y),x1(x1),x2(x2),q(q){}
};

bool cmp(Q p,Q q){
  if(p.y!=q.y)return p.y<q.y;
  return p.q<q.q;
}

int main(void){

  int m;
  cin >> m;
  
  vector<int>xv;
  vector<segment>v(m);
  vector<Q>query;
  
  for(int i=0;i<m;i++){
    cin >> v[i].a.x >> v[i].a.y >> v[i].b.x >> v[i].b.y;

    if(v[i].a.x>v[i].b.x)swap(v[i].a.x,v[i].b.x);
    if(v[i].a.y>v[i].b.y)swap(v[i].a.y,v[i].b.y);

    xv.push_back(v[i].a.x);
    xv.push_back(v[i].b.x);
  }
  
  sort(all(xv));

  for(int i=0;i<m;i++){
    if(v[i].a.x==v[i].b.x){
      query.push_back(Q(v[i].a.y,v[i].a.x,-1,PUSH));
      query.push_back(Q(v[i].b.y,v[i].b.x,-1,POP));
    }
    if(v[i].a.y==v[i].b.y){
      query.push_back(Q(v[i].a.y,v[i].a.x,v[i].b.x,SUM));
    } 
  }
  
  sort(all(query),cmp);

  xv.erase(unique(all(xv)),xv.end());

  for(int i=0;i<xv.size();i++)xmp[xv[i]]=i;

  init(xv.size());
  
  int cnt=0;
  for(int i=0;i<query.size();i++){
    if(query[i].q==PUSH)update(xmp[query[i].x1],1);
    if(query[i].q==POP)update(xmp[query[i].x1],0);
    if(query[i].q==SUM)cnt+=sum(xmp[query[i].x1],xmp[query[i].x2]+1);
  }
  cout << cnt << endl;
  
  return 0;
}