SymPy で Matrix() を使ったベクトル・行列の定義と演算，および勾配 (grad)，発散 (div)，回転 (rot) を定義してベクトル解析。

SymPy の import

In [1]:

from sympy.abc import *
from sympy import *

`Matrix()` を使ったベクトルの定義

Matrix() を使った3次元ベクトル $\boldsymbol{a} = (a_1, a_2, a_3)$ の定義例。成分は実数とします。

In [2]:

# 成分に使う変数の宣言
var('a1:4', real=True) # var('a1 a2 a3') と同じ

a = Matrix([a1, a2, a3])
a

Out[2]:

$\displaystyle \left[\begin{matrix}a_{1}\\a_{2}\\a_{3}\end{matrix}\right]$

$\boldsymbol{b} = (b_1, b_2, b_3), \ \boldsymbol{c} = (c_1, c_2, c_3)$ も同様に定義。以下のような書き方もできるようだ。

In [3]:

# Matrix() の中に変数の宣言を入れる方法
b = Matrix(var('b1:4', real=True))
c = Matrix(var('c1:4', real=True))

display(b)
display(c)

$\displaystyle \left[\begin{matrix}b_{1}\\b_{2}\\b_{3}\end{matrix}\right]$

$\displaystyle \left[\begin{matrix}c_{1}\\c_{2}\\c_{3}\end{matrix}\right]$

ベクトルの演算

たし算，引き算 `+, -`

ベクトルの足算引き算をしてみると，成分同士の足算引き算になっていることがわかります。

In [4]:

a + b

Out[4]:

$\displaystyle \left[\begin{matrix}a_{1} + b_{1}\\a_{2} + b_{2}\\a_{3} + b_{3}\end{matrix}\right]$

In [5]:

a - b

Out[5]:

$\displaystyle \left[\begin{matrix}a_{1} -b_{1}\\a_{2} -b_{2}\\a_{3} -b_{3}\end{matrix}\right]$

ベクトルの定数倍（スカラー倍） `*`

ベクトルの定数倍（スカラー倍）は，各成分にそれぞれ定数（スカラー）をかけたものになっています。

In [6]:

var('C')

C*a

Out[6]:

$\displaystyle \left[\begin{matrix}C a_{1}\\C a_{2}\\C a_{3}\end{matrix}\right]$

ベクトルの内積 `.dot()`

ベクトルの内積は $\boldsymbol{a}\cdot\boldsymbol{b} = $a.dot(b) です。

In [7]:

a.dot(b)

Out[7]:

$\displaystyle a_{1} b_{1} + a_{2} b_{2} + a_{3} b_{3}$

ベクトルの大きさ `.norm()`

ベクトル $\boldsymbol{a}$ の大きさ $|\boldsymbol{a}|$ は
$|\boldsymbol{a}| \equiv \sqrt{\boldsymbol{a}\cdot\boldsymbol{a}} =$ a.norm()

ベクトル $\boldsymbol{a}$ の大きさ $|\boldsymbol{a}|$ を計算してみます。

In [8]:

a.norm()

Out[8]:

$\displaystyle \sqrt{a_{1}^{2} + a_{2}^{2} + a_{3}^{2}}$

ベクトルの外積 `.cross()`

ベクトルの外積は $\boldsymbol{a}\times\boldsymbol{b} = $a.cross(b) です。

In [9]:

a.cross(b)

Out[9]:

$\displaystyle \left[\begin{matrix}a_{2} b_{3} -a_{3} b_{2}\\-a_{1} b_{3} + a_{3} b_{1}\\a_{1} b_{2} -a_{2} b_{1}\end{matrix}\right]$

外積は「反交換」

$\boldsymbol{a}\times\boldsymbol{b} {\color{red}{\neq}} \boldsymbol{b}\times\boldsymbol{a}$ であること，具体的には以下のように「反交換」であることを確認します。
$$\boldsymbol{a}\times\boldsymbol{b} = -\boldsymbol{b}\times\boldsymbol{a}$$

$\boldsymbol{a}\times\boldsymbol{b} + \boldsymbol{b}\times\boldsymbol{a} = \boldsymbol{0}$ となることを示せばいいです。

In [10]:

a.cross(b) + b.cross(a)

Out[10]:

$\displaystyle \left[\begin{matrix}0\\0\\0\end{matrix}\right]$

あるいは比較演算子 == を使って a.cross(b) == -b.cross(a) が True であることを示してもよいです。

In [11]:

a.cross(b) == -b.cross(a)

Out[11]:

True

ベクトルの三重積

スカラー三重積

3つのベクトル $\boldsymbol{a}, \boldsymbol{b}, \boldsymbol{c}$ によるスカラー三重積は

$$\boldsymbol{a} \cdot (\boldsymbol{b}\times\boldsymbol{c})$$

であり，これは結局はスカラーとなります。スカラー三重積には，以下のような公式があります。

