以下 Tutorial Chapter 2 に従います。

Windows の操作

対話的操作

ヘルプ

数値計算


1 + 5/2

(1 + (5/2 * 3)) / (1/7 + 7/9)^2

^ -- 指数演算子

1234^123

100!

fact, ! -- 階乗

isprime(123456789)

TRUE|FALSE -- 真偽値

isprime -- 素数判定
           Miller-Rabin test (確率的素数判定法の一種)
           FALSE -> 確実に合成数
           TRUE -> 素数、もしくは(「ランダム」に選んだ10個の底に対する)強偽素数
           ref. numlib::proveprime

あまり長い桁の数値を扱うと…


a := 2^(2^13)+1
length(a)
isprime(a)


b := 2^(2^17)+1
length(b)
isprime(b)

タスクマネージャーでプロセスを見てみる

ifactor -- 素因数分解


ifactor(123456789)
3^2 * 3607 * 3803

a := 2^(2^5)+1; length(a); ifactor(a)

b := 2^(2^7)+1
length(b)
ifactor(b)

"Session" -> "STOP Kernel"

「厳密な」数値計算


sqrt -- 平方根


sqrt(56)
sqrt(14)^4


limit -- 極限値を求める
infinity -- 正の無限大
exp -- 指数関数
ln -- 自然対数
E -- 自然対数の底


limit((1 + 1/n)^n, n = infinity)
ln(1/exp(1))

「非厳密な」数値計算


float -- 浮動小数点形式への変換(近似値計算)
DIGITS -- 浮動小数点形式への有効桁数(初期値は10)


float(sqrt(56))
DIGITS; float(67473/6728)

DIGITS := 100; float(67473/6728); DIGITS := 10


(1 * (5/2 * 3)) / (1/7 + 7/9)^2    <- 1
(1.0 * (5/2 * 3)) / (1/7 + 7/9)^2  <- 1.0


2/3 * sin(2); 0.6666666 * sin(2)


float(2/3 * sin(2)); 0.6666666 * float(sin(2)) <- ≠

「厳密(exact)な値(無限精度)」と「非厳密(inexact)な値(有限精度)」の区別


sqrt(56.0); sin(3.14)

cos(PI); ln(E)

DIGITS := 100; float(PI); float(E); DIGITS := 10


Excecise 2.2:
√27 - 2√3 と cos(π/8) の「厳密な」値は?

√27 - 2√3 -> 3√3 - 2√3 -> √3


sqrt(27)-2 sqrt(3)
float(%)


cos(A/2) = √{(1+cos(A))/2}

cos(π/8)  -> cos(π/4/2)
           -> √{(1+cos(π/4))/2}
           -> √{(1+1/√2)/2}
           -> {√(2+√2)}/2

cos(PI/8)
float(%)

複素数


虚数単位: I


sqrt(-1), I^2

(1+2*I)*(4+I), (1/2+I)*(0.1+I/2)^3

1/(sqrt(2)+I)


必ず x+yI という形で答を返してくれるわけではない


z := 1/(sqrt(2) + I)
rectform(z)
Re(z)
Im(z)
abs(z)
arg(z)
conjugate(z)


rectform: x+yI という形に整頓した式を返す
Re: 実数部 (x+yI -> x)
Im: 虚数部 (x+yI -> y)
abs: 絶対値
arg: 偏角
conjugate: 共役 (x+yI -> x-yI)

数式処理


reset()

例: x と y の多項式を定義して f と置く

f := y^2 + 4*x + 6*x^2 + 4*x^3 + x^4


例: f を x で微分

diff(f,x)


例: f を y で微分

diff(f,y)


例: f を x で3階微分

diff(diff(diff(f,x),x),x)
diff(f,x,x,x)


D
f' は D(f)の省略記法

例: sin を微分

