« 妖怪ラーメンは手作り感満載 | トップページ | 伏見ででかい釜を見た。 »

2008年12月17日 (水)

古今算法記の問題4をkeisan.casio.jpで解いてみる。

数学セミナー2008年12月号の「現代数学から見た関孝和」に出てた問題。

Seki

をx,y,zについて解くというもの。未知数を消去するとz^1/4の108次式になるそうです。でこれExcelのソルバーで解けないかな?とおもってやったら全然値が違う。あー、桁数が全然足りてないんだ。ということでここはkeisan.casio.jpでやってみよう(これはUPしてもしょうがないかな)。

解くと、30桁で

x=36.8641465377853041287298140705
y=44.3516761876674573034940403248
z=32.5563821163811031295834914427

x=51.4570363183687277311149213375
y=10.293599457590667227853861756
z=49.4144096308403396077899331929

の2つが得られた。解き方はニュートン・ラフソン法。でも微分して逆行列を求めるのは死ぬほどめんどくさい。なのでやっぱりMaximaでやってみた。keisan.casio.jpに対応したプログラムはこちら。ね?めんどくさいでしょ。

for(k=0;k<=1;k=k+1){
  if (k==0) {
    x=36;
    y=44;
    z=32;
  } else {
    x=51;
    y=10;
    z=49;
  }
  for(j=0;j<=100;j=j+1) {
    x2=x-(-(3*sqrt(x)*y^(8/3)*(z^3+y^3-121750))/((12*y^(8/3)-8*x^(5/2))*z^(11/4)+6*x^(5/2)*y^(8/3))-((y^3+x^3-137340)*(4*sqrt(x)*z^(11/4)-3*sqrt(x)*y^(8/3)))/((12*y^(8/3)-8*x^(5/2))*z^(11/4)+6*x^(5/2)*y^(8/3))+(12*sqrt(x)*y^(8/3)*(z^(1/4)+y^(1/3)+sqrt(x)-12)*z^(11/4))/((6*y^(8/3)-4*x^(5/2))*z^(11/4)+3*x^(5/2)*y^(8/3)));
    y2=y-((3*x^(5/2)*y^(2/3)*(z^3+y^3-121750))/((12*y^(8/3)-8*x^(5/2))*z^(11/4)+6*x^(5/2)*y^(8/3))-(12*x^(5/2)*y^(2/3)*(z^(1/4)+y^(1/3)+sqrt(x)-12)*z^(11/4))/((6*y^(8/3)-4*x^(5/2))*z^(11/4)+3*x^(5/2)*y^(8/3))+(3*y^(2/3)*(y^3+x^3-137340)*z^(11/4))/((6*y^(8/3)-4*x^(5/2))*z^(11/4)+3*x^(5/2)*y^(8/3)));
    z2=z-(((3*y^(8/3)-2*x^(5/2))*z^(3/4)*(z^3+y^3-121750))/((6*y^(8/3)-4*x^(5/2))*z^(11/4)+3*x^(5/2)*y^(8/3))+(12*x^(5/2)*y^(8/3)*(z^(1/4)+y^(1/3)+sqrt(x)-12)*z^(3/4))/((6*y^(8/3)-4*x^(5/2))*z^(11/4)+3*x^(5/2)*y^(8/3))-(3*y^(8/3)*(y^3+x^3-137340)*z^(3/4))/((6*y^(8/3)-4*x^(5/2))*z^(11/4)+3*x^(5/2)*y^(8/3)));
    x=x2;
    y=y2;
    z=z2;
  }
println(x2,y2,z2);
}

« 妖怪ラーメンは手作り感満載 | トップページ | 伏見ででかい釜を見た。 »

学問・資格」カテゴリの記事

コメント

コメントを書く

(ウェブ上には掲載しません)

トラックバック

この記事のトラックバックURL:
http://app.cocolog-nifty.com/t/trackback/512682/43450068

この記事へのトラックバック一覧です: 古今算法記の問題4をkeisan.casio.jpで解いてみる。:

« 妖怪ラーメンは手作り感満載 | トップページ | 伏見ででかい釜を見た。 »

最近のコメント

2018年11月
        1 2 3
4 5 6 7 8 9 10
11 12 13 14 15 16 17
18 19 20 21 22 23 24
25 26 27 28 29 30  
フォト
無料ブログはココログ