$$\boldsymbol{a} \cdot (\boldsymbol{b}\times\boldsymbol{c})
= \boldsymbol{b} \cdot (\boldsymbol{c}\times\boldsymbol{a})
= \boldsymbol{c} \cdot (\boldsymbol{a}\times\boldsymbol{b}) $$

○練習：スカラー三重積の公式

$\boldsymbol{a} \cdot (\boldsymbol{b}\times\boldsymbol{c})
= \boldsymbol{b} \cdot (\boldsymbol{c}\times\boldsymbol{a})$ が成り立つことを示せ。

ヒント：

a.dot(b.cross(c)) - b.dot(c.cross(a)) を expand() でかっこをはずして展開してゼロになることを示せばよい。

In [ ]:

同様にして，$\boldsymbol{b} \cdot (\boldsymbol{c}\times\boldsymbol{a})
= \boldsymbol{c} \cdot (\boldsymbol{a}\times\boldsymbol{b}) $ が成り立つことも示せ。

In [ ]:

ベクトル三重積

3つのベクトル $\boldsymbol{a}, \boldsymbol{b}, \boldsymbol{c}$ によるベクトル三重積は

$$\boldsymbol{a} \times (\boldsymbol{b}\times\boldsymbol{c})$$

であり，これは結局はベクトルをつくります。ベクトル三重積には，以下のような公式があります。

$$\boldsymbol{a} \times (\boldsymbol{b}\times\boldsymbol{c})
= (\boldsymbol{a}\cdot\boldsymbol{c})\,\boldsymbol{b} –
(\boldsymbol{a}\cdot\boldsymbol{b})\,\boldsymbol{c} $$

○練習：ベクトル三重積の公式

$\boldsymbol{a} \times (\boldsymbol{b}\times\boldsymbol{c})
= (\boldsymbol{a}\cdot\boldsymbol{c})\,\boldsymbol{b} –
(\boldsymbol{a}\cdot\boldsymbol{b})\,\boldsymbol{c} $ を示せ。

In [ ]:

○練習：ベクトルの内積・外積

$\displaystyle \boldsymbol{a} = \left(-\frac{1}{\sqrt{2}}, \frac{1}{\sqrt{2}}, \sqrt{3} \right), \ \boldsymbol{b} = \left(-1, 1, -\sqrt{2} \right)$ について以下の量を計算しなさい。

$$1. \ |\boldsymbol{a}|, \quad 2. \ |\boldsymbol{b}|, \quad 3. \ \boldsymbol{a}\cdot\boldsymbol{b}, \quad 4. \ \boldsymbol{a}\times\boldsymbol{b}$$

また，2つのベクトル $\boldsymbol{a}, \boldsymbol{b}$ のなす角を $\theta$ としたときの $\cos\theta$ および $\sin\theta$ の値を求めなさい。

まず，ベクトル $\boldsymbol{a}, \boldsymbol{b}$ を Matrix() で定義し…

In [12]:

a = Matrix([-1/sqrt(2), 1/sqrt(2), sqrt(3)])
b = Matrix([-1, 1, -sqrt(2)])

display(a)
display(b)

$\displaystyle \left[\begin{matrix}-\frac{\sqrt{2}}{2}\\\frac{\sqrt{2}}{2}\\\sqrt{3}\end{matrix}\right]$

$\displaystyle \left[\begin{matrix}-1\\1\\-\sqrt{2}\end{matrix}\right]$

$1. \ |\boldsymbol{a}| = $ a.norm() であるから…

In [ ]:

$2. \ |\boldsymbol{b}| = $ b.norm() であるから…

In [ ]:

$3. \ \boldsymbol{a}\cdot\boldsymbol{b} = $ a.dot(b) であるから…

In [ ]:

$4. \ \boldsymbol{a}\times\boldsymbol{b} = $ a.cross(b) であるから…

In [ ]:

$\boldsymbol{a}\cdot\boldsymbol{b} = |\boldsymbol{a}| |\boldsymbol{b}| \cos\theta$ より，$\displaystyle \cos\theta = \frac{\boldsymbol{a}\cdot\boldsymbol{b}}{|\boldsymbol{a}| |\boldsymbol{b}|} = \cdots$

In [13]:

uhen1 = a.dot(b)/(a.norm()*b.norm())
Eq(cos(theta), uhen1)

Out[13]:

$\displaystyle \cos{\left(\theta \right)} = -\frac{\sqrt{6}}{4} + \frac{\sqrt{2}}{4}$

ちなみに $\theta = $ acos(uhen1) であるから，2つのベクトル $\boldsymbol{a}, \boldsymbol{b}$ のなす角 $\theta$ を求めると…

In [14]:

th1 = acos(uhen1)
th1

Out[14]:

$\displaystyle \frac{7 \pi}{12}$

$|\boldsymbol{a}\times\boldsymbol{b}| = |\boldsymbol{a}| |\boldsymbol{b}| \sin\theta$ より $\displaystyle \sin\theta = \frac{|\boldsymbol{a}\times\boldsymbol{b}|}{|\boldsymbol{a}| |\boldsymbol{b}|} = \cdots$

