266びぼーろく;maximaでの二分法

数値計算でニュートン法を使おうと思ったら全然動かない。

今回は変数も少ないので二分法に挑戦。

maximaにはパッケージがないので、自分で作るしか無い。

面倒なのでgoogle先生に教えを請うと、すぐに答えが返ってきた。

まず、予め解きたい方程式をfとして与える。
次に解の範囲をgnuplot 等で調べ、解xの範囲をa,bとして与える。
さらに求めたい解の精度をrで与える。

以下はそのコピペ。

bs(f,a,b,r):=block(
     f(x):=f,
     c:a,
     for i:1 while abs(subst(c,x,f(x)))>r do(
     c:(a+b)/2,
     if subst(c,x,f(x))=0 then return(float(c)),
     if (subst(a,x,f(x)))*(subst(c,x,f(x)))<0 b:c="" p="" then="">     else a:c
     ),
     return(float(c))
     );

うまく動いたので感動。ありがとうございました。

多変数の二分法ってあるんだろうか?


さて、あまりに寒いので同僚とおでんを食べに行った。せりが美味しかった。
ついでに二次会で海鮮丼を食べていたので激写。よく食べるなぁ。

さらに家に帰ると室温が5度を切っていたので、コタツを出した。ぬくぬくで幸せ。