So-net無料ブログ作成

式が消えた!! また、余計な真似をしてくれやがって [ひとこと言わねば]

何か、So-netブログ、変なことをしてくれたようだにゃ。

編集画面から過去に書いた数学の原稿を記事を見ると、式が表示されない!!
なぜかわからないけれど、FireFoxからだと、原稿から式が消えて見ることをできない!!



これじゃ〜、保存した原稿の修正などの編集ができないケロよ。
毎日、10万字以内エラーで苦しめられているというのに、どこまで、ネムネコを傷めつけるつもりだにゃ、このブログは。


幸い、Chomeもどきのブラウザだと、式を見ることができる。



So-netブログは、オレに、
「2つのブラウザを同時に使って、ブログの数学の記事を書き、そして、編集しろ」
って言いたいケロか。
アンマリだにゃ(涙)

ネムネコがこのあまりにひどい仕打ちにブチッと切れて、怒りに任せて、「こんなブログは、永遠に、この世から消し去ってやる」と、「ブログを削除する」をクリックしたら、どうするつもりだにゃ。



上の赤いところをクリックするだけで、「ねこ騙し数学」というブログをこの世から消し去ることができる!!んだケロよ。

質・量ともに、「ねこ騙し数学」に取って代わる数学のブログなんてどこにもないケロよ。だから、日本の将来にとって、とんでもなく大きな禍根を残すことになると思うにゃ。

ただでも、So-netブログは使いづらくて、役立たずのブログなんだから、余計な真似をして、これ以上、ネムネコの足を引っ張る真似はやめて欲しいケロ。



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

スプレッドシートを公開(4月16日)!! [ひとこと言わねば]

例によりまして、今日4月16日の数学の記事に使用したスプレッドシートを公開したにゃ。



「こんな簡単なものを公開して、いったい、何の役に立つのか」という、強い思いはあったけれど、清水の舞台から飛び降りる思いで、公開したにゃ。



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

今日のアニソン、東方から『燐火〜根ノ國〜』 [今日のアニソン]

今日のアニソンは、東方から『燐火〜根ノ國〜』です。


我ら(化け)ネコ界のスーパーアイドルである「お燐」の曲だにゃ。だから、紹介しないわけにはいかないケロ!!


上の2つの曲は、ネムネコがよく使っている




とは違って、『廃獄ララバイ』という曲だにゃ。



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

境界値問題の前進積分による解法1 [数値解析]

境界値問題の前進積分による解法

 

微分方程式の境界値問題を数値的に解くには、差分法などを使い離散化した連立方程式を解く必要があった。しかし、2階の境界値問題は、初期値問題のように、前進積分法によって解ける場合がある。その方法を紹介する。

 

問題 非線形境界値問題

  

h=0.2として、前進積分法で解け。

【解】

  

とすると、

  

h=0.2だから、

  

という漸化式が得られる。

問題の条件よりy(0)=0だからy₀=0

ここで、y₁=0.3と仮定すると、(2)式よりy₂

  

以下同様に、

  

が得られる。

1回目のy₁の推測値を後で使うので、と書くことにする。

y(1)=1だから、

  

y₁=0.3だと、y₅y(1)=1を超過するので、y₁0.3より少し減らし、y₁=0.25とし再計算する。

  

したがって、

  

ここで、補間を使って、y₁を次のように推測する。

  

この値を用いて、また、再計算する

  

(解答終)

 

以下に、表計算ソフトを用いて計算した結果を示す。

 

 

 

あくまで、このように計算できる場合もあるという話!!

 

ということで、これまでのように、連立方程式を解く方法でこの問題を解くことにする。

 

微分方程式を差分で近似すると、

  

いつも同じ方法だと芸がないということで、今回は、この連立2次方程式をニュートン法で解いた。

やっぱ、ニュートン法は、収束が早いね。

 

以下に、そのプログラムを示す。

 

parameter(n=10)
real y(0:n)
real a(n),b(n),c(n),d(n)
real w(n)