In [15]:

uhen2 = (a.cross(b)).norm()/(a.norm()*b.norm())
Eq(sin(theta), uhen2)

Out[15]:

$\displaystyle \sin{\left(\theta \right)} = \frac{\sqrt{2} \cdot \left(1 + \sqrt{3}\right)}{4}$

$\theta = $ asin(uhen2) であるから，2つのベクトル $\boldsymbol{a}, \boldsymbol{b}$ のなす角 $\theta$ を求めると…

In [16]:

asin(uhen2)

Out[16]:

$\displaystyle \frac{5 \pi}{12}$

acos(uhen1) から求めた角度 $\theta$ と asin(uhen2) から求めた角度 $\theta$ は異なっているが，これはどのように考えたらよいか。

ヒント： $\arccos x$ と $\arcsin x$ の定義域はどちらも $-1 \leq x \leq 1$ だが，値域は？

`Matrix()` を使った行列の定義

Matrix() を使った行列

$$ A = \begin{pmatrix}
a_{11} & a_{12} \\
a_{21} & a_{22}
\end{pmatrix}$$

の定義例。

In [17]:

var('a11:13, a21:23')
A = Matrix([[a11, a12], 
            [a21, a22]])
A

Out[17]:

$\displaystyle \left[\begin{matrix}a_{11} & a_{12}\\a_{21} & a_{22}\end{matrix}\right]$

以下のような書き方もできるようだ。

In [18]:

# Matrix() の中に変数の宣言を入れる方法
B = Matrix([var('b11:13'), 
            var('b21:23')])
B

Out[18]:

$\displaystyle \left[\begin{matrix}b_{11} & b_{12}\\b_{21} & b_{22}\end{matrix}\right]$

行列の演算

たし算，引き算 `+, -`

In [19]:

A + B

Out[19]:

$\displaystyle \left[\begin{matrix}a_{11} + b_{11} & a_{12} + b_{12}\\a_{21} + b_{21} & a_{22} + b_{22}\end{matrix}\right]$

In [20]:

A - B

Out[20]:

$\displaystyle \left[\begin{matrix}a_{11} -b_{11} & a_{12} -b_{12}\\a_{21} -b_{21} & a_{22} -b_{22}\end{matrix}\right]$

行列の積 `*`

In [21]:

A * B

Out[21]:

$\displaystyle \left[\begin{matrix}a_{11} b_{11} + a_{12} b_{21} & a_{11} b_{12} + a_{12} b_{22}\\a_{21} b_{11} + a_{22} b_{21} & a_{21} b_{12} + a_{22} b_{22}\end{matrix}\right]$

In [22]:

B * A

Out[22]:

$\displaystyle \left[\begin{matrix}a_{11} b_{11} + a_{21} b_{12} & a_{12} b_{11} + a_{22} b_{12}\\a_{11} b_{21} + a_{21} b_{22} & a_{12} b_{21} + a_{22} b_{22}\end{matrix}\right]$

行列の積に関しては交換則は成り立ちません。$A B \neq B A$
比較演算子 != を使って「左辺」と「右辺」等しくないことを確認します。

In [23]:

A * B != B * A

Out[23]:

True

転置行列 `.transpose()`

行列 $A$ の転置行列は $A^{T} = $ transpose(A) または A.transpose()

In [24]:

Out[24]:

$\displaystyle \left[\begin{matrix}a_{11} & a_{12}\\a_{21} & a_{22}\end{matrix}\right]$

In [25]:

A.transpose()

Out[25]:

$\displaystyle \left[\begin{matrix}a_{11} & a_{21}\\a_{12} & a_{22}\end{matrix}\right]$

逆行列 `.inv()`

行列 $A$ の逆行列 $\displaystyle A^{-1}$ は A**(-1) または A.inv()

In [26]:

A.inv()

Out[26]:

$\displaystyle \left[\begin{matrix}\frac{a_{22}}{a_{11} a_{22} -a_{12} a_{21}} & -\frac{a_{12}}{a_{11} a_{22} -a_{12} a_{21}}\\-\frac{a_{21}}{a_{11} a_{22} -a_{12} a_{21}} & \frac{a_{11}}{a_{11} a_{22} -a_{12} a_{21}}\end{matrix}\right]$

$A^{-1} A$ や $A A^{-1}$ が単位行列になることを確認します。

In [27]:

A**(-1) * A

Out[27]:

$\displaystyle \left[\begin{matrix}\frac{a_{11} a_{22}}{a_{11} a_{22} -a_{12} a_{21}} -\frac{a_{12} a_{21}}{a_{11} a_{22} -a_{12} a_{21}} & 0\\0 & \frac{a_{11} a_{22}}{a_{11} a_{22} -a_{12} a_{21}} -\frac{a_{12} a_{21}}{a_{11} a_{22} -a_{12} a_{21}}\end{matrix}\right]$

In [28]:

