Hatena::ブログ(Diary)

純粋関数型雑記帳 このページをアンテナに追加 RSSフィード Twitter

このページはHaskellを愛でるページです。
日刊形式でHaskellなどについての記事をだらだらと掲載しとります。
2004 | 07 | 08 | 09 | 10 | 11 | 12 |
2005 | 01 | 03 | 04 | 05 | 06 | 07 | 08 | 09 | 10 | 11 | 12 |
2006 | 01 | 02 | 03 | 04 | 06 | 07 | 08 | 09 | 11 |
2007 | 03 | 04 | 05 | 07 | 08 | 09 | 12 |
2008 | 02 | 04 | 05 | 06 | 07 | 08 | 09 | 10 | 11 | 12 |
2009 | 03 | 05 | 06 | 09 | 10 | 11 | 12 |
2010 | 02 | 03 | 04 | 05 | 06 | 07 | 08 | 09 | 10 | 12 |
2011 | 01 | 02 | 05 |
本体サイト

2005年03月08日(火) 複素数と平面幾何

この前の文章がやねうらおさんのところで紹介されたもんで、

予想以上の反響になった。

何はともあれ、いろいろな人に読んでいただけて有難いもんです。

内容のほうが今ひとつ煮え切らないような感じになってしまっていたので

もっと批判的なご意見もいただくかと思ったが、

同じように感じているというものが多く、

嬉しくもあり意外でもあった。

まぁ、なにかもっと真髄のようなものをつかめたら

続きでもかきます。

[]複素数による平面幾何プログラミング 複素数による平面幾何プログラミングを含むブックマーク 複素数による平面幾何プログラミングのブックマークコメント

ICPCが近いので、ライブラリを作らないといかん…。

というわけで、平面幾何。

平面幾何を扱うのにcomplexを使うのはもはや常識(?)

なのかもしれないが、そこをちょっとおさらいしてみることにする。

最初このテクニックを見たとき、

自分でstruct point{ int x,y; }; とか、

演算子の定義とかを書く必要がないのか、と思ったのだが

私の考えが浅かった。複素数で点を扱うと、

複素数平面上で平面幾何を扱うことになり、

便利な関数がたくさん利用できるということになる。

さらに、回転が簡単にかけたりするので、問題を回転させて

より簡単な問題に帰着したりできるかもしれない。


typedef complex<double> pt;

pt x,y;
pt z=x+y;
double c;
pt w=c*z;

このようなことが何もせずに使える。

  • 二点間の距離
abs(x-y);

ベクトルvのx軸となす角は、

arg(v);

ベクトルuとvのなす角は

arg(u/v);

と書ける。


ベクトルuをベクトルvに移す写像fを

pt f=polar(abs(u)/abs(v),arg(u/v));

簡潔に書ける。

ベクトルu,vの内積

real(u*conj(v));

で求めることができる。

(成分を計算することにより簡単に確かめられる)

ベクトルu,vに対して行列[u,v]の行列式を

imag(conj(u)*v);

で求めることができる。

(同上)

  • 三角形の符号付面積
double area(pt a,pt b,pt c)
{
  return imag(conj(b-a)*(c-a))*0.5;
}

点aを原点に平行移動させて、行列式より計算できる。

  • 二直線の交点

直線 a+su,b+tv (a,b,u,v はベクトル、s,tはスカラ変数)

の交点。(s,tを返す)

pair<double,double> cross_lines(pt a,pt u,pt b,pt v)
{
  double det=imag(conj(u)*v);
  if (abs(det)<epsilon) throw "交差しない";
  return make_pair(
    -imag(conj(v)*(b-a))/det,
    -imag(conj(u)*(b-a))/det);
}

証明は…成分を計算すると確かに正しいのだが、

式の意味するところは良く分からんです…。

交差しないときに例外を投げてしまっているが、

これはあんまり良くないので、実際に使うときは

うまいインターフェースを考えよう。

  • 円と直線の交点

直線 x=a+su と円 |x-c|=r の交点を求める。

(sを二つ返す)

pair<double,double> cross_circle_and_line(pt a,pt u,pt c,double r)
{
  pt ac=a-c;
  double b=real(u*conj(ac));
  double d=b*b-norm(u)*(norm(ac)-r*r);
  if (d<0) throw "交差しない";
  return make_pair(
    (-b+sqrt(d))/norm(u),
    (-b-sqrt(d))/norm(u));
}

単純に直線と円の式を解くと上の式は得られる。

  • 二円の交点

円 |x-c|=r と |x-d|=s の交点を求める。

pair<pt,pt> cross_circles(pt c,double r,pt d,double s)
{
  double di=abs(c-d);
  if (di>r+s || di<abs(r-s))
    throw "交差しない";
  double l=((r*r-s*s)/di+di)/2.0;
  double n=sqrt(r*r-l*l);
  pt e=c+(d-c)*l/di;
  pt p=(d-c)*n/di*pt(0,-1);
  return make_pair(e+p,e-p);
}

二つの円の交点と、二つの円の中心のなす凧形を考えれば

式が導出できるが、ややこしいので省略。

複素数ならではなことはできていないが、

短くする効果は十分だろう。

yaneuraoyaneurao 2005/03/09 09:28 こうして見ると複素数での平面幾何って面白いですネ。でも、3次元上の幾何には使えないので、どっちかって言うと、行列演算をtemplateで実装したものがあったほうがいいような..?

Nazo9xNazo9x 2005/03/09 17:42 使っていないので詳細はわからないけど。こんなのあるよ。
http://www.page.sannet.ne.jp/d_takahashi/boost/ublas/index.html

tanakhtanakh 2005/03/10 03:36 三次元以上の幾何は正直問題自体が難しいです…。
行列周りのライブラリがあれば確かに便利でしょうけど。
boostにこんなのがあったんですか。
boostは標準じゃないので、コンテストで使えないのが痛い…
逆に使えるとC++があまりにも有利になってしまいますが。

daismdaism 2005/03/11 20:07 検討違いならごめんなさい
http://vision.kuee.kyoto-u.ac.jp/~hiroaki/numerical/

tanakhtanakh 2005/03/11 23:54 これも行列が扱えるんですかね。いろいろあるんですね。

トラックバック - http://d.hatena.ne.jp/tanakh/20050308
Connection: close