「ライブラリ」の編集履歴（バックアップ）一覧に戻る

ライブラリ - (2006/05/23 (火) 02:35:14) の編集履歴(バックアップ)

---

平面幾何計算

class Double2D{
    double x;
    double y;
    static final double EPS = 1.0E-9;
    Double2D(double x,double y){
        this.x = x;
        this.y = y;
    }
    @Override
    public boolean equals(Object o){
        Double2D p = (Double2D)o;
        return x==p.x&&y==p.y;
    }
    public double squareDist(Double2D p){
        return (x-p.x)*(x-p.x) + (y-p.y)*(y-p.y);
    }
    public double dist(Double2D p){
        return Math.sqrt(this.squareDist(p));
    }
    // thisとpを通る半径rの円の中心
    public Double2D[] centerOfPassingCircle1(Double2D p,double r){
        if(this.dist(p)>r) return null; // 境界条件はテキトー
        double a = p.x - x;
        double b = p.y - y;
        double d = this.dist(p);
        double qx0 = x / 2 + p.x / 2 + Math.sqrt(r - d*d/4)*(-b)/d;
        double qy0 = y / 2 + p.y / 2 + Math.sqrt(r - d*d/4)*a/d;
        double qx1 = x / 2 + p.x / 2 - Math.sqrt(r - d*d/4)*(-b)/d;
        double qy1 = y / 2 + p.y / 2 - Math.sqrt(r - d*d/4)*a/d;
        Double2D[] ret = new Double2D[2];
        ret[0] = new Double2D(qx0,qy0);
        ret[1] = new Double2D(qx1,qy1);
        return ret;
    }
    // 三点の作る三角形の面積
    public static double squareOfTriangle(Double2D p1,Double2D p2,Double2D p3){
        double x1 = p1.x-p2.x;
        double x2 = p3.x-p2.x;
        double y1 = p1.y-p2.y;
        double y2 = p3.y-p2.y;
        double s = (x1*y2-x2*y1)/2;
        return Math.abs(s);
    }
    // 多角形の内部にあるかどうか 自己ループはないものとする
    public static boolean inPolygon(ArrayList<Double2D> poly){
        double acc = 0.0;
        int n = poly.size();
        for(int i=0;i<n-1;i++){
            acc += arg(poly.get(i),poly.get(i+1));
        }
        acc += arg(poly.get(n-1),poly.get(0));
        
        if(Math.abs(acc)<EPS) return false;
        else return true;
    }
    // -Pi < theta < Pi の範囲で偏角を返す　複素数でいう Arg(p2/p1)
    public static double arg(Double2D p1,Double2D p2){
        double ax = p1.x;
        double ay = p1.y;
        double bx = p2.x;
        double by = p2.y;
        double theta = Math.atan2((ax*by-ay*bx)/(ax*ax+ay*ay),(ax*bx+ay*by)/(ax*ax+ay*ay));
        return theta;
    }
    // 外接円の中心 Online Judge 検証なし 簡易検証あり
    public static Double2D centerOfCircumscribingCircle(Double2D p1, Double2D p2, Double2D p3){
        double x1 = p1.x, y1 = p1.y, x2 = p2.x, y2 = p2.y, x3 = p3.x, y3 = p3.y;
        double det = (x1*(y2-y3)+x2*(y3-y1)+x3*(y1-y2));
        if(Math.abs(det)<EPS) return null;
        double sx = 
            -y1*(x2*x2+y2*y2-x3*x3-y3*y3)
            -y2*(x3*x3+y3*y3-x1*x1-y1*y1)
            -y3*(x1*x1+y1*y1-x2*x2-y2*y2);
        double sy =
             x1*(x2*x2+y2*y2-x3*x3-y3*y3)
            +x2*(x3*x3+y3*y3-x1*x1-y1*y1)
            +x3*(x1*x1+y1*y1-x2*x2-y2*y2);
        return new Double2D(0.5*sx/det, 0.5*sy/det);
    }
}

