#freeze
[[FrontPage]]

#contents

2009/11/14からのアクセス回数 &counter;

このページで紹介したsageノートブックは、以下のURLから参照することができます。

http://www.sagenb.org/home/pub/1048/

** 変数の宣言と値の代入 [#l61101df]

sageの式で変数として認識させるには変数をvar関数で宣言しなくてなりません。

変数の宣言は、

 var('変数名')
複数の変数を宣言する場合には、スペースを空けて指定します。 宣言される変数を参照したりする場合には、
 x, y = var('x y')
のようにpython変数x, yに宣言したsage変数を代入します。

変数の値を指定して、関数の返す値を調べるには、

 f1(x=1.5)
のように式を表す変数にカッコを付けて、x=1.5と値を指定するとその時の関数の値が表示されます。

sage入力：
#pre{{
x = var('x')
f1 = (x - 1)*(x^2 - 1)
f1(x=1.5)
}}
   	
&color(blue){0.625000000000000};

変数をx, yと複数存在する場合の例をです。f5に\(f5(x,y)=x2+y\)を定義し、f5(x=2, y=3)でx=2, y=3の時の値を出力します。

sage入力：
#pre{{
var('x y w')
f5=x^2 + y
f5(x=2, y=3)
}}
&color(blue){7};

** 規則の代入 [#u0590089]
規則の代入には、subs_exprを使います。

sage入力：
#pre{{
f6=x^2 + 1; f6
}}
&color(blue){x^2 + 1};

sage入力：
#pre{{
f6.subs_expr(x^2 == w)
}}
&color(blue){w + 1};
	
** 関数の定義 [#z365a917]
sageで関数を定義する方法が２つあります。
- pythonの関数を定義する
- lambda式を使って関数を定義する

です。

*** python関数の定義 [#kdb65a09]
手続き的な関数を定義する場合には、python関数を使うと便利です。

pythonでは

#pre{{
def 関数名(引数のリスト):
    関数の本体
}}

のように定義します。

以下にx3を返す関数cubeをpythonの関数を使って定義しました。

sage入力：
#pre{{
def cube(x):
    return x^3

cube(2)
}}
&color(blue){8};

#pre{{
x, y = var('x y')
cube(x+y)
}}
&color(blue){(x + y)^3};

*** lambda関数の定義 [#cb427343]
つぎにlambda式を使って関数を定義する方法を示します。 x3のように式で表現できる場合には、lambda関数が便利です。
 lambda 変数のリスト: 式
のように定義します。 先ほどと同じx3を返す関数の例を以下に示します。

sage入力：
#pre{{
f = lambda x: x^3
f(2)
}}
&color(blue){8};

#pre{{
f(x+y)
}}
&color(blue){(x + y)^3};

** 方程式の解法 [#a78e7648]
数式処理で助かる機能は、方程式の解法でしょう。

方程式の解法は、solve関数を使って以下のように行います。

 solve(解放したい関数また関数のリスト, 解を得る変数)

例として、\(f(x)=ax+b\)の一次方程式を解いてみます。結果から\(x=−b/a\) が得られたので、\(a=2,b=1\)を代入して解の値を取得します。

sage入力：
#pre{{
# 方程式の解法
var('x a b')
f = a*x + b

sln = solve(f, x); sln
}}
&color(blue){[x == -b/a]};

#pre{{
x = -b/a
x(a=2, b=1)
}}
&color(blue){-1/2};

これだと、解から式をもう一度入力しなくてはなりません。そこで、a,bも変数として関数を定義し、solveでa=2,b=1をセットするようにします。こうすれば簡単に解の値が確認できます。

sage入力：
#pre{{
var('x a b')
g = lambda x, a, b : a*x + b

solve(g(x, a, b), x)
}}
&color(blue){[x == -b/a]};

#pre{{
solve(g(x, a=2, b=1), x)
}}
&color(blue){[x == (-1/2)]};

解の値がほしいときには、find_root関数で数値解析で値を計算します。
 find_root(関数, 計算する変数の開始値、終了値)
