半径rの円があって、図のグレー部分の面積Sが分かっています。

半径rの円があって、図のグレー部分の面積Sが分かっています。
このとき、xの長さを求めたいのですが、どうすればいいでしょうか？

ちなみにプログラムで計算するので、三角関数と逆三角関数が計算可能です。

No.1 ベストアンサー 回答者： himajin100000 回答日時： 2010/06/29 03:31

(図を後で添付)

弧ABと弦ABに囲まれた部分の面積がS
原点をO、半径rとする円がある。
∠AOD = p(単位:ラジアン)とする。
なお、CD = xとおき、OD = r - x, Math.Cos(p) = (r - x) / r…A
であるとする。

このとき、弧ACと弦ADに囲まれた部分の面積は S / 2…式1
また、この部分の面積は
扇形ACOの面積 Math.Pow(r,2) * Math.PI * p / (2 * Math.PI) = Math.Pow(r,2) * p / 2 から
三角形AODの面積
1/2 * OD * DA
= 1 /2 * r * Math.Cos(p) * r * Math.Sin(p)
= 1 /2 * Math.Pow(r,2) * Math.Cos(p) * Math.Sin(p)
を引いたものといえ、
Math.Pow(r,2) * p / 2 - 1 /2 * Math.Pow(r,2) * Math.Cos(p) * Math.Sin(p)…式2
前者の算出法と後者の算出法の誤差を求めることを考え（式2 - 式1)、この関数をf(p)とする。
f(p) = Math.Pow(r,2) * p / 2 - 1 /2 * Math.Pow(r,2) * Math.Cos(p) * Math.Sin(p) - S / 2

f(p)が0になるようなpを求めるためにニュートン法を使う。
f(p)のpによる微分fdiff(p)は
Math.Pow(r,2) / 2 - 1 /2 * Math.Pow(r,2) * (Math.Cos(p) * Math.Cos(p) - Math.Sin(p) * Math.Sin(p)) = Math.Pow(r,2) / 2 - 1 /2 * Math.Pow(r,2) * Math.Cos(2 * p)

でニュートン法をそのまま適用。pを0から始めると、

pさえ解ってしまえば
Math.Cos(p) = r - x / r
しかないわけで後は簡単。

という感じでやってみた。ニュートン法のわかりやすい説明って中々なく、説明省略orz

#上記の式ではつい整数同士の除算とか放置しているけど、そのへんを考えて書いてみた物。
値はp = Math.PI / 3を想定してr,とSを与えた例。大体あってますね。
pを0から始めるとfdiffが0とかになるのか、NaNが答になっちゃうのでMath.PI / 2からスタートです。

http://ideone.com/okI58

この回答へのお礼

コード付きの丁寧な回答ありがとうございます。
xからSを出すプログラムで検算してみたところ、精度ばっちりでした！

お礼日時：2010/06/29 21:57