SSブログ

ネムネコ、桁落ちに会う!! [ひとこと言わねば]

数日後の記事に使おうと、Fotranでパタパタとプログラムを作ってみた。
関数f(x)の微分f'(x)を使う箇所があるので、「微分なんていちいちしてられるか!!」と思い、この部分は中心差分を使って、微分係数f'(x)を求めることにした。
そうしたら、真面目に微分を計算したものと計算結果が違いやがった!!

  

中心差分というのは、f(x)の微分を上のように計算する方法。h=1.0×10⁻⁶くらいにとれば、真面目に微分を計算したものとほとんど一致するはずだ(註を参照)と考えたのだが――打切誤差という観点からすると、この予想はまったく正しい――、これが合わないんだケロ。h=1.0×10⁻⁴あたりで精度がもっとも良くなり、それ以上hを小さくとると、精度が逆に落ちてくるんだケロ。

「何故だろう?」としばし考えたところ、どうやら、桁落ち情報落ち、積み残し――以下、これらを総称して桁落ちと呼ぶことにする――が発生し、それで精度が落ちているみたいだね。

かれこれ1年ほど前に、数値計算のところで、丸め誤差、打切誤差、桁落ちなどを記事に書き、その時、桁落ちを発生させる計算例を示したけれど、あれは意図的に作り出した計算例で、ネムネコは、これまでに、自分の計算で、桁落ちで計算結果が大きく異なるという目に会ったことがなかった。それだけに、この計算結果を目の当たりにした際の衝撃は大きく、「何故(なぜ)に・・・」と、しばらく思考停止状態に陥ってしまうほどであった。
プログラムのどこかが間違っているというのならば、ソースコードを見なおせばすぐにわかるけれど、プログラム自体はどこも間違っていないんだから、そりゃ〜、すぐに、この不可解な現象の原因特定ができないケロ。

真面目に微分を計算した計算結果


微分を中心差分で計算した結果


計算結果が4.11111ホニャララになるところが約4.068になっており、計算結果が1%以上も異なってしまう。微分はもっとひどくて、なんと、最大約10%の誤差、アンダー・エスティメイトになるんだよね〜。しかも皮肉なもので、微分を真面目に計算した奴は20回の反復で収束しないのに対して、微分を中心差分で代用したものは反復回数18回で収束する――微分係数を過小評価するので、その結果、修正量が大きくなり、計算が早く収束する!!――。微分計算で桁落ちさせ、この部分の有効精度を落とした方が収束が早いのだから、もう笑うしかない。
果たして、どちらが、数値解法として優れているのやら(苦笑い)


う・ふ・ふ・ふ 毎日だれかに
う・ふ・ふ・ふ 見られることが
ビタミンになる

ではあるが、ネムネコの場合、これだけでは足りないにゃ。


これくらいしないと・・・。

さて、この事象の根本的解決策は、単精度の計算を倍精度にする以外にないね。計算が単純すぎるので、計算方法の工夫によって桁落ちを解消するということは、どうやら、できそうにない。
プログラムを作るのが楽なので、Fortranでプログラムを作ったのだけれど、使用する変数の全ての型宣言をするのならば、C言語の方が楽だな。しかも、精度的にも、FotranよりCの方が優れているし。
30行に満たないオモチャみたいなプログラムだから、最初から作りなおしたところで大変じゃないけれど――プログラムを書き、コンパイル、バグ取り、計算を実行させてたとしても、この全作業は5分以内に終わるはず。10分なんて時間は、絶対に、かからない――、プログラムを作った時間が無駄になってしまったケロよ。

でも、実際、「桁落ちで本当に精度が落ちてしまう」ということを知ることができて良い経験になった。これまで、そんなことに遭遇したことがなかったから。



微分を中心差分で計算した結果(C言語で倍精度計算)



表計算ソフトで計算した結果(真面目に微分計算)



プログラムで解いたものは、3列目が修正量、対して、表計算ソフトの3列目は推測値と解との誤差で違うのだけれど、倍精度で計算すれば、桁落ちが解消されていることがわかる。

ここしばらくC言語なんて使っていなかったから、「C言語の出力コマンドprintfの書式ってどうだったっけ」と、一瞬、考えこんでしまったケロよ。

C言語は、UNIXというOSを開発するために作られた、限りなくアセンブラに近いプログラミング言語だから、入出力系のコマンドは弱いね〜。

【補足】

