aaaa - suchikeisan @ ウィキ

>

#1

• 1

restart;
with(plots):
y:=x->4*sqrt(1-x^2);
a:=0;
b:=1;
s:=evalf(int(y(x),x=a..b));
(1)
>
NumIntTpz:=proc(f,N,a,b)
local S, h, xi, i;
h:=(b-a)/N:
S:=f(a)/2;
for i from 1 to N-1 do
xi:=a+i*h;
S:=S+f(xi);
end do:
S:=S+f(b)/2;
return S*h;
end proc:
dX2:=[]:
>
for
i from 1 to nops(N) do
dX2:=[op(dX2),abs(s-NumIntTpz(y,N[i],a,b))];
end do:
>
lgp2:=loglogplot(N,dX2,color=blue):
>

#2

• 1

restart;
Answer5:=0.78541*(23.129-23.001);
Answer4:=0.7854*(23.12-23.00);
Answer3:=0.785*(23.1-23.0);
Answer2:=0.78*(23-23);
(2)
>

#2

• 2

数値計算の「つぼ」となるのは「誤差」と「安定性」である。
この「誤差」には「丸め誤差」と「打ち切り誤差」が存在する。
浮動小数点数間の算術演算では少なくともεmの相対誤差が生じる。
(εmとは浮動小数点数1.0を加えたときに1.0と異なる結果になるような最小の浮動小数点数のことである)
この誤差のことを「丸め誤差」という。
真の答えと現実の計算による答えとの食い違いのことを「打ち切り誤差」という。
これはすべてのプログラマの制御下に存在している。
これらの誤差についてきちんと心得ておくことが数値計算のつぼであるといえる。
それに加えて、解が見つかる保証である安定性が重要となってくる。
しかし、安定性は問題や解法、初期値に依存するということに注意しなければならない。
>

#3

• 1

restart;
with(LinearAlgebra):
RANK:=Matrix(1,1,0],[1,0,1],[0,1,1);
(3)
>

#3

• 2

大阪のトラックは各月に半分が名古屋に行き、残りは域内で使用されるとあるので
翌月には
大阪：15台
名古屋：15台
東京：0台
となっているはずである。
>

#3

• 3

>
restart;
with(LinearAlgebra):
RANK:=Matrix(1,1,0],[1,0,1],[0,1,1);
(4)
>
TR:=Transpose(RANK);
(5)
>
V:=[]:
for i from 1 to 3 do
s:=0;
for j from 1 to 3 do
s:=s+TR[j][i];
end do:
V:=[op(V),1/s];
end do:
V;
VA:=DiagonalMatrix(Vector(V),3,3);
PAGE:=TR.VA;
(6)
>
e,v:=evalf(Eigenvectors(PAGE)):
>
maxnum:=0: maxi:=1: s:=0:
for i from maxi to 3 do
if (maxnum < abs(e[i])) then
maxnum:=abs(e[maxi]);
maxi:=i;
end if;
end do:
PAGERANK:=evalf(Column(v,maxi)/add(Column(v,maxi)[i],i=1..3));
(7)
>
ans:=Vector(1/3,3):
N:=100: eps:=10^(-8):
for i from 1 to N do
d:=abs(abs(ans)-abs(PAGE.ans));
delta:=add(d[j],j=1..3);
if (delta<eps) then
break;
end if;
ans:=PAGE.ans;
end do:
evalf(ans);  
(8)
>

#4

二次曲線
restart;
with(LinearAlgebra):
with(plots):
with(plottools):
X:=[0,1,3,6,12];
Y:=[50,45,37,27,15];
y1:=[seq(point([X[i],Y[i]],symbolsize=30),i=1..5)]:
display(y1);
>
DM:=proc(F,X)
local A, n, m, i, j;
n:=5;
m:=3;
A:=Matrix(1..n,1..m);
for i from 1 to n do
for j from 1 to m do
A[i,j]:=F[j](X[i]);
end do:
end do:
return A;
end proc:
>
solveNE:=proc(A,X,Y)
return MatrixInverse(Transpose(A).A).Transpose(A).Vector(Y);
end proc:
>
n:=3:
F:=[x->1,x->x,x->x^2];
A:=DM(F,X);
a:=solveNE(A,X,Y);
(9)
>
f:=unapply(add(a[i]*F[i](x),i=1..3),x);
f1:=plot(f(x),x=X[1]..X[5]):
display(y1,f1);
>

#4

直線
restart;
with(LinearAlgebra):
with(plots):
with(plottools):
X:=[0,1,3,6,12];
Y:=[50,45,37,27,15];
y1:=[seq(point([X[i],Y[i]],symbolsize=30),i=1..5)]:
display(y1);
>
DM:=proc(F,X)
local A, n, m, i, j;
n:=5;
m:=2;
A:=Matrix(1..n,1..m);
for i from 1 to n do
for j from 1 to m do
A[i,j]:=F[j](X[i]);
end do:
end do:
return A;
end proc:
>
solveNE:=proc(A,X,Y)
return MatrixInverse(Transpose(A).A).Transpose(A).Vector(Y);
end proc:
>
n:=2:
F:=[x->1,x->x];
A:=DM(F,X);
a:=solveNE(A,X,Y);
(10)
>
f:=unapply(add(a[i]*F[i](x),i=1..2),x);
f1:=plot(f(x),x=X[1]..X[5]):
display(y1,f1);
>