先ほどの関数g(x):a=2,b=1の解をx=−1,1 の範囲で数値解析します。

sage入力：
#pre{{
find_root(g(x, a=2, b=1), -1, 1)
}}
&color(blue){-0.5};

*** 多項式の解 [#u06bdee2]
一次方程式ではあまりに簡単なので、多項式の解を求めてみます。
$$
f(x)=x^3−2x+4 
$$

のプロット図を以下に表示します。x軸と交差しているのはx=−2 です。-1と1付近に極値があります。

sage入力：
#pre{{
var('x')
f = x^3-2*x+4
plot(f, -2.5, 2.5)
}}
&ref(sage0.png);

solveの結果、\(x=−2,x=1−i,x=1+i \)を得ました。
sage入力：
#pre{{
solve(f, x)
}}
&color(blue){[x == (-I + 1), x == (I + 1), x == -2]};

sage入力：
#pre{{
factor(f)
}}
&color(blue){(x + 2)*(x^2 - 2*x + 2)};

plotで可視化することでfind_rootでの計算範囲を大まかに把握することができます。

図から-3から3の間に解があることが見て取れたので、この範囲で数値解を求めます。
sage入力：
#pre{{
find_root(f, -3, 3)
}}
&color(blue){-1.9999999999999949};

*** 連立方程式の解 [#q37df93e]
同様にsolveに関数に式のリストを与えることで、連立方程式の解を得ることができます。

$$
\begin{eqnarray}
x^2 - 2y &=& 2 \\
x - y &=& 1
\end{eqnarray}
$$

の解を例に以下に示します。

sage入力：
#pre{{
var('x y')
f = [y == 1/2*x^2 - 1, y == x - 1]
solve(f, x, y)
}}
&color(blue){[ [x == 2, y == 1], [x == 0, y == -1] ]};

上記関数のプロット結果は、以下のとおりです。

sage入力：
#pre{{
p1 = plot(lambda x : 1/2*x^2 - 1, (x, -1, 3))
p2 = plot(lambda x : x - 1, (x, -1, 3))
(p1+p2).show()
}}
&ref(sage0-1.png);

もう一つ、

$$
\begin{eqnarray}
x^2 + y^2 &=& 1 \\
y &=& x - 1
\end{eqnarray}
$$

の解を例に以下に示します。

sage入力：
#pre{{
var('x y')
f = [x^2 + y^2 == 1, y == x - 1]
solve(f, x, y)
}}
&color(blue){[ [x == 1, y == 0], [x == 0, y == -1] ]};

この関数をimplicit_plotを使って関係式の形をそのまま使ってプロットしてみます。
円のなので、aspect_ratio=1としました。

sage入力：
#pre{{
var('x y')
p1 = implicit_plot(x^2 + y^2 == 1, (x, -1, 1), (y, -1, 1))
p2 = plot(lambda x : x - 1, (x, -1, 1))
(p1+p2).show(aspect_ratio=1)
}}
&ref(sage0-2.png);

** 導関数 [#v361bc68]
高校で習った、導関数をsageを使って導いてみましょう。

関数\(f(x)=\frac{1}{2} x^3\)とします。

平均変化率は、\(g(x)=hf(x+h)−f(x)\)を以下のように定義します。

sage入力：
#pre{{
# 導関数
def f(x):
    return (x^3)/2

var('x, h')
g=(f(x + h) - f(x))/h; g
}}
&color(blue){1/2*((h + x)^3 - x^3)/h};

g(x)を展開し、整理すると

sage入力：
#pre{{
g1= g.rational_simplify(); g1
}}
&color(blue){1/2*h^2 + 3/2*h*x + 3/2*x^2};

となります。

ここで、h→0 の極値を取ったものが求める導関数\(f′(x)\)です。

sage入力：
#pre{{
limit(g1, h=0)
}}
&color(blue){3/2*x^2};

結果は、f(x)をxで微分したものと同じになります。

sage入力：
#pre{{
diff(f(x), x)
}}
&color(blue){3/2*x^2};