# 上記の結果を簡単化
simplify(_)

Out[28]:

$\displaystyle \left[\begin{matrix}1 & 0\\0 & 1\end{matrix}\right]$

In [29]:

# together で一緒に（通分）
together(A * (A.inv()))

Out[29]:

$\displaystyle \left[\begin{matrix}1 & 0\\0 & 1\end{matrix}\right]$

行列式 `.det()`

$A$ の行列式 $\det A$ は A.det() または det(A)

In [30]:

A.det()

Out[30]:

$\displaystyle a_{11} a_{22} – a_{12} a_{21}$

In [31]:

det(A)

Out[31]:

$\displaystyle a_{11} a_{22} – a_{12} a_{21}$

$\det A = \det A^T$ および $\displaystyle \det A^{-1} = \frac{1}{\det A}$ であることを確認します。

In [32]:

det(A) == det(A.transpose())

Out[32]:

True

In [33]:

det(A.inv()) == 1/det(A)

Out[33]:

True

○練習：行列の積の行列式

$\displaystyle \det\, (A B) = (\det A) \cdot (\det B)$ を確認する。

In [34]:

display(A)
display(B)

$\displaystyle \left[\begin{matrix}a_{11} & a_{12}\\a_{21} & a_{22}\end{matrix}\right]$

$\displaystyle \left[\begin{matrix}b_{11} & b_{12}\\b_{21} & b_{22}\end{matrix}\right]$

In [ ]:

○練習：3次元の回転行列

3次元空間での $z$ 軸周りの角度 $\phi$ の回転を表す回転行列は以下のように書けます。

$$R(\phi) =
\begin{pmatrix}
\cos \phi & -\sin \phi & 0\\
\sin \phi & \cos \phi & 0\\
0 & 0 & 1
\end{pmatrix}
$$

In [35]:

def R(phi):
    return Matrix([[cos(phi), -sin(phi), 0],
                   [sin(phi), cos(phi),  0], 
                   [0,        0,         1]])

$\displaystyle 1.\ R(\phi)$ が直交行列であることを示せ。

直交行列とは，転置行列が逆行列になっている行列のことですから，以下を示せばよいです。

$$R(\phi) \,R^T(\phi) = I$$

ここで，$I$ は単位行列。

In [ ]:

ちなみに，$R(\phi)$ の逆行列は $R(-\phi)$ であることも確認できます。

In [ ]:

位置ベクトル $\boldsymbol{r} = (x, y, z)$ をこの回転行列で直交変換してみます。

In [36]:

# 位置ベクトルの定義

r = Matrix(var('x y z', real = True))
r

Out[36]:

$\displaystyle \left[\begin{matrix}x\\y\\z\end{matrix}\right]$

直交変換されたベクトルの大きさは，元のベクトル $\boldsymbol{r}$ の大きさと同じであることを確かめます。

まず，元のベクトル $\boldsymbol{r}$ の大きさ（の2乗）は…

In [37]:

r.dot(r)

Out[37]:

$\displaystyle x^{2} + y^{2} + z^{2}$

直交変換されたベクトルは…

In [38]:

r1 = R(phi) * r
r1

Out[38]:

$\displaystyle \left[\begin{matrix}x \cos{\left(\phi \right)} – y \sin{\left(\phi \right)}\\x \sin{\left(\phi \right)} + y \cos{\left(\phi \right)}\\z\end{matrix}\right]$

このベクトルの大きさ（の2乗）は…

In [39]:

(r1.dot(r1)).simplify()

Out[39]:

$\displaystyle x^{2} + y^{2} + z^{2}$

ベクトル $\boldsymbol{r}$ の大きさ $|\boldsymbol{r}|$ は r.norm() ですが，なぜこれを使わずに大きさの2乗 $\boldsymbol{r}\cdot\boldsymbol{r}$ を使ったか，わかりますか？

In [40]:

r.norm()

Out[40]:

$\displaystyle \sqrt{x^{2} + y^{2} + z^{2}}$

In [41]:

r1.norm()

Out[41]:

$\displaystyle \sqrt{z^{2} + \left|{x \sin{\left(\phi \right)} + y \cos{\left(\phi \right)}}\right|^{2} + \left|{x \cos{\left(\phi \right)} – y \sin{\left(\phi \right)}}\right|^{2}}$

$\displaystyle 2. \ R(\beta) \, R(\alpha) = R(\alpha+\beta)$ を示せ。

これは，まず $z$ 軸のまわりに角度 $\alpha$ だけ回転させ，引き続き角度 $\beta$ だけ回転させる変換は，$z$ 軸のまわりに一度に角度 $\alpha+\beta$ 回転させる変換と同じだということを表します。

In [42]:

# 変数の宣言
var('alpha beta')

R(beta) * R(alpha)

Out[42]:

