ということで、Bézier曲線で囲まれた領域の面積とか(2)

ということで、Bézier曲線で囲まれた領域の面積とか - らんだむな記憶で面積計算のための式を導いた。が検証はしていなかった。
まぁ、計算とか苦手なんで間違っている可能性は大なのだが！

すでに準備はできている。rubyとBézierと(2) - らんだむな記憶で十分にBézier曲線で円を近似できることを感覚的に見たので、半径1の円を近似して、その面積を求めると、円周率に十分近い値が求まるはずだ。真円よりはちょっと膨らむということなのだから、円周率に近いけど、ちょっと大きい値、になるはずだ。
さて、rubyで書こう。お料理番組ではそう言った次の瞬間にカメラが横に移動して完成したボールが出てくる。勿論、ここでも同様に書いた後のコードを貼り付ける。しかも、ちょっと富豪的プログラミング(念のため引用:富豪的プログラミング)してみたw 今回はUIとか全然関係ないけどね！インスタンスなんかどんどん生成すればいいんだよ！GCが回収するから！(とかすると結構コストが嵩んで遅くなるんだけどね。破壊的メソッドを使ってインスタンスを使いまわしたほうが速いと思う。)

#! /usr/bin/ruby -Ku

K = 4.0 * (Math.sqrt(2) - 1) / 3

class Point
    attr_reader :x, :y

    def initialize(x, y)
        @x, @y = x, y
    end

    def coord
        return @x, @y
    end

    def +(rhs)
        return Point.new(@x + rhs.x, @y + rhs.y)
    end

    def *(scalar)
        return Point.new(@x * scalar, @y * scalar)
    end

    def cross_product(rhs)
        @x * rhs.y - rhs.x * @y
    end
end

class Bezier
    def initialize(p0, p1, p2, p3)
        @p0, @p1, @p2, @p3 = p0, p1, p2, p3
    end

    def pos(t)
        return (@p0 * (1-t)**3 + @p1 * 3 * (1-t)**2 * t + @p2 * 3 * (1-t) * t**2 + @p3 * t**3).coord
    end

    # 面積計算のためのある量の20倍の値を返す。直接的に面積を与えるものではない。
    def twenty_times_area_value_of_segment
        6 * @p0.cross_product(@p1) + 3 * @p0.cross_product(@p2) + @p0.cross_product(@p3) + 3 * @p1.cross_product(@p2) +  3 * @p1.cross_product(@p3) + 6 * @p2.cross_product(@p3)
    end
end

class BezPath
    def initialize
        @bez_segments = []
    end

    def add(bez)
        @bez_segments << bez

        self
    end

    def area
        twenty_times_area = 0
        @bez_segments.each do |seg|
            twenty_times_area += seg.twenty_times_area_value_of_segment
        end

        Float(twenty_times_area) / 20
    end
end


################################################################################

b1 = Bezier.new(Point.new(1,  0), Point.new(1,   K), Point.new(K,   1), Point.new(0,  1))
b2 = Bezier.new(Point.new(0,  1), Point.new(-K,  1), Point.new(-1,  K), Point.new(-1, 0))
b3 = Bezier.new(Point.new(-1, 0), Point.new(-1, -K), Point.new(-K, -1), Point.new(0, -1))
b4 = Bezier.new(Point.new(0, -1), Point.new(K,  -1), Point.new(1,  -K), Point.new(1,  0))

print "#{BezPath.new.add(b1).add(b2).add(b3).add(b4).area}\n"

Bézierの座標計算で Point×数値の順にしないとならんのは残念。数値×Pointも実現するには、数値型のクラスに一時的にPoint型変数をも引数として受け取る形での四則演算子のオーバーロードをかまさんとならないが、それをするには特異メソッドや特異クラスのような "黒魔術" すれすれな感じになって嫌だからやらない。やりすぎるとドラッグと同じで廃人になりそうだw

Rubyの黒魔術でも読んで

"ムスカ".balse!

したらいいと思うw

なお、面積計算の結果は、

$ ruby bezier_test.rb
3.1424723326565074

だった。円周率は3もとい、3.1415以下略であることを考えると、一応期待する結果かなということで、たぶん元の面積計算式も恐らくあっていたのだろう。
なお、$\pi$の連分数第一近似は$\frac{22}{7} = 3.142857...$なのだが、案外こっちに近い値が出たなぁ、とかかんとか。いや、ただの偶然だけど。