class Line{
    double a;
    double b;
    double c;
    static final double EPS = 1.0E-9;
    Line(double a, double b, double c) {
        this.a = a;
        this.b = b;
        this.c = c;
    }
    Line(Double2D p1,Double2D p2){
        this.a = p1.y-p2.y;
        this.b = -p1.x+p2.x;
        this.c = -p1.x*(p1.y-p2.y)+p1.y*(p1.x-p2.x);
    }
    // 線分と線分の交点 ない場合は null
    // 2線分が重なりをもつときは考えていない
    public Double2D crossPoint(Line l){
        // 平行
        if(Math.abs(this.a*l.b-this.b*l.a)<EPS){
            return null;
        }
        else{
            double x = (-l.b*this.c+this.b*l.c)/(this.a*l.b-this.b*l.a);
            double y = (l.a*this.c-this.a*l.c)/(this.a*l.b-this.b*l.a);
            return new Double2D(x,y);            
        }
    }
}

class Circle{
    double x;
    double y;
    double r;
    static final double EPS = 1.0E-9;
    Circle(double x, double y, double r) {
        this.x = x;
        this.y = y;
        this.r = r;
    }
    public double centerDist(Circle c){
        return Math.sqrt((x-c.x)*(x-c.x)+(y-c.y)*(y-c.y));
    }
    // 2円の交点
    public Double2D[] crossPoint1(Circle c){
        double d = centerDist(c);
        if(r+c.r > d) return null; // 境界条件はテキトー
        double d1 = (d*d + r*r - c.r*c.r)/(2*d);
        double ex = x + (c.x-x)*d1/d;
        double ey = y + (c.y-y)*d1/d;
        double l = Math.sqrt(r*r-d1*d1);
        double lx = -(c.y-y)*l/d;
        double ly = (c.x-x)*l/d;
        Double2D[] ret = new Double2D[2];
        ret[0] = new Double2D(ex+lx, ey+ly);
        ret[1] = new Double2D(ex-lx, ey-ly);
        return ret;
    }
    // 1点を通る接線
    public Line[] tangentLine(Double2D p){
        Line[] ret = new Line[2];
        if(Math.abs((x-p.x)*(x-p.x)-r*r)<EPS){
            ret[0] = new Line(1, 0, -p.x);
            double m = ((y-p.y)*(y-p.y)-r*r)/(2*(x-p.x)*(y-p.y));
            ret[1] = new Line(m, -1, -m*p.x+p.y);
        }
        else{
            double m1 = ((x-p.x)*(y-p.y)+r*Math.sqrt((x-p.x)*(x-p.x)+(y-p.y)*(y-p.y)-r*r))/((x-p.x)*(x-p.x)-r*r);
            double m2 = ((x-p.x)*(y-p.y)-r*Math.sqrt((x-p.x)*(x-p.x)+(y-p.y)*(y-p.y)-r*r))/((x-p.x)*(x-p.x)-r*r);
            ret[0] = new Line(m1, -1, -m1*p.x+p.y);
            ret[1] = new Line(m2, -1, -m2*p.x+p.y);
        }
        return ret;
    }
}

class LineSegment{
    Double2D p1;
    Double2D p2;
    static final double EPS = 1.0E-9;
    LineSegment(double x1, double y1, double x2, double y2){
        p1 = new Double2D(x1, y1);
        p2 = new Double2D(x2, y2);
    }
    // 2線分が重なりをもつときは考えていない
    public Double2D crossPoint(LineSegment ls){
        Double2D p3 = ls.p1;
        Double2D p4 = ls.p2;
        if(fcompare(Math.max(p1.x, p2.x), Math.min(p3.x, p4.x))<0 ||
           fcompare(Math.max(p3.x, p4.x), Math.min(p1.x, p2.x))<0 ||
           fcompare(Math.max(p1.y, p2.y), Math.min(p3.y, p4.y))<0 ||
           fcompare(Math.max(p3.y, p4.y), Math.min(p1.y, p2.y))<0) return null;
        
        double dx00 = p2.x - p1.x;
        double dy00 = p2.y - p1.y;
        double dx01 = p3.x - p1.x;
        double dy01 = p3.y - p1.y;
        double dx02 = p4.x - p1.x;
        double dy02 = p4.y - p1.y;
        
        double dx10 = p4.x - p3.x;
        double dy10 = p4.y - p3.y;
        double dx11 = p1.x - p3.x;
        double dy11 = p1.y - p3.y;
        double dx12 = p2.x - p3.x;
        double dy12 = p2.y - p3.y;

        // 平行
        if(Math.abs(dx00*dy10-dy00*dx10)<EPS){
            if(p1.x==p3.x&&p1.y==p3.y) return p1;
            if(p1.x==p4.x&&p1.y==p4.y) return p1;
            if(p2.x==p3.x&&p2.y==p3.y) return p2;
            if(p2.x==p4.x&&p2.y==p4.y) return p2;
        }
        
        if(fcompare((dx00*dy01-dx01*dy00)*(dx00*dy02-dx02*dy00), 0)<=0 &&
           fcompare((dx10*dy11-dx11*dy10)*(dx10*dy12-dx12*dy10), 0)<=0){
            Line l1 = new Line(p1, p2);
            Line l2 = new Line(p3, p4);
            return l1.crossPoint(l2);
        }
        else return null;
    }
    private int fcompare(double a, double b){
        if(Math.abs(a-b)<EPS) return 0;
        if(a>b) return 1;
        else return -1;
    }

}