$\displaystyle \left[\begin{matrix}- \sin{\left(\alpha \right)} \sin{\left(\beta \right)} + \cos{\left(\alpha \right)} \cos{\left(\beta \right)} & – \sin{\left(\alpha \right)} \cos{\left(\beta \right)} – \sin{\left(\beta \right)} \cos{\left(\alpha \right)} & 0\\\sin{\left(\alpha \right)} \cos{\left(\beta \right)} + \sin{\left(\beta \right)} \cos{\left(\alpha \right)} & – \sin{\left(\alpha \right)} \sin{\left(\beta \right)} + \cos{\left(\alpha \right)} \cos{\left(\beta \right)} & 0\\0 & 0 & 1\end{matrix}\right]$

In [ ]:

ベクトル解析

SymPy にはベクトル解析のためのモジュール sympy.vector があるようですが，簡単なので，以下のようにその都度定義してやってみます。

スカラー場の勾配（`grad()` の定義）

スカラー場 $\varphi(\boldsymbol{r})$ の勾配を与える関数 grad() を定義する。

$$\nabla \varphi \equiv \mbox{grad}\, \varphi = \left(\frac{\partial \varphi}{\partial x}, \frac{\partial \varphi}{\partial y},\frac{\partial \varphi}{\partial z}\right)$$

In [43]:

def grad(phi):
    return Matrix([diff(phi,x), diff(phi,y), diff(phi,z)])

スカラー場 $\varphi(x,y,z)$ の勾配 $\nabla\varphi$ を確認します。

In [44]:

# varphi を x, y, z の関数として宣言
varphi = Function('varphi')(x,y,z)
grad(varphi)

Out[44]:

$\displaystyle \left[\begin{matrix}\frac{\partial}{\partial x} \varphi{\left(x,y,z \right)}\\\frac{\partial}{\partial y} \varphi{\left(x,y,z \right)}\\\frac{\partial}{\partial z} \varphi{\left(x,y,z \right)}\end{matrix}\right]$

勾配の計算例

位置ベクトル $$\boldsymbol{r} = (x, y, z)$$ の大きさ
$$r = |\boldsymbol{r}| = \sqrt{x^2 + y^2 + z^2}$$ について，以下の計算をしてみます。

$$\nabla r = \frac{\partial r}{\partial x} \boldsymbol{i}
+\frac{\partial r}{\partial y} \boldsymbol{j}
+\frac{\partial r}{\partial z} \boldsymbol{k} $$

In [45]:

# 位置ベクトルの定義

rvec = Matrix(var('x y z', real = True))
display(rvec)

r = rvec.norm()
display(r)

$\displaystyle \left[\begin{matrix}x\\y\\z\end{matrix}\right]$

$\displaystyle \sqrt{x^{2} + y^{2} + z^{2}}$

$\nabla r$ を計算します。

In [46]:

grad(r)

Out[46]:

$\displaystyle \left[\begin{matrix}\frac{x}{\sqrt{x^{2} + y^{2} + z^{2}}}\\\frac{y}{\sqrt{x^{2} + y^{2} + z^{2}}}\\\frac{z}{\sqrt{x^{2} + y^{2} + z^{2}}}\end{matrix}\right]$

上の結果が

$$\nabla r =\frac{x}{r} \,\boldsymbol{i} + \frac{y}{r} \,\boldsymbol{j} + \frac{z}{r} \,\boldsymbol{k} = \frac{\boldsymbol{r}}{r}$$

となっていることを確認します。

In [47]:

grad(r) == rvec/r

Out[47]:

True

また，電磁気学では以下のような計算をする必要も出てくるので，やっておきます。

$$\nabla \left(\frac{1}{r}\right) = – \frac{\boldsymbol{r}}{r^3}$$

In [48]:

grad(1/r)

Out[48]:

$\displaystyle \left[\begin{matrix}- \frac{x}{\left(x^{2} + y^{2} + z^{2}\right)^{\frac{3}{2}}}\\- \frac{y}{\left(x^{2} + y^{2} + z^{2}\right)^{\frac{3}{2}}}\\- \frac{z}{\left(x^{2} + y^{2} + z^{2}\right)^{\frac{3}{2}}}\end{matrix}\right]$

In [49]:

grad(1/r) == -rvec/r**3

Out[49]:

True

ベクトル場の発散 (`div()` の定義)

ベクトル場の発散を定義する。

$$\nabla\cdot\boldsymbol{v} \equiv \mbox{div}\, \boldsymbol{v} = \frac{\partial v_x}{\partial x} + \frac{\partial v_y}{\partial y} + \frac{\partial v_z}{\partial z}$$

In [50]:

def div(v):
    return diff(v[0], x) + diff(v[1], y) + diff(v[2], z)

ベクトル場 $\boldsymbol{B}$ の発散

$$\nabla\cdot\boldsymbol{B} = \mbox{div}\, \boldsymbol{B} =
\frac{\partial B_x}{\partial x} + \frac{\partial B_y}{\partial y} +\frac{\partial B_x}{\partial x}$$

を確認します。

In [51]:

# ベクトル場 B の定義
B_x = Function('B_x')(x,y,z)
B_y = Function('B_y')(x,y,z)
B_z = Function('B_z')(x,y,z)