sin'
sin'(x)
D(sin)
D(sin')


diff と D の違い ->荒っぽく言うと定義域が異なる

diff は式に対して式を返す
ex. x+x^2 -> 2*x+1
D は関数に作用して関数を返す(微分作用素)
ex. (x -> x+x^2) -> (x -> 2*x+1)

現時点ではあまり差がはっきりしないでしょう
  これは数学の知識・学習度の問題
  (紙と鉛筆と頭を使った)数学のトレーニングが大事

余談
  プログラミング・実装的観点からの違い
  diff は builtin 関数
  D はライブラリ関数

  expose ... procedure (など)のソースコードを表示する

  expose(diff)
  expose(D)

  数学で定義・定理を積み重ねていくのと同じこと

  expose(sin)
  expose(abs)

  ...


例: f を x=0..1 まで定積分

int(f, x=0..1)


例: f を x で不定積分

int(f,x)


例: f を y で不定積分

int(f,y)


必ずうまいこと「綺麗な」形で結果が帰ってくるわけではない
  「初等関数」で表現できるか否か
int(1/sqrt(x^2-1), x)
int(1/sqrt(x^3-1), x) <- 楕円積分(でよかったっけ)


例: 不定積分・定積分

g := 1/(exp(x^2)+1)
h := int(g, x)
diff(h, x)
s := int(g, x=0..1)
float(s)
plotfunc2d(g, x=0..1)


MuPAD で定義されている関数一覧
  see Quick Reference p11 (Chapter 8)


例: 各々の関数についての特別な値を知っている

cos(0)
sin(PI/2)
exp(0)
ln(1)


例: 特別な値でないときはそのままの形を返す

sqrt(2)
exp(1)
sin(x+y)


例: 式を展開する(expand)
  「展開」の定義は…
    (たとえば)和積を積和に変形してくれる


(x+1)*(x-1)
expand(%)

exp(x+y)
expand(%)

sin(x+y)
expand(%)

tan(x+3*PI/2)
expand(%)


例: 式を正規化する(normal)
   「正規化」の定義は…
    (たとえば)既約分数の形に整頓してくれる

normal((x^2-1)/(x+1))
f := x/(1+x) - 2/(1-x): g := normal(f)
normal(x^2/(x+y) - y^2/(x+y))


例: 部分分数展開(partfrac)

partfrac(g)
partfrac(g, x)

f := x^2/(x^2-y^2)
partfrac(f,x)
partfrac(f,y)


例: 簡約(simplify)
    「より単純な」式へ変形してくれる
    かなり複雑なことをしてくれる(詳細は略)

(exp(x)-1)/(exp(x/2)+1)
simplify(%)

simplify(sin(x)^2+cos(x)^2)


例: 巾乗根についての簡約(radsimp)

f := sqrt(4+2*sqrt(3))
radsimp(%)

x := 1/2 + sqrt(23/108);
y := x^(1/3)+1/3/x^(1/3);
z := y^3 - y;
radsimp(z)


例: 因数分解 (factor)

reset()
factor(x^3 + 3*x^2 + 3*x + 1)
factor(2*x*y - 2*x - 2*y + x^2 + y^2)
factor(x^2/(x+y) - z^2/(x+y))

# Exercise 2.3:

expand((x^2+y)^5)
factor(%)


例: 式の書換え (rewrite)

   rewrite(f, t)   t = tan | sincos | exp | ln |...


rewrite(cot(x), sincos)
rewrite(sin(x), tan)
rewrite(exp(I*x), sincos)
rewrite(tan(x), exp)


例: 極限値を求める (limit)

limit(sin(x)/x, x=0);

# Exercise 2.5:

f := sin(x)/x;
limit(f, x=0);
plotfunc2d(f, x=-10..10)

  特異点 x=0 の処理もうまくやってくれる(ときもある)

f := (1-cos(x))/x;
limit(f, x=0);
plotfunc2d(f, x=-PI..PI)

f := ln(x);
limit(f, x=0, Right);
plotfunc2d(f, x=0..1)

f := x^sin(x);
limit(f, x=0);
plotfunc2d(f, x=0..PI/2)

f := (1+1/x)^x;
limit(f, x=infinity);
plotfunc2d(f, x=0..1000)

f := ln(x)/exp(x);
limit(f, x=infinity);
plotfunc2d(f, x=0..10)

f := x^ln(x);
limit(f, x=0, Right);
plotfunc2d(f, x=0..5)

f := (1+PI/x)^x;
limit(f, x=infinity);
plotfunc2d(f, x=0..1000);
float(exp(PI))

f := 2/(1+exp(-1/x));
limit(f, x=0, Left);
plotfunc2d(f, x=-3..0)

f := sin(x)^(1/x); ← x in [-PI/2, PI/2]
limit(f, x=0);
limit(f, x=0, Right);
limit(f, x=0, Left);
plotfunc2d(f, x=0..PI/2)


例: 関数定義 (->)
    プログラミングとも関連(詳細は略)
    仮引数

F := x -> x^2
F(x), F(y), F(a+b), F(a+b), F'(x)

G := y -> y^2

F と G は「同じ」関数
本当は ■ -> ■^2 と書きたい:-)


