読者です 読者をやめる 読者になる 読者になる

領域の内側か外側かを調べる

はじめに

複素関数の定理で閉領域内なら2πi、領域外なら0を返すものがあります。この定理を使うことで指定した点が領域内にあるかどうかを判定できます。参考書は「複素関数攻略の一本道,森北出版」です。

2πiの定理

コーシーの積分定理はwikipediaによると、以下の通りです。

D を単連結領域とし、f(z) は D 上で正則である複素関数とするとき、C を D 内にある長さを持つ閉曲線とすると、
f:id:wildpie:20141019145106p:plain

この定理から「2πiの定理」が成り立つことが分かります。
f:id:wildpie:20141019145758p:plain
z0がC外なら関数はC内で正則になるので、積分値はコーシーの積分定理により0となります。
一方、z0がC内なら
f:id:wildpie:20141019150942p:plainとおくと、
f:id:wildpie:20141019151010p:plain
より、積分値は2πiになります。
つまり「2πiの定理」を評価すれば点が閉領域内なのか領域外なのかを調べることができます。

実装

実装をするためには積分の計算をする必要があります。積分の実装は難しいので、参考書では"G.S. Gipson, Use of the residue theorem in locating points within an arbitrary multiply-connected region, 1986."で提案されている式変換を使っています。
これは閉曲線CをN個のセグメントに分けて計算する手法です。
f:id:wildpie:20141019152537p:plain
より、
f:id:wildpie:20141019152617p:plain
が成り立ちます。
f:id:wildpie:20141019153914p:plain
θはz-z0が実軸となす角です。r1をz0を始点として、セグメントの開始位置までのベクトルとすると、
f:id:wildpie:20141019153302p:plainより
f:id:wildpie:20141019153307p:plain
となります。また、閉曲線の進む向きは外積で求めることができます。

実行例

適当に決めた点の内側にカーソルがあるかを判定します。マウスのMoveイベントで毎回評価しています。
f:id:wildpie:20141019154902p:plain
ソースコード

private enum Result
{
    IN,
    OUT,
};

List<PointF> _points = new List<PointF>
        {
            new PointF(100f, 49.0f),
            new PointF(59f, 77.0f),
            new PointF(57.0f, 178.0f),
            new PointF(122.0f, 207.0f),
            new PointF(191.0f, 167.0f),
            new PointF(197.0f, 105.0f),
            new PointF(135.0f, 77.0f),
        };

public Form1()
{
    InitializeComponent();

    var bitmap = new Bitmap(pictureBox1.Size.Width, pictureBox1.Size.Height);
    pictureBox1.Image = bitmap;

    // 点と線を表示
    var g = Graphics.FromImage(bitmap);
    g.DrawPolygon(Pens.Blue, _points.Select(v => new Point((int)v.X, (int)v.Y)).ToArray());
    foreach (var p in _points)
    {
        //g.FillRectangle(Brushes.Black, new Rectangle((int)p.X, (int)p.Y, 10, 10));
        g.FillEllipse(Brushes.Black, new Rectangle((int)p.X-5, (int)p.Y-5, 10, 10));
    }
}

private Result InOrOut(List<PointF> points, PointF z0)
{
    // 内積
    Func<PointF, PointF, double> dp 
        = (r1, r2) => r1.X * r2.X + r1.Y * r2.Y;

    // ノルム
    Func<PointF, double> norm 
        = r => Math.Sqrt(r.X * r.X + r.Y * r.Y);

    // 一個ずらした要素をまとめる
    var pack = points.Zip(points.Skip(1), (z1, z2) => new {Z1=z1, Z2=z2}).ToList();
    pack.Add(new { Z1=points.Last(), Z2=points.First() });

    double integral = 0.0;
    foreach (var p in pack)
    {
        var r1 = new PointF(p.Z1.X - z0.X, p.Z1.Y - z0.Y);
        var r2 = new PointF(p.Z2.X - z0.X, p.Z2.Y - z0.Y);

        // 反時計回り+ or 時計回り-
        int sign = Math.Sign(r1.X*r2.Y - r1.Y*r2.X);

        integral += sign * Math.Acos(dp(r1, r2) / (norm(r1) * norm(r2)));
    };

    // 積分結果がだいたい2πiなら内側
    const double epsilon = 0.1;
    return Math.Abs(Math.Abs(integral / (2.0 * Math.PI)) - 1.0) < epsilon ? Result.IN : Result.OUT;
}

private void pictureBox1_MouseMove(object sender, MouseEventArgs e)
{
    var r = InOrOut(_points, new PointF(e.X, e.Y));
    this.label1.Text = string.Format("({0}, {1}) : {2}", e.X, e.Y, r);
}

複素関数の定理がプログラムに使えるのはおもしろいですね。関係ないですが「複素関数攻略の一本道」は途中式を省略せずに書いてあるので、数学の苦手は私でも分かった気になります。そういう本がふえてくれないかなー。