Bvec = Matrix([B_x, B_y, B_z])
Bvec

Out[51]:

$\displaystyle \left[\begin{matrix}B_{x}{\left(x,y,z \right)}\\B_{y}{\left(x,y,z \right)}\\B_{z}{\left(x,y,z \right)}\end{matrix}\right]$

In [52]:

# B の発散
div(Bvec)

Out[52]:

$\displaystyle \frac{\partial}{\partial x} B_{x}{\left(x,y,z \right)} + \frac{\partial}{\partial y} B_{y}{\left(x,y,z \right)} + \frac{\partial}{\partial z} B_{z}{\left(x,y,z \right)}$

発散の計算例

$$\nabla\cdot \boldsymbol{r} = \frac{\partial x}{\partial x} + \frac{\partial y}{\partial y} +\frac{\partial z}{\partial z} = 1+1+1 = 3$$

（答えは空間の次元の数）

In [53]:

# 位置ベクトルが定義されていることを確認
rvec

Out[53]:

$\displaystyle \left[\begin{matrix}x\\y\\z\end{matrix}\right]$

In [54]:

div(rvec)

Out[54]:

$\displaystyle 3$

ベクトル場の回転 (`rot()` の定義)

ベクトルの回転を以下のようにして定義します。

\begin{eqnarray}
\nabla\times\boldsymbol{v} &\equiv& \mbox{rot}\, \boldsymbol{v} \\
&=&\left( \frac{\partial}{\partial y} v_z – \frac{\partial}{\partial z} v_y, \frac{\partial}{\partial z} v_x – \frac{\partial}{\partial x} v_z, \frac{\partial}{\partial x} v_y – \frac{\partial}{\partial y} v_x\right)
\end{eqnarray}

In [55]:

def rot(v):
    vx = v[0]
    vy = v[1]
    vz = v[2]
    return Matrix(
            [diff(vz, y) - diff(vy, z), 
             diff(vx, z) - diff(vz, x), 
             diff(vy, x) - diff(vx, y)]
            )

ベクトル場 $\boldsymbol{A}$ の回転 $\nabla\times\boldsymbol{A} = \mbox{rot}\, \boldsymbol{A} $ を確認します。

In [56]:

# ベクトル場 A の定義
A_x = Function('A_x')(x,y,z)
A_y = Function('A_y')(x,y,z)
A_z = Function('A_z')(x,y,z)

Avec = Matrix([A_x, A_y, A_z])
Avec

Out[56]:

$\displaystyle \left[\begin{matrix}A_{x}{\left(x,y,z \right)}\\A_{y}{\left(x,y,z \right)}\\A_{z}{\left(x,y,z \right)}\end{matrix}\right]$

In [57]:

rot(Avec)

Out[57]:

$\displaystyle \left[\begin{matrix}- \frac{\partial}{\partial z} A_{y}{\left(x,y,z \right)} + \frac{\partial}{\partial y} A_{z}{\left(x,y,z \right)}\\\frac{\partial}{\partial z} A_{x}{\left(x,y,z \right)} – \frac{\partial}{\partial x} A_{z}{\left(x,y,z \right)}\\- \frac{\partial}{\partial y} A_{x}{\left(x,y,z \right)} + \frac{\partial}{\partial x} A_{y}{\left(x,y,z \right)}\end{matrix}\right]$

回転の計算例

$$\nabla\times \left(\frac{\boldsymbol{r}}{r}\right)$$

In [58]:

rvec/r

Out[58]:

$\displaystyle \left[\begin{matrix}\frac{x}{\sqrt{x^{2} + y^{2} + z^{2}}}\\\frac{y}{\sqrt{x^{2} + y^{2} + z^{2}}}\\\frac{z}{\sqrt{x^{2} + y^{2} + z^{2}}}\end{matrix}\right]$

In [59]:

rot(rvec/r)

Out[59]:

$\displaystyle \left[\begin{matrix}0\\0\\0\end{matrix}\right]$

○練習：発散や回転の計算

$\displaystyle 1. \ \nabla\cdot\left(\frac{\boldsymbol{r}}{r}\right)$

In [ ]:

$\displaystyle 2. \ \nabla\cdot\left(\frac{\boldsymbol{r}}{r^n}\right)$

In [ ]:

$\displaystyle 3. \ \nabla\times\left(\frac{\boldsymbol{r}}{r}\right)$

In [ ]:

$\displaystyle 4. \nabla\times\left(\frac{\boldsymbol{r}}{r^n}\right)$

In [ ]:

勾配 (grad) の発散 (div)：ラプラシアンの定義

スカラー場 $\varphi$ の勾配 $\nabla \varphi$ の発散 $$\nabla\cdot(\nabla \varphi) = \nabla^2 \varphi$$
ここで

$$\nabla^2 \equiv \frac{\partial^2}{\partial x^2} + \frac{\partial^2}{\partial y^2} +\frac{\partial^2}{\partial z^2}$$

