今回はfortranを使って行列の積を求めていきたいと思います。
紙の上では計算できても、いざコンピューター上で計算させようとすると非常に混乱したので参考にしてもらえたらと思います。
行列の積
行列の書き方
行列は以下のようにして書くことができます。
real(8) :: amat(3,3) data amat /2.,1.,0., 1.,2.,1., 0.,1.,2./
この場合、amatという行列を(3,3)で定義したあと(一行目)、実際に値を入れています(二行目)。
数式で書くと
$$
amat = \begin{pmatrix}2&1&0\\1&2&1\\0&1&2\end{pmatrix}
$$
という式を作ったことになります。
行列の表示方法
行列は表示するにもひと工夫必要です。
do i=1,3 write(*,'(3F15.7)')amat(i,1),amat(i,2),amat(i,3) end do
と繰り返し処理を行って、列ごとに数字を入れていきます。
行列の積のコード
では行列の定義・表示方法が分かったところで、行列の積を求めるコードを書いていきましょう。
このような形になります。
program matrix implicit none real(8) :: amat(3,3),bmat(3,3),cmat(3,3) data amat /2.,1.,0., 1.,2.,1., 0.,1.,2./ data bmat /1.,1.,3., 2.,3.,4., 7.,2.,5./ integer :: i,j,k do i=1,3 do j=1,3 do k=1,3 cmat(i,j) = cmat(i,j)+amat(i,k)*bmat(k,j) end do end do end do write(*,*) "cmat=" do i=1,3 write(*,'(3F15.7)')cmat(i,1),cmat(i,2),cmat(i,3) end do stop end program matrix
繰り返し処理部分
これだけだとよくわからないと思うので、解説していきたいと思います。
まず、
do i=1,3 do j=1,3 do k=1,3 cmat(i,j) = cmat(i,j)+amat(i,k)*bmat(k,j) end do end do end do
で、amatとbmatの掛け算をしています。
何をしているのか数式で見てみましょう。
まず、次のような行列を考えます。
$$
A = \begin{pmatrix}a_{11}&a_{12}&a_{13}\\ a_{21}&a_{22}&a_{23} \\ a_{31}&a_{32}&a_{33} \end{pmatrix}
$$
$$
B = \begin{pmatrix}b_{11}&b_{12}&b_{13}\\ b_{21}&b_{22}&b_{23} \\ b_{31}&b_{32}&b_{33} \end{pmatrix}
$$
ここでAとBの積Cを考えます。つまり、
$$
C = A\times B$$
全部書くと大変なのでCの一行一列目つまり、(11)成分だけを考えると、
$$
c_{11} = a_{11}b_{11} + a_{12}b_{21} + a_{13}b_{31}
$$
となります。
ここで、内側の数字と外側の数字に注目します。 (\(a_{12}b_{21}だったら、22が内側、11が外側の数字\))
すると
- 内側の数字は常に同じで、同成分内で1~3まで変化
- 外側の数字はその成分がどこに位置しているかを表す
ということが分かります。
なので、これをfortranで書くと
do k=1,3 cmat(i,j) = cmat(i,j)+amat(i,k)*bmat(k,j) end do
となります。同じ成分の中(i,j)で内側の数字kが1~3まで変化していて、それらを足しています。
さらにこれをすべての行と列で行っていくので、
do i=1,3 do j=1,3 do k=1,3 cmat(i,j) = cmat(i,j)+amat(i,k)*bmat(k,j) end do end do end do
となります。
matmul()
このように大変な思いをして積を計算してきたわけですが、もっと簡単に計算することができます。
それが組み込み関数matmul()です。
使い方としては
matmul(amat,bmat)
とするだけでamatとbmatの積が計算できます。
非常に簡単ですね!!