情報落ち、積み残しの例 通常のコンピュータの計算(浮動小数点を用いた計算)では、百億+1=百億となるので、1を加えたという情報が落ちてしまう、消失する。このように、絶対値の大きな数と絶対値の数を加えた時、絶対値の小さな数が計算結果に反映されない現象を情報落ちや積み残しなどという。
金銭計算、大きな整数同士の四則演算に特化したプログラム言語、COBOLが現在も使われ続けるのは、整数同士の演算に限っては、このような情報落ち、積み残しが起きにくいため。COBOLならば、1000兆+1=1000兆1と、正しい値を出してくれる。お金が絡む計算だから、1円でも計算結果があわないと、大変なことになるからね〜。

桁落ちの例 たとえば、近い数字同士で引き算をすると、計算の有効精度が大きく落ちる、このような現象を桁落ちという。100.0001−100.0000=0.0001 7桁の精度をもっていたものが、この引き算を行うと、精度はなんと1桁に落ちてしまう!!
 ――説明を容易にするために、ここでは、7桁、1桁の精度と書いた。――
それだけで済めばいいが・・・


単精度のFortranの場合、さらに、丸め誤差の影響がこれに加わり、100.0001−100.0000という引き算の結果は上のようになってしまうのであった。


【註】

  

h=0.000001として、ネムネコのPCにある電卓で計算すると、

  

となるんですがね〜。打切誤差が高精度な電卓の計算の丸め誤差の範囲に収まり、厳密な36とドンピシャ一致する。

しかし、Fortranなどのプログラミング言語でこの計算をすると、この値は36とかなり違う値が出るんだケロ。「誤差はあっても、せいぜい、0.1%程度だろう」という、ネムネコの甘い期待は木っ端微塵に打ち砕かれてしまった。

それにしても、電卓って、ホント、高精度だよな〜。

 


nice!(2)  コメント(0) 
共通テーマ:音楽

今日のクラシック、グローフェ作曲『グランド・キャニオン』 [今日のクラシック]

今日のクラシックは、グローフェ作曲『グランド・キャニオン』です。


日本のコアなクラシック・ファンの間では、ホルストの『惑星』とともに、この曲は、セミ・クラシックの曲とされ、当たり前のクラシックとして扱ってもらえない曲であります。ですから、通を気取るクラシック・ファンの前では、口が避けても、『グランド・キャニオン』が好きだとか、『グランド・キャニオン』を名曲だと言ってはいけません。言ったら最後、「コイツは本格的なクラシック音楽を理解できない」という、生涯、消すことのできない烙印を押されてしまいますから(笑)。

写実的、絵画的な音楽で、いわゆる、絶対音楽とされる曲から外れる曲ですが、この曲は間違いなく名曲だと思いますがね〜。


また、クラシック・マニアの間には、
 ドイツ・オーストリア音楽>フランス音楽>ロシア音楽>・・・>イギリス音楽>アメリカ音楽
という、暗黙の、序列、格付けがあります。
ホルストはイギリス、グローフェはアメリカの作曲家ですから、この観点からしても、『グランド・キャニオン』はセミ・クラシックに分類されてしまうのです(^^ゞ

『グランド・キャニオン』と比較すると知名度はかなり落ちますが、グローフェには、『ミシシッピ組曲』という曲もあります。


そして、グローフェは、ガーシュインの『ラプソディー・イン・ブルー』の編曲者として、名を残す作曲家でもあります。



nice!(0)  コメント(0) 
共通テーマ:音楽

今日のアニソン、「ハイスクール・フリート」から『High Free Spirits』 [今日のアニソン]

今日のアニソンは、アニメ「ハイスクール・フリート」から『High Free Spirits』です。


たぶん、この曲はこれまで取り上げたことがないと思うのですが・・・。
さらに、このアニメのED曲を♪



nice!(0)  コメント(4) 
共通テーマ:音楽

第8回 単射、全射、逆写像 [集合論入門]

第8回 単射、全射、逆写像

 

§1 単射と全射

 

tanzen-001.pngを写像とする。

任意のtanzen-002.pngに対して

  tanzen-003.png

であるとき、fは単射であるという。

単射の定義には、(1)の対偶をとった次のものを使ってもよい。

任意のtanzen-002.pngに対して

  tanzen-004.png

 

問 次のことを示せ。

(1) f(x)=x³で定義される写像f:R→Rは単射である。

(2) f(x)=x²で定義されるf:R→Rは単射でない。

【解】

(1) f(x₁)=f(x₂)とすると、

  

これが成立するのは、

  tanzen-006.png

よって、

  tanzen-007.png

ゆえに、この写像は単射である。

 

(2) f(1)=f(−1)=1だから、この写像は単射でない。

(解答終)

 

tanzen-001.pngを写像とする。

任意のy∈Yに対して、f(x)=yを満たすx∈Xが存在するとき、すなわち、

  tanzen-008.png

であるとき、f全射という。

また、fが全射かつ単射であるとき、f全単射であるという。

 

問2 次のことを示せ。

(1) f(x)=2x+3で与えられる写像f:R→Rは全射である。

(2) f(x)=x²で与えられる写像f:R→Rは全射でない。

【解】

(1) 任意のy∈Rに対して

  tanzen-009.png

よって、この写像は全射である。

 

(2) y=−1とすると、

  tanzen-010.png

を満たす実数xは存在しない。

したがって、この写像は全射ではない。

(解答終)

 

(1)のf(x)=2x+3で与えられるf:R→Rは単射でもあるので、全単射である。

 

問3 f(x)=x²で与えられる写像f:X→Yがあるとする。

(1) Xを実数全体の集合Rとし、Y=y∈Ry≧0}とすると、fは全射になることを示せ。