はラプラス演算子（ラプラシアン）。

In [60]:

varphi = Function('varphi')(x, y, z)
varphi

Out[60]:

$\displaystyle \varphi{\left(x,y,z \right)}$

In [61]:

div(grad(varphi))

Out[61]:

$\displaystyle \frac{\partial^{2}}{\partial x^{2}} \varphi{\left(x,y,z \right)} + \frac{\partial^{2}}{\partial y^{2}} \varphi{\left(x,y,z \right)} + \frac{\partial^{2}}{\partial z^{2}} \varphi{\left(x,y,z \right)}$

スカラーに作用するラプラシアンの定義

In [62]:

# スカラー関数に作用するラプラス演算子の定義
def laplacian(phi):
    # スカラー関数 phi のみに適用
    return div(grad(phi))

ラプラス演算子の計算例

$$\nabla^2 \left(\frac{1}{r}\right) = 0 \quad (r \neq 0)$$

を確認します。

In [63]:

laplacian(1/r).simplify()

Out[63]:

$\displaystyle 0$

勾配 (grad) の回転 (rot)：恒等式

スカラー場 $\varphi$ の勾配 $\nabla \varphi$ の回転 $\nabla\times(\nabla \varphi) $。

恒等的に
$$\nabla\times(\nabla \varphi) = \boldsymbol{0}$$である。

In [64]:

grad(varphi)

Out[64]:

In [65]:

rot(grad(varphi))

Out[65]:

$\displaystyle \left[\begin{matrix}0\\0\\0\end{matrix}\right]$

回転 (rot) の発散 (div)：恒等式

ベクトル場 $\boldsymbol{a}$ の回転 $\nabla \times\boldsymbol{a} $ の発散 $\nabla\cdot(\nabla \times\boldsymbol{a}) $

恒等的に
\begin{eqnarray}
\nabla\cdot(\nabla \times\boldsymbol{a}) &=&
0 \end{eqnarray}
である。

まずはベクトル場 $\boldsymbol{a}$ の定義：

In [66]:

# ベクトル場 a の定義
a_x = Function('a_x')(x,y,z)
a_y = Function('a_y')(x,y,z)
a_z = Function('a_z')(x,y,z)

avec = Matrix([a_x, a_y, a_z])
avec

Out[66]:

$\displaystyle \left[\begin{matrix}a_{x}{\left(x,y,z \right)}\\a_{y}{\left(x,y,z \right)}\\a_{z}{\left(x,y,z \right)}\end{matrix}\right]$

$\nabla\times \boldsymbol{a}$

In [67]:

rot(avec)

Out[67]:

$\displaystyle \left[\begin{matrix}- \frac{\partial}{\partial z} a_{y}{\left(x,y,z \right)} + \frac{\partial}{\partial y} a_{z}{\left(x,y,z \right)}\\\frac{\partial}{\partial z} a_{x}{\left(x,y,z \right)} – \frac{\partial}{\partial x} a_{z}{\left(x,y,z \right)}\\- \frac{\partial}{\partial y} a_{x}{\left(x,y,z \right)} + \frac{\partial}{\partial x} a_{y}{\left(x,y,z \right)}\end{matrix}\right]$

$\nabla\cdot\left(\nabla\times \boldsymbol{a}\right)$

In [68]:

div(rot(avec))

Out[68]:

$\displaystyle 0$

回転 (rot) の回転 (rot)

ベクトル場 $\boldsymbol{a}$ の回転 $\nabla \times\boldsymbol{a} $ の回転 $\nabla\times(\nabla \times\boldsymbol{a}) $

$$ \nabla\times(\nabla \times\boldsymbol{a}) = \nabla\left( \nabla\cdot\boldsymbol{a} \right) – \nabla^2 \boldsymbol{a} $$

が成り立つことを確認する。

ベクトルに作用するラプラシアンの定義

In [69]:

# ベクトルに作用する Laplacian を定義
def Laplacian(v):
    return Matrix([laplacian(v[0]),
                   laplacian(v[1]),
                   laplacian(v[2])])

In [70]:

rot(rot(avec))

Out[70]:

$\displaystyle \left[\begin{matrix}- \frac{\partial^{2}}{\partial y^{2}} a_{x}{\left(x,y,z \right)} – \frac{\partial^{2}}{\partial z^{2}} a_{x}{\left(x,y,z \right)} + \frac{\partial^{2}}{\partial y\partial x} a_{y}{\left(x,y,z \right)} + \frac{\partial^{2}}{\partial z\partial x} a_{z}{\left(x,y,z \right)}\\- \frac{\partial^{2}}{\partial x^{2}} a_{y}{\left(x,y,z \right)} – \frac{\partial^{2}}{\partial z^{2}} a_{y}{\left(x,y,z \right)} + \frac{\partial^{2}}{\partial y\partial x} a_{x}{\left(x,y,z \right)} + \frac{\partial^{2}}{\partial z\partial y} a_{z}{\left(x,y,z \right)}\\- \frac{\partial^{2}}{\partial x^{2}} a_{z}{\left(x,y,z \right)} – \frac{\partial^{2}}{\partial y^{2}} a_{z}{\left(x,y,z \right)} + \frac{\partial^{2}}{\partial z\partial x} a_{x}{\left(x,y,z \right)} + \frac{\partial^{2}}{\partial z\partial y} a_{y}{\left(x,y,z \right)}\end{matrix}\right]$