eps=1.e-6

! 初期化
y=0

! 境界条件
y(0)=0.; y(n)=1.0

h=1./n

! ニュートン法
do k=1,10
    do i=1,n-1
        a(i)=1.
        b(i)=2.*(h*h*y(i)-1)
        c(i)=1.
        d(i)=y(i-1)+(h*h*y(i)-2.)*y(i)+y(i+1)
    end do
   
    call tdma(a,b,c,d,n-1)

    e_max = 0.
    do i=1,n-1
        y(i)=y(i)-d(i)
        e_max=amax1(e_max,abs(d(i)))
    end do
!    write(*,*) k, e_max
    if (e_max.lt.eps) exit
end do

! 結果の出力
write(*,*) ' k      x       y'
do i=0,n
    write(*,100) i, i*h, y(i)
end do

100 format(i3,1x,f8.5,1x,f8.5)

end

!  TDMA
subroutine tdma(a,b,c,d,n)
real a(n), b(n), c(n),d(n)
do i=1,n-1
    ratio=a(i+1)/b(i)
    b(i+1)=b(i+1)-ratio*c(i)
    d(i+1)=d(i+1)-ratio*d(i)
end do
d(n)=d(n)/b(n)
do i=n-1,1,-1
    d(i)=(d(i)-c(i)*d(i+1))/b(i)
end do

end

 

計算結果は、次の通り。

 

Hisen-Result.png

 

ほとんど直線的に変化しているということもあるのだけれど、問題の解答のような粗い計算でも、小数点3位くらいまではこの数値解と一致しており、結構、いい感じで解けていることがわかる。

 

問題の微分方程式

  

は、一見、解くのは簡単そうに見えるけれど・・・。

解けると思うヒトは、是非、チャレンジして欲しい。

そして、もし、運良く解くことができたら――正確には、三角関数や指数、対数関数などの初等関数を用いて、この解を表わせたら――、教えて欲しいにゃ。

 

 


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

オレって実はスゴイ(笑) Part2 [ひとこと言わねば]

ネムネコが、昨夜、閃いた(思いついた)、4次のルンゲ・クッタ法を用いた2階常微分方程式の境界値問題の解法は、第一種の境界条件と第2種の境界条件の両方が使用されている場合にも、適用が可能なようだね。

  

h=0.1で、上の境界値問題を解いた場合の計算結果は次の通り。

 

Nemuneko-ho.png

 

参考までに、h=0.01とし、差分を用いた結果も以下に示す。

 

 

h=0.01とルンゲ・クッタ法の10倍計算格子の間隔を細かくとっているものより、心もち、高精度で計算できていることがわかる。

 

 

ではあるが、

混合条件にまで拡張できるか、どうやって拡張したらよいかについてはまったくの白紙状態。

大体、ネムネコは、数値計算に関しては素人だにゃ。

化けネコではあるが、ネコがここまでできるってこと自体が既に驚異だにゃ。

そう思わないかい?




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

今日のアニソン、「名探偵ホームズ」から『グッバイ・スウィートハート』 [今日のアニソン]

今日のアニソンは、映画「名探偵ホームズ」から『グッバイ・スウィートハート』です。


ネムネコは、この曲↑、すごく好きなんだケロ。


ED曲、いいよね。そう思わないケロか?


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

差分法による2階常微分方程式の境界値問題の数値解法3 [数値解析]

差分法による2階常微分方程式の境界値問題の数値解法3

 

前回、次の微分方程式

  

x=1における第2種の境界条件(ノイマンの境界条件)

  

を後退差分

  

と近似し計算したところ、うまくいかなかった。(下図参照)

 

kyoukaichi-graph-002.png

 

これは、

(1)式を差分方程式にするときに用いた

  

であるのに対し、x=1における境界条件の近似で用いた

  

の誤差がO(h²)と、(2)、(3)式に対して誤差が大きく、これが数値解の精度を悪化させたためである。

したがって、この問題を解決するためには、x=1における境界条件の近似を(4)より高精度にする必要がある。

