4章 数値計算とコンピュータ

§1 簡単なプログラム
プログラム:program ギリシャ語「公然と書く」の意から
BASIC:Beginner's All-purpose Symbolic Instruction Code の略とはできすぎな話
§2 いろいろなアルゴリズム

二分法:Bisection method
台形公式:Trapezium rule trap(罠)とかtrapeze(空中ブランコ)とは関係ないか
平方の和とか約数とか最大公約数はそのまま関数が用意されているし

(%i1) load(functs);

(%i2) wa(M):=nusum(k^2,k,1,M);
Result

(%i3) wa(10);
Result

(%i4) wa(M);
Result

(%i5) yak(A):=divisors(A);
Result

(%i6) yak(100);
Result

台形公式をMaximaでプログラミングしてみました。
関数と区間と分割数を入力して面積を台形公式で計算する。
収束は遅い方です。

(%i7) trap(f,a,b,n):=block(
    f(x):=f,
    return(float((b-a)/n/2*(a+b+2*nusum(subst((b-a)/n*k,x,f(x)),k,0,n-1))))
);

Result

(%i8) for n:1 thru 10 do print(trap(x^2,0,1,n));
Result

2分法で方程式の解を近似してみると,思ったより収束は早い。

(%i9) 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 then b:c
        else a:c
    ),
    return(float(c))
);

Result

(%i10) for i:0 thru -10 step -1 do print(bs(x^2-2,1,2,10^i));
Result

inserted by FC2 system