例: 数式 (expression)
    数学系で一般的に使われる「式」と考えてよい
    きちんとした定義は略、例示のみ
    ":=" と "=" の違いに注意!

0
1
1/2
PI
x+1
x+y
x^3+y+z
sin(x)+ln(x)
x+y+z = 1
x+y+z > 1
1 = 2
x <> 0
3 <= 1
x -> x^2


例: 集合 (set, {})
    数学系で一般的に使われる「集合」と考えてよい
      要素の重複に意味はない
      要素の順序に意味はない

{}   空集合

{1,2,3,4,a,b,1,2,a}
{a,b,1,2,3,4} 
{x+y+z=1, x-y-z>1}


例: 等式系・不等式系を解く (solve)
solve(式)
solve(式, 変数)
solve(式集合, 変数集合)
  
  例: 連立一次方程式を解く
  e := { x+y=a, x-a*y=b };
  u := { x, y };
  solve(e,u)
  
  例: xについて解く
  solve(x^4 - 5*x^2 + 6*x = 2, x)
  solve(x^2-2*x+2 = 0, x)
  
  例: (暗黙にxについて)解く
  solve(x^2 = 3)
  solve(x^2-2*x+2 = 0)
  
  例: xについて解く
  solve(sin(x*PI/7) = 0, x)
  
  例: xについて解く
  solve(cos(PI*x) = 0)← 簡約が足らなそう
  normal(%)
  
  例: xについて解く(不等式なので x の範囲が求まる)
  solve(x^2 > 5, x)
  
  例: xについて解く(不等式なので x の範囲が求まる)
  solve(x^2 <> 7, x) == solve(x^2 =  7, x) の補集合
  
  例: 二次方程式の解の公式
      piecewise ... 条件分け
      solve(a*x^2+b*x+c=0, x)


例: 総和 Σ (sum)
sum(i, i=0..100) ← ガウスの伝記によく出てくる計算
sum(j, j=0..100) ← i でも j でも同じ

sum(i, i=0..n)
factor(%)
sum(i^2, i=0..n)
factor(%)
sum(i^3, i=0..n)
factor(%)
sum(i^k, i=0..n) ← ここまでの一般化はちょっと無理そう

# Excersice 2.7

sum(k^2+k+1, k=1..n);
factor(%)

sum((2*k-3)/((k+1)*(k+2)*(k+3)), k=1..infinity)

sum(k/(k-1)^2/(k+1)^2, k=2..infinity)


例: 累積 Π (product)

product(i, i=1..10)
10!
product(i, i=1..n) ←  Γ(n+1) ... ガンマ関数
                       n! =  Γ(n+1) 

product((1-1/i), i=2..10) 
product((1-1/i), i=2..100)

おもしろくもない計算だが 1/2 * 2/3 * 3/4 * ... (n-1)/n = 1/n なので
product((1-1/i), i=2..n) == 1/n になる(n>1)

product((1-1/i), i=2..n) ← ???

行列


例: 行列の生成 (matrix)
    詳細は ?matrix
A := matrix([[1,2],[a,4]])
B := matrix([[y,3],[z,5]])
A
B

例: 行列演算 (+, -, *)
A+B
A-B
A*B

例: 逆行列の生成・演算 (^-1, /)
A^-1
A^(-1)
1/A
A*B^(-1)
A/B
matrix([[x,y,z],[x,y,z],[z,y,x]])
%^(-1)

例: 行列式の計算 (linalg::det)
    LINear ALGebra DETerminant
    "::" の意味は後述予定(プログラム関連)
linalg::det(A)

例: 行ベクトル (n * 1 行列)
matrix([1,x])   角括弧の数の違いに注意

例: 線形方程式 Av=b を解く (A の逆行列を利用)
v := A^(-1) * b

例: 行列の各要素に対する数式操作
map(v, normal)

例: 検算
A * %
map(%, normal) == b