そこで、

  

[0,1]の外に仮想の点を設け、さらに、微分方程式

  

x=1でも成立すると仮定する。

すると、x=1における差分方程式は

   

になる。

x=1における境界条件に中心差分を用いると

  

これを(5)に代入すると、

  

これと、k=1n−1に関する差分方程式

  

を合わせると、未知数の個数はn個、対して、方程式の数もnとなるので、次の連立方程式を解くことによって、を全て求めることができる。

  

この方針に従ってプログラムを書き換え、n=10h=0.1の場合について計算した結果は次の通り。

前回と異なり、良好な計算結果が得られていることがわかるだろう。

 

Kyoukaichi-Graph-001.png

 

なお、本計算に使用したプログラムは、(1)より一般の

  

を解くものであり、以下に、計算結果とプログラムを示す。

 

Keisan-Kekka-001.png 

 

! 差分法を用いて y''+f(x)y'+g(x)y=h(x) の境界値問題(ノイマン条件)を解くプログラム
! ただし、固有値問題は解けない!!
parameter (ns=100)
real x(0:ns), y(0:ns)

x = 0.; y =0; ! 変数の初期化

n=10 ! nは分割数

a=0.; b=1. ! 境界のx座標
x(0)=a; x(n)=b
y(0)=1. ! 境界条件

call Solver(x,y,n)

! 計算結果の出力
write(*,*) ' i      x         y      Exact'
do i=0, n, 1
    write(*,100) i, x(i), y(i), exact(x(i))
end do

100 format(i3,1x,f9.5,1x,f9.5,1x,f9.5)
end

function exact(x)
e=exp(1.)
exact=2/(2-e)*exp(-x)- e/(2-e)*exp(-2.*x)
end

function f(x)
f=3.
end

function g(x)
g=2.
end

function h(x)
h=0.
end

! ここより下はいじると危険
! ブラックボックスとして使うべし


subroutine Solver(x,y,n)
real a(n),b(n),c(n),d(n)
real x(0:n), y(0:n)

n1=n-1
dx=(x(n)-x(0))/n
dx2=dx*dx

do i=1,n1
    x(i)=x(0)+i*dx
end do

! 差分法によって得られる連立方程式の係数の計算
do i=1,n1
    a(i)=1-dx/2.*f(x(i))
    b(i)=-(2-dx2*g(x(i)))
    c(i)=1+dx/2.*f(x(i))
    d(i)=dx2*h(x(i))
end do

! 境界条件
    d(1)=d(1)-a(1)*y(0) ! x=0 ディリクレ条件
!    x=1 dy/dx=0 ノイマン条件
    a(n)=2.
    b(n)=(dx2*g(x(n))-2.)
    c(n)=0.
    d(n)=dx2*h(x(n))

call tdma(a,b,c,d,n) ! TDMAで連立方程式を解く

do i=1,n
    y(i)=d(i) ! 計算結果をセット
end do

end

!  TDMA
subroutine tdma(a,b,c,d,n)
real a(n), b(n), c(n),d(n)
do i=1,n-1
    ratio=a(i+1)/b(i)
    b(i+1)=b(i+1)-ratio*c(i)
    d(i+1)=d(i+1)-ratio*d(i)
end do
d(n)=d(n)/b(n)
do i=n-1,1,-1
    d(i)=(d(i)-c(i)*d(i+1))/b(i)
end do

end

 

x=0.3近傍、y=0付近での数値解と厳密解との差が多少気になるので、n=100h=0.01で計算した結果を以下に示す。

 

Keisan-Kekka-002.png

 

h=0.1のとき、小数点2桁目に、h=0.01の場合は、小数点4桁目で、数値解と厳密解が異なっているが、数値計算の精度がこの程度だから、これはしょうがない。

 

 


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

オレって実はスゴイ(笑) [ひとこと言わねば]

ルンゲ・クッタ法(4次)で境界値問題

  