(2) X=x∈Rx≧0}、Y=y∈Ry≧0}とすると、fは全単射になることを示せ。

【解】

(1) 任意のy≧0に対して、とおくと、x∈Xであり、

  

となるので、fは全射である。

 

(2)は略

(解答終)

 

定理1 写像について次のことが成り立つ。

(1) を満たす写像が存在すれば、fは全射である。

(2) を満たす写像が存在すれば、fは単射である。

(3) かつを満たす写像ghが存在するならばfは全単射である。また、g=hである。

【証明】

(1) 任意のy∈Yに対して、とすれば、

  tanzen-011.png

よって、fは全射である。

 

(2) だから

  

とすると、

  tanzen-012.png

よって、

  tanzen-016.png

ゆえに、fは単射である。

 

(3) (1)、(2)よりfが全単射。

fは全単射だからy=f(x)を満たすx∈Xがただ一つ存在する。

一方、任意のyに対して

  

が成り立つので、

  

(証明終)

 

定理2 XYZを集合、tanzen-013.pngを写像とする。このとき、次のことが成り立つ。

(1) tanzen-014.pngが単射ならば、fは単射である。

(2) tanzen-014.pngが全射ならば、gは全射である。

【証明】

(1) f(x₁)=f(x₂)とすると、

  tanzen-015.png

tanzen-014.pngは単射なので、x₁=x₂である。

したがって、

  tanzen-016.png

よって、fは単射である。

 

(2) tanzen-014.pngは全射なので、任意のz∈Zに対して、あるx∈Xがあって、

  

そこで、y=f(x)∈Yとおくと、

  

よって、gは全射である。

(証明終)

 

問3 写像tanzen-013.pngであるfgともに全単射ならば、XからZへの合成写像tanzen-014.pngは全単射であることを示せ。

【解】

仮定より、任意のz∈Zに対してz=g(y)であるy∈Yが存在し、また、y=f(x)であるx∈Xも存在し、

  tanzen-018.png

よって、tanzen-014.pngは全射である。

また、とすると、fは単射なので、

  

また、gも単射なので、

  

よって、tanzen-014.pngは単射である。

したがって、tanzen-014.pngは全単射である。

(解答終)

 

この逆、

tanzen-014.pngが全単射ならば、fgは全単射である

は、一般に、成立しないので、注意。

 

§2 逆写像

 

写像が全単射ならば、任意のy∈Yに対して、あるx∈Xが存在してy=f(x)である。一方、このようなx∈Xは各yに対してただ一つである。したがって、任意のy∈Yに対して、y=f(x)を満たすx∈Xがただ1つ定まり、写像を定義することができる。この写像gも全単射であり、を満たす。このような写像f逆写像といい、f⁻¹で表す。

 

問4 写像が全単射のならば、次のことが成り立つことを示せ。

  tanzen-019.png

【解】

(1) 任意のx∈Xに対して

  tanzen-020.png

(2) 任意のy∈Yに対して

  tanzen-021.png

 

定理 とする。このとき、次のことが成り立つ。

  tanzen-019.png

【証明】

(1) 任意のx∈Xに対して

  

とおけば、

  tanzen-025.png

よって、

  tanzen-023.png

 

(2)、さらに、とすると、

  

である。

(証明終)

 


nice!(0)  コメント(0) 

この広告は前回の更新から一定期間経過したブログに表示されています。更新すると自動で解除されます。