In [71]:

grad(div(avec)) - Laplacian(avec)

Out[71]:

In [72]:

rot(rot(avec)) == grad(div(avec)) - Laplacian(avec)

Out[72]:

True

○練習：ベクトルの外積の発散

以下の式が成り立つことを示せ。

$$\nabla\cdot\left( \boldsymbol{E}\times\boldsymbol{B}\right)
= \left(\nabla\times\boldsymbol{E} \right)\cdot\boldsymbol{B}
– \boldsymbol{E}\cdot\left(\nabla\times\boldsymbol{B} \right)$$

In [73]:

E_x = Function('E_x')(x,y,z)
E_y = Function('E_y')(x,y,z)
E_z = Function('E_z')(x,y,z)
B_x = Function('B_x')(x,y,z)
B_y = Function('B_y')(x,y,z)
B_z = Function('B_z')(x,y,z)

E = Matrix([E_x, E_y, E_z])
B = Matrix([B_x, B_y, B_z])

display(E)
display(B)

$\displaystyle \left[\begin{matrix}E_{x}{\left(x,y,z \right)}\\E_{y}{\left(x,y,z \right)}\\E_{z}{\left(x,y,z \right)}\end{matrix}\right]$

$\displaystyle \left[\begin{matrix}B_{x}{\left(x,y,z \right)}\\B_{y}{\left(x,y,z \right)}\\B_{z}{\left(x,y,z \right)}\end{matrix}\right]$

In [ ]:

○練習：スカラーとベクトルの積の発散

以下の式が成り立つことを示せ。

$$\nabla\cdot\left( \varphi\,\boldsymbol{D}\right) =
\left(\nabla\varphi\right)\cdot\boldsymbol{D} +\varphi\left(\nabla\cdot\boldsymbol{D}\right)$$

In [74]:

varphi = Function('varphi')(x, y, z)
display(varphi)

D_x = Function('D_x')(x,y,z)
D_y = Function('D_y')(x,y,z)
D_z = Function('D_z')(x,y,z)

D = Matrix([D_x, D_y, D_z])
display(D)

$\displaystyle \varphi{\left(x,y,z \right)}$

$\displaystyle \left[\begin{matrix}D_{x}{\left(x,y,z \right)}\\D_{y}{\left(x,y,z \right)}\\D_{z}{\left(x,y,z \right)}\end{matrix}\right]$

In [ ]:

5. SymPy でベクトル・行列・ベクトル解析

SymPy の import

Matrix() を使ったベクトルの定義

ベクトルの演算

たし算，引き算 +, -

ベクトルの定数倍（スカラー倍） *

ベクトルの内積 .dot()

ベクトルの大きさ .norm()

ベクトルの外積 .cross()

外積は「反交換」

ベクトルの三重積

スカラー三重積

○練習：スカラー三重積の公式

ベクトル三重積

○練習：ベクトル三重積の公式

○練習：ベクトルの内積・外積

Matrix() を使った行列の定義

行列の演算

たし算，引き算 +, -

行列の積 *

転置行列 .transpose()

逆行列 .inv()

行列式 .det()

○練習：行列の積の行列式

○練習：3次元の回転行列

ベクトル解析

スカラー場の勾配（grad() の定義）

勾配の計算例

ベクトル場の発散 (div() の定義)

発散 の計算例

ベクトル場の回転 (rot() の定義)

回転の計算例

○練習：発散や回転の計算

勾配 (grad) の発散 (div)：ラプラシアンの定義

スカラーに作用するラプラシアンの定義

ラプラス演算子の計算例

勾配 (grad) の回転 (rot)：恒等式

回転 (rot) の発散 (div)：恒等式

回転 (rot) の回転 (rot)

ベクトルに作用するラプラシアンの定義

○練習：ベクトルの外積の発散

○練習：スカラーとベクトルの積の発散

このページの目次

このセクションの構成

`Matrix()` を使ったベクトルの定義

たし算，引き算 `+, -`

ベクトルの定数倍（スカラー倍） `*`

ベクトルの内積 `.dot()`

ベクトルの大きさ `.norm()`

ベクトルの外積 `.cross()`

`Matrix()` を使った行列の定義

たし算，引き算 `+, -`

行列の積 `*`

転置行列 `.transpose()`

逆行列 `.inv()`

行列式 `.det()`

スカラー場の勾配（`grad()` の定義）

ベクトル場の発散 (`div()` の定義)

発散の計算例

ベクトル場の回転 (`rot()` の定義)