が簡単に解けてしまう(・・?

 

ルンゲ・クッタを3回使うだけで、何故か分からないけれど、簡単に解けてしまったケロよ。

この解法を思いついたネムネコがすごいってことか(笑)。

 

 

何でこんなに簡単に解けるだろう。不思議だな〜。



すげぇな、4次のルンゲ・クッタ法って。



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

微分方程式の境界条件についてのよもやま話 [ひとこと言わねば]

微分方程式の境界条件についてのよもやま話

 

たとえば、次のような2階の常微分方程式があるとする。

  

この微分方程式の境界条件が

  

と与えられるものをディリクレの境界条件第1種の境界条件)、

   

と1階微分の形で与えられたものをノイマンの境界条件第2種の境界条件)などという。

 

今日の数学の記事で書いたけれど、ノイマンの境界条件をうまくプログラムに入れるのは、結構、厄介。

単に面倒なだけではなく、数値解はノイマン条件に、結構、敏感に反応するので、とんでもない結果をもたらすことがあるにゃ。うまく入れないと下の図のように数値計算と厳密解が大きく食い違ったりする。

 

 

ノイマンの境界条件をプログラムにうまく入れることが如何に大変か、

この件についてddt³さんがきっと何か素敵な記事を書いてくれると思うにゃ。みんな、期待するにゃ。

 

微分方程式の境界条件は、ディリクレ、ノイマンの境界条件の他に、混合型の境界条件第3種の境界条件)という更に厄介なものがある。

非定常1次元の熱伝導方程式

  

の場合は、

  

  

といった奴だにゃ。

(1)でも、結構、大変なのに、(2)なんかは涙がチョチョギ出るほど大変だと思うね。

(2)は、プログラムの部分的修正ではダメで、計算に使うアルゴリズム、使用する計算法などを全面的に見直し、プログラムを一から作り直さないといけないかもしれない。

そして、

新たに作ったプログラムによる計算結果が現実と大きく食い違うなんて可愛い話ではすまず、計算方法――計算に使うアルゴリズム、スキームなど――によっては、数値解が簡単に発散したり強く振動したりして、収束解が得られないなんて事態に遭遇するかもしれない(>_<)

 

⑨の定常解、つまり、

  

である、

  

の場合ですら、4次の代数方程式を解かなければ、この微分方程式は解くことができないんだから、(2)の条件を入れて解くことが如何に大変か容易に想像できるだろう。

 

差分法で解こうとする場合、最も簡単な方法でも、境界条件(2)から

  

という4次方程式が出てくる!!

単体の方程式(3)だけならば、ニュートン法や2分法を用いて、この近似解を求めることができるだろうが、

  

というn−1本の線形の連立方程式(4)に非線形の方程式(3)を同時に解かないといけないのだから、これはそんなに簡単な話ではない。

ひょっとしたら、(3)と(4)からなる連立方程式系は実数解をもたないかもしれない(^^

 

 

4次方程式(3)は上のグラフのような状態になっており、解けないかもしれない(^^ゞ

 

 

2種や3種の境界条件をうまく計算プログラムに入れるってのは、それほど、大変なんだケロ。

そして、多くの数値解析、数値計算の(入門的な?)本には、こういった条件の取り扱い方については、ほとんど書いていないのであった。

 

英語で書かれた「Computer Analisys」といったようなタイトルの1万円くらいの、結構、専門的な英語の本であっても、この手の話はほとんど出ていないと思うよ。

 

ってことは、ネムネコが、いま、書いている記事には、1万円以上の価値があるってこと(・・?

これらの記事を読んだヒト、一人から300円ではなく千円くらいのお賽銭をもらってもバチは当たらないに違いない(笑)。



カッパネットのように「お値段異常」ではなく、適正価格だと思うにゃ。



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

今日のアニソン、東方から『幻想河童行進曲』 [今日のアニソン]

今日のアニソンは、東方から『幻想河童行進曲』です。


河童の「にとり」は、可愛いにゃ。



「詐欺」つながりで、この曲も埋め込んでおこう!!



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