** 微分 [#gc56b461]
高校で習う微分の公式は、

- \(\frac{d}{dx} c = 0\)
- \(\frac{d}{dx} x^n = n x^{x-1}\)

です。

sageの微分diff関数は、以下のように使います。

 diff(関数, 微分する変数)
上記公式のsageで結果は、以下のようになります。

sage入力：
#pre{{
# 微分
x, c, n = var('x c n')
f = function('f', x)
g = function('g', x)

print diff(c, x)
print diff(x^n, x)
}}
&color(blue){0};

&color(blue){n*x^(n - 1)};

sageで変数xの関数fを以下のように定義し微分すると、\(f′\)は以下のようにD[0](f)(x)となります。

#pre{{
x = var('x')
f = function('f', x)
}}

これで変数xの関数fを定義したことになります。

sage入力：
#pre{{
diff(f, x)
}}
&color(blue){D[0](f)(x)};

このD[0](f)を\(f′\)と読み替えると、

sage入力：
#pre{{
print diff(c*f, x)
print diff(f+g, x)
print diff(f*g,x)
print diff(f/g, x).simplify_full()
print diff(1/g,x)
print diff(f(g), x)
}}
&color(blue){c*D[0](f)(x)};

&color(blue){D[0](f)(x) + D[0](g)(x)};

&color(blue){D[0](f)(x)*g(x) + f(x)*D[0](g)(x)};

&color(blue){(D[0](f)(x)*g(x) - D[0](g)(x)*f(x))/g(x)^2};

&color(blue){-D[0](g)(x)/g(x)^2};

&color(blue){D[0](f)(g(x))*D[0](g)(x)};

は、以下のようようになり、微分の各公式を表しています。
- \( (cf(x))′=cf′(x)\)
- \( (f(x)+g(x))′=f′(x)+g′(x)\)
- \( (f(x)g(x))′=f′(x)g(x)+f(x)g′(x)\)
- \( f(x)/g(x)=(f(x)′g(x)−g′(x)f(x))/g(x)^2 \)
- \( (1/g(x))′=−g′(x)/g(x)^2\)
- \( (f(g(x)))′=f′(g(x))g′(x)\)

*** いろんな関数の微分 [#l9d4ca27]
以下の関数を微分した結果を示します。
- \( (x^2+1)^3\)
- \( sin(x)^3\)
- \( sin^{-1}(x)\)
- \( e^x \)
- \( log(x)\)

sage入力：
#pre{{
# いろんな関数の微分
print diff((x^2+1)^3).factor()
print diff(sin(x)^3, x)
print diff(arcsin(x), x)
print diff(exp(x), x)
print diff(log(x), x)
}}
&color(blue){6*(x^2 + 1)^2*x};

&color(blue){3*sin(x)^2*cos(x)};

&color(blue){1/sqrt(-x^2 + 1)};

&color(blue){e^x};

&color(blue){1/x};

*** 微分の応用 [#ob786c66]
微分の応用として、接線の方程式とその法線を計算してみましょう。\( f(x)=e^{−x^2}\)の接線は、

$$
y−y_1=f′(x_1)(x−x_1) 
$$

となります。 また、法線は、
$$
y−y_1=\frac{1}{f′(x1)}(x−x_1) 
$$
となります。

微分した式をgとし、x=1での接線と法線を計算し、プロットしてみます。

sage入力：
#pre{{
f = exp(-(x^2)); f
}}
&color(blue){e^(-x^2)};

#pre{{
g = diff(f, x); g
}}
&color(blue){-2*x*e^(-x^2)};

#pre{{
m = g(x=1); m
}}
&color(blue){-2*e^(-1)};

#pre{{
p1 = plot(f, (x, -2, 2))
p2 = plot(m*(x-1)+f(1), (x, -2, 2))
p3 = plot(-(x-1)/m+f(1), (x, -2, 2))
pt = point([1, f(1)])
(p1+p2+p3+pt).show(aspect_ratio=1, ymin=-1, ymax=2)
}}
&ref(sage0-3.png);