高速Kruskal with 高速MFSet

O(NlogN)のMFsetを用いたKruskalのアルゴリズム．

class KruskalSolver{
    int n;
    double[][] dist;
    double len;
    private PriorityQueue<Edge> q;
    KruskalSolver(int n,double dist[][]){
        this.n = n;
        this.dist = dist;
        this.q = new PriorityQueue<Edge>();
        this.len = 0.0;
        
        for(int i=0;i<n;i++){
            for(int j=i+1;j<n;j++){
                if(dist[i][j]<Double.MAX_VALUE){
                    q.offer(new Edge(i, j, dist[i][j]));
                }
            }
        }
        
        //Kruskal's algorithm
        MFSet mfs = new MFSet(n);
        while(!q.isEmpty()){
            Edge e = q.poll();
            int i = e.p1;
            int j = e.p2;
            if(mfs.find(i)!=mfs.find(j)){
                mfs.merge(mfs.find(i),mfs.find(j));
                len += e.length;
            }
        }
    }
}

class Edge implements Comparable<Edge>{
    int p1;
    int p2;
    double length;
    Edge(int i,int j,double l){
        p1 = i;
        p2 = j;
        length = l;
    }
    public int compareTo(Edge e){
        if(length>e.length) return 1;
        else if(length==e.length) return 0;
        else return -1;
    }
}

class MFSet{
    int n;
    private int[] setsize;
    private int[] firstelem;
    private int[] set;
    private int[] next;
    MFSet(int n){
        this.n = n;
        setsize = new int[n];
        Arrays.fill(setsize, 1);
        firstelem = new int[n];
        for(int i=0;i<n;i++) firstelem[i] = i;
        set = new int[n];
        for(int i=0;i<n;i++) set[i] = i;
        next = new int[n];
        for(int i=0;i<n;i++) next[i] = -1;
    }
    int find(int x){
        return set[x];
    }
    void merge(int a, int b){
        if(setsize[a]<setsize[b]){
            merge(b, a);
            return;
        }
        int i = firstelem[b];
        while(next[i]>=0){
            set[i] = a;
            i = next[i];
        }
        set[i] = a;
        next[i] = firstelem[a];
        firstelem[a] = firstelem[b];
        setsize[a] = setsize[a] + setsize[b];
        setsize[b] = 0;
        firstelem[b] = -1;
    }
}

Kruskal

O(N^2)のMFSetを用いたKruskalのアルゴリズム．

class KruskalSolver{
    int n;
    double[][] dist;
    double len;
    private PriorityQueue<Edge> q;
    int[] set;
    KruskalSolver(int n,double dist[][]){
        this.n = n;
        this.dist = dist;
        this.len = 0.0;
        this.q = new PriorityQueue<Edge>();
        this.set = new int[n];
        for(int i=0;i<n;i++) set[i] = i;
        
        for(int i=0;i<n;i++){
            for(int j=i+1;j<n;j++){
                if(dist[i][j]<Double.MAX_VALUE){
                    q.offer(new Edge(i, j, dist[i][j]));
                }
            }
        }
        
        //Kruskal's algorithm
        while(!q.isEmpty()){
            Edge e = q.poll();
            int i = e.p1;
            int j = e.p2;
            int si = set[i];
            int sj = set[j];
            if(si!=sj){
                len += e.length;
                for(int k=0;k<n;k++){
                    if(set[k]==sj) set[k] = si;
                }
            }
        }
    }
}

class Edge implements Comparable<Edge>{
    int p1;
    int p2;
    double length;
    Edge(int i,int j,double l){
        p1 = i;
        p2 = j;
        length = l;
    }
    public int compareTo(Edge e){
        if(length>e.length) return 1;
        else if(length==e.length) return 0;
        else return -1;
    }
}