** 積分 [#p3cac054]
sageでは積分は、integrate関数で計算します。

#pre{{
integrate(被積分関数, 積分変数)
また
integrate(被積分関数, 積分変数, 積分範囲)
}}
とすると定積分を計算します。

\(\int xdx\)を計算した後、その微分を求めています。

sage入力：
#pre{{
f = integrate(x); f
}}
&color(blue){1/2*x^2};

#pre{{
diff(f)
}}
&color(blue){1};

以下の関数を積分した結果示します。

- \( \int cdx \)定数c
- \( \int \frac{1}{x}dx\)
- \( \int e^x dx\)
- \( \int log(x)dx \)
- \( \int sin(x)dx \)
- \( \int \frac{f'(x)}{f(x)}dx\)

sage入力：
#pre{{
# いろんな関数の積分
x = var('x')
c = var('c')
f = function('f', x)
g = function('g', x)

print integrate(c, x)
print integrate(1/x, x)
print integrate(exp(x), x)
print integrate(log(x), x)
print integrate(sin(x), x)
print integrate(diff(f, x)/f(x), x)
}}
&color(blue){c*x};

&color(blue){log(x)};

&color(blue){e^x};

&color(blue){x*log(x) - x};

&color(blue){-cos(x)};

&color(blue){log(f(x))};

*** 定積分 [#oa2182e4]
次に定積分\( \int_0^{\pi/2}sin(x)dx \) の計算を示します。

sage入力：
#pre{{
integrate(sin(x), x, 0, pi/2)
}}
&color(blue){1};

もう一つ
$$
\int_0^3 x^2 + 2 sind(x) dx
$$

の例を示し、

sage入力：
#pre{{
integrate(x^2+2*sin(x), x, 0, 3)
}}
&color(blue){-2*cos(3) + 11};

上記結果を数値とすると以下のようになります。
#pre{{
N(_)
}}
&color(blue){12.9799849932009};

*** 数値積分 [#l2f0a04e]
数値積分は、numerical_integral関数を使って計算します。

 numerical_integral(被関分館数, 積分範囲)

で戻り値は、積分結果と誤差を返します。

以下は、\(sin(x\))を0から\(\pi\)で数値積分した結果です。

sage入力：
#pre{{
sigma, error = numerical_integral(sin(x),  0, pi); (sigma, error)
}}
&color(blue){(1.9999999999999998, 2.2204460492503128e-14)};

*** 積分の応用 [#hbbb8e1c]
積分の応用例として、サイクロイドの計算をしてみます。

サイクロイドは、tをパラメータとして、以下のような式で表されます。

$$
\begin{eqnarray}
x &=& 2(t-sin(t)) \\
y &=& 2(1-cos(t))
\end{eqnarray}
$$  

この曲線の長さを計算すると、

$$
\int_{0}^{2\pi}\sqrt{\left(\frac{dx}{dt}\right)^2+\left(\frac{dy}{dt}\right)^2}dt
$$

となります。解析解では、この結果は16となります。

sage入力：
#pre{{
t = var('t')
x = 2*(t-sin(t))
y = 2*(1-cos(t))
parametric_plot([x, y], (t, 0, 2*pi))
}}
&ref(sage0-4.png);

#pre{{
cycloid = sqrt(diff(x,t)^2+diff(y,t)^2)
}}

残念ながら、定積分は結果がでませんでした。

#pre{{
integrate(cycloid, t, 0, 2*pi)
}}
&color(blue){integrate(sqrt(4*(cos(t) - 1)^2 + 4*sin(t)^2), t, 0, 2*pi)};

数値積分で、期待した結果となりました。

#pre{{
numerical_integral(cycloid, 0, 2*pi)
}}
&color(blue){(15.999999999999998, 1.7763568394002502e-13)};


** コメント [#s902e4bb]
#vote(おもしろかった[11],そうでもない[0],わかりずらい[0])
#vote(おもしろかった[11],そうでもない[0],わかりずらい[1])

皆様のご意見、ご希望をお待ちしております。
#comment_kcaptcha
トップ   新規 一覧 単語検索 最終更新   ヘルプ   最終更新のRSS
SmartDoc