.NET5 のサンプル： MatrixClass.dll

MatrixClass.dll には Matrix と Vector クラスが含まれています。

Vector クラスは実数のベクトル、Matrix クラスは実数の行列をモデル化したクラスです。

この DLL のソリューションエクスプローラは次の画像のようになっています。

クラスの定義

名前空間

MatrixClass

クラス

• Vector　実数のベクトル
• Matrix　実数の行列

Vector クラス

コンストラクタ

• Vector(int)　パラメータで指定した長さの Vector オブジェクトを構築する。
• Vector(double[])　配列を元に Vector オブジェクトを構築する。

プロパティ

• Length　ベクトルの要素数 (読みのみ)
• this[int]　指定した位置の要素の値 (読み書き可能)

メソッド

• Zero()　要素をすべてゼロクリアする。
• Inv()　要素の符号をすべて反転する。
• ToArray()　配列に変換する（メソッドの戻り値として返す)
• ToString()　文字列表現を返す。
• Clone()　クローンを作成してその参照を返す。
• Equals(object)　オーバーライドメソッド。パラメータと比較して等しければ true、そうでなければ false を返す。
• Dot(Vector, Vector)　スタティックメソッド。ベクトルの内積の値を返す。

演算子

• operator + (Vector, double[])　Vector と実数配列の和を求める。結果は実数配列。
• operator + (Vector, Vector)　Vector と Vector の和を求める。結果は Vector。
• operator – (Vector, double[])　Vector と実数配列の差を求める。結果は実数配列。
• operator – (Vector, Vector)　Vector と Vector の差を求める。結果は Vector。
• operator + (double, Vector)　スカラーと Vector の掛け算を行う。結果は Vector。

Matrix クラス

コンストラクタ

• Matrix(int, int)　行と列を指定して Matrix オブジェクトを構築する。
• Matrix(double[ , ])　2次元配列を元にして Matrix オブジェクトを構築する。
• Matrix(double[][])　2次元ジャグ配列 (配列の配列) を元に Matrix オブジェクトを構築する。

プロパティ

• Size　int[] で行数、列数を返す。
• this[int, int]　行、列で指定した要素の値。(読み書き可能)

メソッド

• Get(int, int)　行、列で指定した要素の値を返す。
• Set(int, int, double)　　行、列で指定した要素の値を変更する。
• AppendRow(double[])　行を追加する。追加後の行数を返す。
• GetRow(int)　指定した行を double[] で返す。
• GetColumn(int)　指定した列を double[] で返す。
• Zero()　行列の要素をすべてゼロクリアする。
• Inv()　行列の要素の符号をすべて反転する。
• Clone()　この行列オブジェクトのクローンを作る。(新しい Matrix オブジェクトを返す)
• ToString()　オーバーライドメソッド。文字列表現を返す。
• Equals(object)　Matrix オブジェクトの要素がすべて等しいか判別する。等しい場合は true を返す。

演算子

• operator + (Matrix, Matrix)　Matrix オブジェクトどうしの加算。結果は Matrix オブジェクト。
• operator – (Matrix, Matrix)　Matrix オブジェクトどうしの減算。結果は Matrix オブジェクト。
• operator * (Matrix, Matrix)　Matrix オブジェクトどうしの乗算。結果は Matrix オブジェクト。
• operator * (double, Matrix)　スカラーとMatrix オブジェクトの乗算。結果は Matrix オブジェクト。
• operator * (Matrix, Vector)　Matrix オブジェクトと Vector オブジェクトの乗算。結果は Vector オブジェクト。

ソースプログラム

using System;
using System.Text;
using System.Collections.Generic;

namespace MatrixClass
{
    /// <summary>
    /// ベクトルクラス (Vector)
    /// </summary>
    public class Vector
    {
        protected List<double> vector;

        /// <summary>
        /// 空のベクトルを作る。
        /// </summary>
        /// <param name="n"></param>
        public Vector(int n)
        {
            vector = new List<double>();
            for (int i = 0; i < n; i++)
            {
                vector.Add(0.0);
            }
        }

        /// <summary>
        /// 配列を元にベクトルを作る。
        /// </summary>
        /// <param name="a">配列</param>
        public Vector(double[] a)
        {
            vector = new List<double>(a);
        }

        /// <summary>
        /// 要素の数
        /// </summary>
        public int Length
        {
            get
            {
                return this.vector.Count;
            }
        }

        /// <summary>
        /// ゼロクリアする。
        /// </summary>
        public void Zero()
        {
            vector.ForEach((x) => { x = 0.0; });
        }

        /// <summary>
        /// 符号を反転する。
        /// </summary>
        public void Inv()
        {
            vector.ForEach((x) => { x = -x; });
        }

        /// <summary>
        /// 配列に変換する。
        /// </summary>
        /// <returns></returns>
        public double[] ToArray()
        {
            return vector.ToArray();
        }

        /// <summary>
        /// 文字列表現を返す。
        /// </summary>
        /// <returns></returns>
        public override string ToString()
        {
            var sb = new StringBuilder();
            foreach (double x in this.vector)
            {
                sb.Append(x);
                sb.Append(" ");
            }
            return sb.ToString().Trim();
        }

        /// <summary>
        /// クローンを作る。
        /// </summary>
        /// <returns>Vector のクローン</returns>
        public Vector Clone()
        {
            var v = new Vector(this.Length);
            for (int i = 0; i < this.Length; i++)
            {
                v[i] = this.vector[i];
            }
            return v;
        }

        /// <summary>
        /// ベクトルの要素
        /// </summary>
        /// <param name="n">要素のインデックス</param>
        /// <returns>要素の値</returns>
        public double this[int n]
        {
            get
            {
                return this.vector[n];
            }

            set
            {
                this.vector[n] = value;
            }
        }

        /// <summary>
        /// ベクトルの加算 (Vector+Array)
        /// </summary>
        /// <param name="v">Vector オブジェクト</param>
        /// <param name="a">配列</param>
        /// <returns></returns>
        public static double[] operator + (Vector v, double[] a)
        {
            var v2 = v.Clone();
            for (int i = 0; i < a.Length; i++)
            {
                v2[i] = v[i] + a[i];
            }
            return v2.ToArray();
        }

        /// <summary>
        /// ベクトルの加算 (Vector+Vector)
        /// </summary>
        /// <param name="v">Vector オブジェクト</param>
        /// <param name="a">Vector オブジェクト</param>
        /// <returns></returns>
        public static Vector operator + (Vector v, Vector a)
        {
            var v2 = v.Clone();
            for (int i = 0; i < a.Length; i++)
            {
                v2[i] = v[i] + a[i];
            }
            return v2;
        }

        /// <summary>
        /// ベクトルの減算 (Vector-Array)
        /// </summary>
        /// <param name="v">Vector オブジェクト</param>
        /// <param name="a">配列</param>
        /// <returns></returns>
        public static double[] operator -(Vector v, double[] a)
        {
            var v2 = v.Clone();
            for (int i = 0; i < a.Length; i++)
            {
                v2[i] = v[i] - a[i];
            }
            return v2.ToArray();
        }

        /// <summary>
        /// ベクトルの減算 (Vector-Vector)
        /// </summary>
        /// <param name="v">Vector オブジェクト</param>
        /// <param name="a">Vector オブジェクト</param>
        /// <returns></returns>
        public static Vector operator -(Vector v, Vector a)
        {
            var v2 = v.Clone();
            for (int i = 0; i < a.Length; i++)
            {
                v2[i] = v[i] - a[i];
            }
            return v2;
        }

        /// <summary>
        /// スカラーとベクトルの積
        /// </summary>
        /// <param name="r">スカラー</param>
        /// <param name="v">ベクトル</param>
        /// <returns></returns>
        public static Vector operator * (double r, Vector v)
        {
            var v2 = v.Clone();
            for (int i = 0; i < v.Length; i++)
            {
                v2[i] *= r;
            }
            return v2;
        }

        /// <summary>
        /// ベクトルの内積
        /// </summary>
        /// <param name="u">Vector オブジェクト</param>
        /// <param name="v">Vector オブジェクト</param>
        /// <returns>内積</returns>
        public static double Dot(Vector u, Vector v)
        {
            double y = 0.0;
            for (int i = 0; i < u.Length; i++)
            {
                y += u[i] * v[i];
            }
            return y;
        }

        /// <summary>
        /// Vector どうしの比較
        /// </summary>
        /// <param name="obj">比較対象のベクトル</param>
        /// <returns>等しいければ true</returns>
        public override bool Equals(object obj)
        {
            var vec = obj as Vector;
            if (vec.Length != this.Length)
            {
                return false;
            }
            for (int i = 0; i < this.Length; i++)
            {
                if (this[i] != vec[i])
                {
                    return false;
                }
            }
            return true;
        }

        /// <summary>
        /// ハッシュコードを返す。
        /// </summary>
        /// <returns></returns>
        public override int GetHashCode()
        {
            return base.GetHashCode();
        }

    }

    /// <summary>
    /// 行列クラス (Matrix)
    /// </summary>
    public class Matrix
    {
        /// <summary>
        /// 行列オブジェクト
        /// </summary>
        protected List<double[]> matrix;

        /// <summary>
        /// 空の行列を作る。
        /// </summary>
        /// <param name="m">行数</param>
        /// <param name="n">列数</param>
        public Matrix(int m, int n)
        {
            matrix = new List<double[]>();
            for (int i = 0; i < m; i++)
            {
                matrix.Add(new double[n]);
            }
        }

        /// <summary>
        /// 2次元配列から行列を作る。
        /// </summary>
        /// <param name="a">元となる配列</param>
        public Matrix(double[,] a)
        {
            matrix = new List<double[]>();
            int m = a.GetLength(0);
            int n = a.GetLength(1);
            for (int i = 0; i < m; i++)
            {
                var v = new double[n];
                for (int j = 0; j < n; j++)
                {
                    v[j] = a[i, j];
                }
                matrix.Add(v);
            }
        }

        /// <summary>
        /// ジャグ配列で初期化
        /// </summary>
        /// <param name="a"></param>
        public Matrix(double[][] a)
        {
            matrix = new List<double[]>();
            int m = a.Length;
            for (int i = 0; i < m; i++)
            {
                matrix.Add(a[i]);
            }
        }

        /// <summary>
        /// 行列のサイズ [行数、列数]
        /// </summary>
        public int[] Size
        {
            get
            {
                var size = new int[2];
                size[0] = matrix.Count;
                size[1] = matrix[0].Length;
                return size;
            }
        }

        /// <summary>
        /// 行列の要素を得る。
        /// </summary>
        /// <param name="m">行の指定</param>
        /// <param name="n">列の指定</param>
        /// <returns>要素の値</returns>
        public double Get(int m, int n)
        {
            double[] row = matrix[m];
            return row[n];
        }

        /// <summary>
        /// 行列の要素の値を変更する。
        /// </summary>
        /// <param name="m">行の指定</param>
        /// <param name="n">列の指定</param>
        /// <param name="v">要素の値</param>
        public void Set(int m, int n, double v)
        {
            double[] row = matrix[m];
            row[n] = v;
        }

        /// <summary>
        /// 行列の要素
        /// </summary>
        /// <param name="m"></param>
        /// <param name="n"></param>
        /// <returns></returns>
        public double this[int m, int n]
        {
            get
            {
                return Get(m, n);
            }

            set
            {
                Set(m, n, value);
            }
        }

        /// <summary>
        /// 配列を行として追加する。
        /// </summary>
        /// <param name="a">追加する配列</param>
        /// <returns>追加後の行数</returns>
        public int AppendRow(double[] a)
        {
            if (matrix.Count > 0 && matrix[0].Length != a.Length)
                throw new Exception("行列の列数と入力ベクトルの整合性がありません。");

            matrix.Add(a);
            return matrix.Count;
        }

        /// <summary>
        /// 指定した行を配列として取得する。
        /// </summary>
        /// <param name="n"></param>
        /// <returns></returns>
        public double[] GetRow(int n)
        {
            int[] size = this.Size;
            double[] row = new double[size[1]];
            for (int i = 0; i < size[1]; i++)
            {
                row[i] = this[n, i];
            }
            return row;
        }

        /// <summary>
        /// 指定した列を配列として取得する。
        /// </summary>
        /// <param name="m"></param>
        /// <returns></returns>
        public double[] GetColumn(int m)
        {
            int[] size = this.Size;
            double[] col = new double[size[0]];
            for (int i = 0; i < size[0]; i++)
            {
                col[i] = this[i, m];
            }
            return col;
        }

        /// <summary>
        /// 行列の要素をゼロクリアする。
        /// </summary>
        public void Zero()
        {
            for (int i = 0; i < this.Size[0]; i++)
            {
                double[] row = matrix[i];

                for (int j = 0; j < this.Size[1]; j++)
                {
                    row[j] = 0.0;
                }
            }
        }

        /// <summary>
        /// 行列の要素の符号を反転させる。
        /// </summary>
        public void Inv()
        {
            for (int i = 0; i < this.Size[0]; i++)
            {
                double[] row = matrix[i];

                for (int j = 0; j < this.Size[1]; j++)
                {
                    row[j] = -row[j];
                }
            }
        }

        /// <summary>
        /// 行列のクローンを作る。
        /// </summary>
        /// <param name="a">行列</param>
        /// <returns>クローン</returns>
        public Matrix Clone()
        {
            int m = this.Size[0];
            int n = this.Size[1];
            Matrix mx = new Matrix(m, n);
            for (int i = 0; i < m; i++)
            {
                for (int j = 0; j < n; j++)
                {
                    mx[i, j] = this[i, j];
                }
            }
            return mx;
        }

        /// <summary>
        /// 文字列表現を返す。
        /// </summary>
        /// <returns>文字列</returns>
        public override string ToString()
        {
            var sb = new StringBuilder();
            int m = this.Size[0];
            int n = this.Size[1];

            for (int i = 0; i < m; i++)
            {
                for (int j = 0; j < n; j++)
                {
                    sb.Append(Get(i, j).ToString());
                    sb.Append(" ");
                }
                sb.Append("\n");
            }
            return sb.ToString();
        }

        /// <summary>
        /// 行列の加算
        /// </summary>
        /// <param name="a">行列</param>
        /// <param name="b">行列</param>
        /// <returns>行列</returns>
        public static Matrix operator + (Matrix a, Matrix b)
        {
            var size = a.Size;
            Matrix m = a.Clone();
            for (int i = 0; i < size[0]; i++)
            {
                for (int j = 0; j < size[1]; j++)
                {
                    m[i, j] += b[i, j];
                }
            }
            return m;
        }

        /// <summary>
        /// 行列の減算
        /// </summary>
        /// <param name="a">行列</param>
        /// <param name="b">行列</param>
        /// <returns>行列</returns>
        public static Matrix operator -(Matrix a, Matrix b)
        {
            var size = a.Size;
            Matrix m = a.Clone();
            for (int i = 0; i < size[0]; i++)
            {
                for (int j = 0; j < size[1]; j++)
                {
                    m[i, j] -= b[i, j];
                }
            }
            return m;
        }

        /// <summary>
        /// 行列の掛け算
        /// </summary>
        /// <param name="a">行列</param>
        /// <param name="b">行列</param>
        /// <returns>行列</returns>
        public static Matrix operator *(Matrix a, Matrix b)
        {
            var size = a.Size;
            double[,] y = new double[size[0], size[1]];
            for (int i = 0; i < size[0]; i++)
            {
                for (int j = 0; j < size[1]; j++)
                {
                    y[i, j] = 0.0;
                }
            }
            for (int i = 0; i < size[0]; i++)
            {
                for (int j = 0; j < size[1]; j++)
                {
                    double[] row = a.GetRow(i);
                    double[] col = b.GetColumn(j);
                    y[i, j] = InnerProduct(row, col);
                }
            }
            return new Matrix(y);
        }

        public static double InnerProduct(double[] a, double[] b)
        {
            double s = 0.0;
            for (int i = 0; i < a.Length; i++)
            {
                s += a[i] * b[i];
            }
            return s;
        }

        /// <summary>
        /// スカラーと行列の掛け算
        /// </summary>
        /// <param name="r">スカラー</param>
        /// <param name="a">行列</param>
        /// <returns>行列</returns>
        public static Matrix operator * (double r, Matrix a)
        {
            var size = a.Size;
            Matrix m = a.Clone();
            for (int i = 0; i < size[0]; i++)
            {
                for (int j = 0; j < size[1]; j++)
                {
                    m[i, j] = r * a[i, j];
                }
            }
            return m;
        }

        /// <summary>
        /// 行列とベクトルの掛け算
        /// </summary>
        /// <param name="a">行列</param>
        /// <param name="v">ベクトル</param>
        /// <returns>ベクトル</returns>
        public static Vector operator * (Matrix a, Vector v)
        {
            var u = new Vector(v.Length);
            int m = a.Size[0];
            int n = a.Size[1];

            for (int i = 0; i < m; i++)
            {
                u[i] = 0.0;
                for (int j = 0; j < n; j++)
                {
                    u[i] += a[i, j] * v[j];
                }
            }
            return u;
        }

        /// <summary>
        /// 行列の比較
        /// </summary>
        /// <param name="obj">比較対象の行列</param>
        /// <returns>等しければ true</returns>
        public override bool Equals(object obj)
        {
            var mat = obj as Matrix;
            var size = this.Size;
            if (size[0] == mat.Size[0] && size[1] == mat.Size[1])
            {
                for (int i = 0; i < size[0]; i++)
                {
                    for (int j = 0; j < size[1]; j++)
                    {
                        if (this[i, j] != mat[i, j])
                        {
                            return false;
                        }
                    }
                }
            }
            else
            {
                return false;
            }
            return true;
        }

        /// <summary>
        /// ハッシュコードを返す。
        /// </summary>
        /// <returns></returns>
        public override int GetHashCode()
        {
            return base.GetHashCode();
        }
    }
}

使用例

using System;
using MatrixClass;

namespace TestMatrixClass
{
    /// <summary>
    /// MatrixClass のテスト
    /// </summary>
    class Program
    {
        static void Main(string[] args)
        {
            Console.WriteLine("MatrixClass のテスト");
            if (args.Length == 0)
            {
                // Matrix のテスト (1)
                var mat1 = new MatrixClass.Matrix(3, 3);
                mat1.Zero();
                mat1.Set(0, 0, 1.0);
                Console.WriteLine(mat1.Get(0, 0));
                double[,] a = new double[,] { { 1.0, 2.0, 3.0 }, { 4.0, 5.0, 6.0 }, { 7.0, 8.0, 9.0 } };
                var mat2 = new MatrixClass.Matrix(a);
                Console.WriteLine(mat2.Get(1, 1));
                double[][] b = new double[3][];
                b[0] = new double[3] { 1.0, 2.0, 3.0 };
                b[1] = new double[3] { 4.0, 5.0, 6.0 };
                b[2] = new double[3] { 7.0, 8.0, 9.0 };
                var mat3 = new MatrixClass.Matrix(b);
                Console.WriteLine(mat3.Get(2, 2));
                Console.WriteLine(mat3[2, 2]);
                Console.WriteLine(mat3.ToString());
            }
            else
            {
                switch (args[0])
                {
                    case "-v":
                        {
                            // Vector のテスト
                            var v1 = new Vector(3);
                            v1.Zero();
                            var v2 = new Vector(new double[] { 1.0, 2.0, 3.0 });
                            Console.WriteLine(v1.ToString());
                            Console.WriteLine(v2.ToString());
                            v1[0] = -1.0;
                            v1[1] = -2.0;
                            v1[2] = -3.0;
                            Vector v = v1 + v2;
                            Console.WriteLine(v.ToString());
                            Vector u = v1 - v2;
                            Console.WriteLine(u.ToString());
                            Vector w = 3.0 * v1;
                            Console.WriteLine(w);
                            double x = Vector.Dot(v1, v2);
                            Console.WriteLine(x);
                        }
                        break;

                    case "-m":
                        {
                            // Matrix 演算のテスト
                            var m1 = new Matrix(2, 2);
                            m1.Zero();
                            m1[0, 0] = 1.0;
                            m1[1, 1] = 1.0;
                            Console.WriteLine("m1 = \n" + m1.ToString());
                            var m2 = new Matrix(new double[,] { { 1.0, -1.0 }, { 0.0, 1.0 } });
                            Console.WriteLine("m2 = \n" +  m2.ToString());
                            var m3 = new Matrix(new double[][] { new double[] {2.0, 1.0 }, new double[] {-1.0, -3.0 } });
                            Console.WriteLine("m3 = \n" + m3.ToString());
                            var m = m1 + m2;
                            Console.WriteLine("m1 + m2 = \n" + m.ToString());
                            m = m3 - m2;
                            Console.WriteLine("m3 - m2 = \n" + m.ToString());

                        }
                        break;

                    default:
                        break;
                }
            }
        }
    }
}

PowerShell で MatrixClass を使う例

Matrix クラスの使用例

using namespace System.Numerics
if ($PSVersionTable.Platform -eq "Unix") {
    Add-Type -Path "/home/user/lib/NET5/MatrixClass.dll"
}
else {
    Add-Type -Path "d:/lib/NET5/MatrixClass.dll"
}

$Matrix = [MatrixClass.Matrix]
$Vector = [MatrixClass.Vector]

$m0 = $Matrix::new(2, 2)
$m0.Zero()
$m0.Set(1, 1, 1)
$m0[0, 1] = 2
$m0.Get(1, 1)
$m0[0, 1]
$m0.ToString()
$a = ((1, 0), (0, 1))
$m1 = $Matrix::new($a)
$m1.ToString()
$b = (1.5, -2.0), (0, 3.0)
$m2 = $Matrix::new($b)
$m2.ToString()
$mx = $m1 + $m2
$mx.ToString()
$mx = $m1 - $m2
$mx.ToString()
$mx = $m1 * $m2
$mx.ToString()
$v = $Vector::new((2, 2))
$mx = $m1 * $v
$mx.ToString()

実行例

PS D:\workspace\PowerShell\NET5DLL> .\Matrix.ps1
1
2
0 2
0 1

1 0
0 1

1.5 -2
0 3

2.5 -2
0 4

-0.5 2
0 -2

1.5 -2
0 3

2 2

Vector クラスの使用例

using namespace System.Numerics
if ($PSVersionTable.Platform -eq "Unix") {
    Add-Type -Path "/home/user/lib/NET5/MatrixClass.dll"
}
else {
    Add-Type -Path "D:/lib/NET5/MatrixClass.dll"
}

$Vector = [MatrixClass.Vector]
$v0 = $Vector::new(3);
$v1 = $Vector::new((-1, 0, 1))

$v0.Zero()
$v0.ToString()
$v1.ToString()
$v0[1] = 2.0
$v0.ToString()
$v2 = $v0 + $v1
$V2.ToString()
$a = $v2.ToArray()
$a.GetType()
$a
$Vector::dot($v1, $v2)

実行例

PS D:\workspace\PowerShell\NET5DLL> .\Vector.ps1
0 0 0
-1 0 1
0 2 0
-1 2 1

IsPublic IsSerial Name                                     BaseType
-------- -------- ----                                     --------
True     True     Double[]                                 System.Array
-1
2
1
2

 

Julia の勉強：線形代数 (いろいろな行列)

単位行列 (Identity Matrix) など

単位行列は対角成分が 1 で、それ以外はすべて 0 の行列です。

成分がすべてゼロ (zeros()関数を利用) やすべて任意の値 (fill()関数を利用) である行列も簡単に作れます。

乱数行列

Julia では単なる乱数だけでなく、行列の成分が乱数である乱数行列を rand() 関数で簡単に作れます。

ベクトルの回転行列

ベクトルに行列をかけると別のベクトルに変換できます。その中で幾何学的に意味を持たせたベクトルを回転するための回転ベクトルがあります。

転置行列 (Transpose Matrix)

正方行列の対角線で要素の上下を入れ替えた行列のことを転置行列と言います。

Julia では LinearAlgebra.Transpose () 関数で求めることができます。

逆行列 (Inverse Matrix)

元の行列とその逆行列を掛けると単位行列になるような行列のこと。逆行列を持つ正方行列を特に正則行列と言います。逆行列を持つ正方行列を正則行列と言います。

逆行列は、Julia では Base.inv() 関数で求めることができます。

置換行列 (Permutation Matrix)

各行各列にひとつだけ 1 を持ち、残りはすべて 0 であるような行列。この置換行列を他の行列に掛けると要素の置換が起きます。

対角行列 (Diagonal Matrix)

対角成分以外はすべて 0 であるような正方行列を対角行列と言います。

対称行列 (Symmetric Matrix)

対角要素を挟んで要素が対称になっている正方行列を対称行列と言います。

交代行列 (Alternating Matrix)

元の行列の転置行列が元の行列の -1 倍であるような正方行列を交代行列と言います。

直行行列 (Orthogonal Matrix)

元の行列の転置行列が元の行列の逆行列が等しくなるような正方行列を直交行列と言います。回転行列や置換行列は直交行列でもあります。

余因子行列 (Cofactor Matrix)

正方行列のある行 i と列 j を取り除いた小行列の行列式を並べた行列を余因子行列と言います。

三角行列 (Triangular Matrix)

対角成分の上側の要素、あるいは下側の要素がすべて 0 であるような正方行列を三角行列と言います。上側がすべて 0 の場合が下三角行列、下側がすべて 0 である行列が上三角行列と言います。

エルミート行列 (Hermitian Matrix)

エルミート行列とは、元の複素正方行列を複素共役にし、さらに転置した行列が元の行列に等しいような行列のことです。(実数行列である) 対称行列はエルミート行列でもあります。

随伴行列 (Adjugate Matrix)

随伴行列とは複素行列を転置し、複素共役を取った行列です。

ユニタリー行列 (Unitary Matrix)

ユニタリー行列とは、元の複素正方行列を複素共役にし、さらに転置した行列が、元の行列の逆行列に等しくなるような行列のことです。(実数行列である) 直行行列はユニタリー行列でもあります。

いろいろな行列を作る野良パッケージを作ってみました。わかりやすいように関数には日本語の名前を付けてあります。

module Matrixes
using LinearAlgebra

# ゼロ行列
function ゼロ行列(行数, 列数)
    return zeros(行数, 列数)
end

# 単位行列
function 単位行列(行列数)
    E = zeros(行列数, 行列数) + I(行列数)
    return E
end

# 同一値行列
function 同一値行列(値, 行数, 列数)
    サイズ = (行数, 列数)
    fill(値, サイズ)
end

# 浮動小数点乱数行列
function 乱数行列(行数, 列数; 倍率=1.0)
    R = rand(行数, 列数)
    if 倍率 != 1.0
        R = 倍率 * R
    end
    return R
end

# 整数乱数行列
function 整数乱数行列(範囲, 行数, 列数)
    rand(範囲, 行数, 列数)
end

# 対角行列
function 対角行列(ベクトル)
    diagm(ベクトル)
end

# 回転（2次元）
function 回転行列2(角度; ラジアン=true)
    theta = 角度
    if ! ラジアン
        theta = 角度 / 180.0 * π
    end
    return [cos(theta) -sin(theta); sin(theta) cos(theta)]
end

# 回転（3次元）
function 回転行列3(角度; 回転軸='x', ラジアン=true)
    theta = 角度
    if ! ラジアン
        theta = 角度 / 180.0 * π
    end
    if 回転軸 == 'x'
        return [1.0 0.0 0.0; 0.0 cos(theta) -sin(theta); 0.0 sin(theta) cos(theta) ]
    elseif 回転軸 == 'y'
        return [cos(theta) 0.0 sin(theta); 0.0 1.0 0.0; -sin(theta) 0.0 cos(theta)]
    elseif 回転軸 == 'z'
        return [cos(theta) -sin(theta) 0.0; sin(theta) cos[theta] 0.0; 0.0 0.0 1.0]
    else
        throw("回転軸指定エラー")
    end
end

# 拡大・縮小・反転（2次元）= 2次元の対角行列
function 拡大縮小2(対角成分)
    A = zeros(2, 2)
    A[1, 1] = 対角成分[1]
    A[2, 2] = 対角成分[2]
    return A
end

# 拡大・縮小・反転（3次元）= 3次元の対角行列
function 拡大縮小3(対角成分)
    A = zeros(3, 3)
    A[1, 1] = 対角成分[1]
    A[2, 2] = 対角成分[2]
    A[3, 3] = 対角成分[3]
    return A
end

# 平行移動（2次元）
function 平行移動2(移動成分)
    A = zeros(2, 3)
    A[1, 3] = 移動成分[1]
    A[2, 3] = 移動成分[2]
    return A
end

# 平行移動（3次元）
function 平行移動3(移動成分)
    A = zeros(3, 4)
    A[1, 4] = 移動成分[1]
    A[2, 4] = 移動成分[2]
    A[3, 4] = 移動成分[3]
    return A
end

# 行列式
function 行列式(正方行列)
    det(正方行列)
end

# 階数(ランク)
function 階数(正方行列)
    rank(正方行列)
end

function ランク(正方行列)
    rank(正方行列)
end

# トレース
function トレース(正方行列)
    tr(正方行列)
end

# 固有値
function 固有値(正方行列)
    f = eigen(正方行列)
    f.values  # ベクトルの一要素が固有値
end

# 転置行列
function 転置行列(正方行列)
    Transpose(正方行列)
end

# 対称行列
function 対称行列(正方行列; 下側=false)
    if 下側
        Symmetric(正方行列, :L)
    else
        Symmetric(正方行列)
    end
end

# 直交行列 : 回転行列や置換行列は直交行列の一種

# 交代行列 (歪対称行列、Alternative Matrix)
function 交代行列(正方行列)
    A = 正方行列
    n = size(A)[1]
    for i = 1:n
        for j = 1:n
            if i > j
                A[i, j] = -A[j, i]
            elseif i == j
                A[i, i] =  0.0
            else
                # A[i, j] = A[i, j]
            end
        end
    end
    return A
end

# 置換行列
function 置換行列(次元)
    A = zeros(次元, 次元)
    for i = 1:次元
        j = 次元 - i + 1
        A[i, j] = 1.0
    end
    return A
end

# 逆行列
function 逆行列(正方行列)
    inv(正方行列)
end

# 余因子行列 2x2 (Adjugate Matrix)
function 余因子行列2(正方行列)
    A = 正方行列
    return [A[2, 2] -A[1, 2]; -A[2, 1] A[1, 2]]
end

# 余因子行列 3x3 (Adjugate Matrix)
function 余因子行列3(正方行列)
    A = 正方行列
    d11 = det(小行列(A, 1, 1))
    d12 = det(小行列(A, 1, 2))
    d13 = det(小行列(A, 1, 3))
    d21 = det(小行列(A, 2, 1))
    d22 = det(小行列(A, 2, 2))
    d23 = det(小行列(A, 2, 3))
    d31 = det(小行列(A, 3, 1))
    d32 = det(小行列(A, 3, 2))
    d33 = det(小行列(A, 3, 3))
    B = [d11 (-d12) d13; (-d21) d22 (-d23); d31 (-d32) d33]
    return B
end

# 小行列
function 小行列(正方行列, 行, 列)
    n = size(正方行列)[1]
    A = Matrix{Float64}(undef, n - 1, n - 1)
    is = 0; js = 0
    for i = 1:n
        if i == 行
            is += 0
        else
            is += 1
        end
        for j = 1:n
            if j == 列
                js += 0
            else
                js += 1
                if i != 行
                    #@show is, js, i, j
                    A[is, js] = 正方行列[i, j]
                end
            end
        end
        js = 0
    end
    return A
end

# 上三角行列
function 上三角行列(正方行列)
    LinearAlgebra.triu(正方行列)
end

# 下三角行列
function 下三角行列(正方行列)
    LinearAlgebra.tril(正方行列)
end

# エルミート行列
function エルミート行列(複素正方行列)
    LinearAlgebra.Hermitian(複素正方行列)
end

# 随伴行列： 複素行列を転置し、複素共役を取った行列 (A*)
function 随伴行列(複素正方行列)
    A = LinearAlgebra.Transpose(複素正方行列)
    return conj(A)
end

# ユニタリー行列
function ユニタリー行列(複素正方行列)
end

end  # end Matrixes

テストプログラム

# Matrixes のテスト
import Pkg
if !isfile("Project.toml")
    println("Project.toml が見つかりません。")
    exit(1)
end
Pkg.activate(".")
import Matrixes

E = Matrixes.単位行列(3)
println("単位行列(3) = ", E)
Z = Matrixes.ゼロ行列(3, 3)
println("ゼロ行列(3, 3) = ", Z)
D = Matrixes.同一値行列(1.0, 3, 3)
println("同一値行列(1.0, 3, 3) = ", D)
R = Matrixes.乱数行列(3, 3, 倍率=100.0)
println("乱数行列(3, 3, 倍率=100.0) = ", R)
R = Matrixes.整数乱数行列(0:9, 3, 3)
println("整数乱数行列(0:9, 3, 3) = ", R)
A = Matrixes.対角行列([10, 5, -4.9])
println("対角行列([10, 5, -4.9] = ", A)
A = Matrixes.回転行列2(30.0, ラジアン=false)
println("回転行列2(30.0, ラジアン=false) = ", A)
A = Matrixes.回転行列3(30.0, 回転軸='x', ラジアン=false)
println("回転行列3(30.0, 回転軸='x', ラジアン=false) = ", A)
A = Matrixes.拡大縮小2([2.0, 0.5])
println("拡大縮小2([2.0, 0.5]) = ", A)
A = Matrixes.拡大縮小3([2.0, 1.0, 0.5])
println("拡大縮小3([2.0, 1.0, 0.5]) = ", A)
A = Matrixes.平行移動2([2.0, -2.0])
println("平行移動2([2.0, -2.0]) = ", A)
A = Matrixes.平行移動3([2.0, 0.0, -2.0])
println("平行移動3([2.0, 0.0, -2.0]) = ", A)
d = Matrixes.行列式([-1 1; 3 3])
println("行列式([-1 1; 3 3]) = ", d)
r = Matrixes.階数([1 2 3; 0 4 4; -1 6 2])
println("階数([1 2 3; 0 4 4; -1 6 2]) = ", r)
r = Matrixes.ランク([1 2 3; 2 4 6; -1 6 2])
println("階数([1 2 3; 2 4 6; -1 6 2]) = ", r)
t = Matrixes.トレース([1 3 3; 0 1 3; -8 -3 1])
println("トレース([1 3 3; 0 1 3; -8 -3 1]) = ", t)
e = Matrixes.固有値([1 3 3; 0 1 3; -8 -3 1])
println("固有値([1 3 3; 0 1 3; -8 -3 1]) = ", e)
T = Matrixes.転置行列([1 3 3; 0 1 3; -8 -3 1])
println("転置行列([1 3 3; 0 1 3; -8 -3 1]) = ", T)
A = Matrixes.対称行列([1 2 3; 4 5 6; 7 8 9])
println("対称行列([1 2 3; 4 5 6; 7 8 9]) = ", A)
A = Matrixes.交代行列([1 2 3; 4 5 6; 7 8 9])
println("交代行列([1 2 3; 4 5 6; 7 8 9]) = ", A)
T = Matrixes.置換行列(3)
println("置換行列(3) = ", T)
A = Matrixes.逆行列([1 0 2; 3 1 3; -2 4 1])
println("逆行列([1 0 2; 3 1 3; -2 4 1]) = ", A)
A = [1 2 3; 4 5 6; 7 8 9]
B = Matrixes.小行列(A, 1, 1)
println("小行列([1 2 3; 4 5 6; 7 8 9], 1, 1) = ", B)
B = Matrixes.小行列(A, 2, 2)
println("小行列([1 2 3; 4 5 6; 7 8 9], 2, 2) = ", B)
B = Matrixes.小行列(A, 3, 1)
println("小行列([1 2 3; 4 5 6; 7 8 9], 3, 1) = ", B)
B = Matrixes.小行列(A, 1, 3)
println("小行列([1 2 3; 4 5 6; 7 8 9], 1, 3) = ", B)
B = Matrixes.余因子行列2([1 2; 3 4])
println("余因子行列2([1 2; 3 4] = ", B)
B = Matrixes.余因子行列3([1 2 4; 4 5 6; 7 8 9])
println("余因子行列2([1 2 4; 4 5 6; 7 8 9] = ", B)
println("上三角行列([1 2 3; 4 5 6; 7 8 9]) = ", Matrixes.上三角行列([1 2 3; 4 5 6; 7 8 9]))
println("下三角行列([1 2 3; 4 5 6; 7 8 9]) = ", Matrixes.下三角行列([1 2 3; 4 5 6; 7 8 9]))
A = [1 0 2+2im 0 3-3im; 0 4 0 5 0; 6-6im 0 7 0 8+8im; 0 9 0 1 0; 2+2im 0 3-3im 0 4]
@show A
println("随伴行列() = ", Matrixes.随伴行列(A))
println("エルミート行列(A) = ", Matrixes.エルミート行列(A))
#println("ユニタリー行列(A) = ", Matrixes.ユニタリー行列(A))

サンプル１ 回転行列を使って単位ベクトルを1回転したグラフ

# 回転行列
import Pkg
if !isfile("Project.toml")
    println("Project.toml が見つかりません。")
    exit(1)
end
Pkg.activate(".")
import Matrixes
import Plots

座標x = Float64[]
座標y = Float64[]

for 角度 = 0:360
    xy = Matrixes.回転行列2(角度, ラジアン=false) * [1.0, 0.0]
    append!(座標x, xy[1])
    append!(座標y, xy[2])
end
Plots.gr()
Plots.plot(座標x, 座標y, show=true, framestyle=:origin, aspect_ratio = 1; title = "Rotation Matrix", w = 2)
print("> ")
while true
    yield()
end

サンプル１のグラフ表示

 

Julia の勉強：行列とベクトル

Julia は数値計算を強く意識しているため、配列も線形代数を意識したものとなっています。

一般的な配列は Array ですが、特別な配列として Matrix や Vector も定義されています。

数学ではベクトルは単なる次元 (要素数) だけでなく横ベクトルと縦ベクトルがありますが、Julia でも縦か横かを区別されています。

配列

配列の初期化を行うサンプルを次に示します。

# 配列
A = Array{Int32}(undef, 3, 3)
B = [1 2 3; 4 5 6; 7 8 9]
A = B
println(B)
println(A[1, 1], ", ", A[2, 2], ", ", A[3, 3])
v = Array{Int32}(undef, 3)
v[1] = 0
v[2] = 2
v[3] = 5
println(v)
v0 = [0; 2; 5]
println(v == v0)

このサンプルプログラムを実行すると次のようになります。

$ julia array.jl
[1 2 3; 4 5 6; 7 8 9]
1, 5, 9
Int32[0, 2, 5]
true

この結果を見て気になるのが、先頭の要素のインデックスが 1 であるということです。他言語では最初の要素のインデックスは 0 なので注意が必要です。

全要素が同じ配列なら fill を使うとよいです。さらに、その要素が 0 なら zero というのがあります。

julia> zeros(3)
3-element Vector{Float64}:
 0.0
 0.0
 0.0
julia> A = fill(3.0, 2, 2)
2×2 Matrix{Float64}:
 3.0  3.0
 3.0  3.0
julia> println(A)
[3.0 3.0; 3.0 3.0]

その他、すべての要素を 1 にしたり、乱数にしたりと様々な初期化関数が用意されています。

配列リテラルの要素の型を指定するには [ の前に「型」を書いておきます。

(例) Int8[1, 2, 3]

LinearAlgebra

行列やベクトルの計算は Julia だけでも可能ですが、LinearAlgebra (線形代数の意味) という強力なパッケージがあるのでインストールしておきます。

Julia の対話モードで “]” を押してパッケージマネージャを起動して add コマンドを実行します。

pkg> add LinearAlgebra

ベクトル

Julia には Vector という型も用意されています。

これは一次元配列なのですが、数学のベクトルに対応していて横と縦のベクトルがあります。

普通に初期化したベクトルは横ベクトルになります。”,” の代わりに要素の区切りを “;” を使って初期化すると縦ベクトルになります。

横ベクトルは行列の掛け算には使えないので注意が必要です。

ベクトルのサンプル１

# ベクタの計算
using LinearAlgebra
v1 = [2.0, 2.5, 8.0]
v2 = [1.2, 3.5, 7.0]
# ノルム
n =norm(v1)
println(n)
# 足し算, 引き算
v = v1 + v2
println(v)
v = v1 - v2
println(v)
# 内積
a = dot(v1, v2)
println(a)
# 外積
b = cross(v1, v2)
println(b)
# 要素ごとの積
c = v1.*v2
println(c)

実行結果

$ julia vector1.jl
8.616843969807043
[3.2, 6.0, 15.0]
[0.8, -1.0, 1.0]
67.15
[-10.5, -4.4, 4.0]
[2.4, 8.75, 56.0]

ベクトルのサンプル２

# ベクトル (同じ型の一次元配列)
v = Vector{Float64}(undef, 3)
v[1] = 1.0
v[2] = 2.0
v[3] = 4.0
println(v)
# nothing を許すベクトル
v = Vector{Union{Nothing, Float64}}(nothing, 3)
v[1] = 4.0
v[2] = nothing
v[3] = 1.0
println(v)
v = [1; 2; 3]

実行結果

$ julia vector2.jl
[1.0, 2.0, 4.0]
Union{Nothing, Float64}[4.0, nothing, 1.0]
[1, 2, 3]

要素の追加や取り出しなど

ベクトルに要素を追加するには append!() 関数を使います。！が付いていることに注意してください。

julia> a = []
Any[]
julia> append!(a, 1)
1-element Vector{Any}:
 1
julia> append!(a, 2)
2-element Vector{Any}:
 1
 2

ベクトルの先頭に挿入する場合は、pushfirst!() 関数を使います。

julia> pushfirst!(a, -1)
3-element Vector{Any}:
 -1
  1
  2

最後の要素を取り出す場合は、pop!()、先頭の要素を取り出す場合は、popfirst!() 関数を使います。現在のベクトルの長さは length() 関数で取得できます。

julia> a = [1, 2, 3, 4, 5]
5-element Vector{Int64}:
 1
 2
 3
 4
 5
julia> pop!(a)
5
julia> popfirst!(a)
1
julia> length(a)
3
julia> a
3-element Vector{Int64}:
 2
 3
 4

行列

Julia には Vector 同様に Matrix という型も用意されています。これらはどちらも Array をベースにしています。

行列のサンプル１

# 行列の定義と計算
using LinearAlgebra  # add LinearAlgebra が必要

# 3 x 3 行列の定義
m1 = [5 1 4; 3 7 9; 0 1 2]
println(m1)
# 単位行列
e = Matrix{Int64}(I, 3, 3)
println(e)
# 行列 x 行列
m2 = m1 * e
println(m2)
# 行列 x ベクトル
a1 = [2, 2, 2]
m3 = m2 * a1
println(m3)
# 行列式
d = det(m1)
println(d)
# 逆行列
mi = inv(m1)
println(mi)

実行結果

$ julia matrix1.jl
[5 1 4; 3 7 9; 0 1 2]
[1 0 0; 0 1 0; 0 0 1]
[5 1 4; 3 7 9; 0 1 2]
[20, 38, 6]
31.0
[0.16129032258064516 0.06451612903225806 -0.6129032258064516; -0.19354838709677422 0.3225806451612903 -1.064516129032258; 0.09677419354838711 -0.16129032258064516 1.032258064516129]

行列のサンプル２

# 行列 (同じ型の二次元配列)
A = Matrix{Float64}(undef, 3, 3)
A[1, 1] = 1.0
A[1, 2] = 0.0
A[1, 3] = 0.0
A[2, 1] = 0.0
A[2, 2] = 2.0
A[2, 3] = 0.0
A[3, 1] = 0.0
A[3, 2] = 0.0
A[3, 3] = 3.0
println(A)
# nothing を許す行列
B = Matrix{Union{Nothing, Float64}}(nothing, 2, 2)
B[1, 1] = 10.0
B[1, 2] = nothing
B[2, 1] = nothing
B[2, 2] = 20.0

実行結果

$ julia matrix2.jl
[1.0 0.0 0.0; 0.0 2.0 0.0; 0.0 0.0 3.0]
Union{Nothing, Float64}[10.0 nothing; nothing 20.0]

行列のサンプル３

部分配列も簡単に作れます。元の行列で部分行列の範囲を指定するとその部分行列が取得できます。また、行だけ取り出す場合は A[2, :]などのようにします。ここで “:” は全体を意味します。

julia> A = [1 2 3 4; 2 3 4 5; 3 4 5 6; 4 5 6 7]
4×4 Matrix{Int64}:
 1  2  3  4
 2  3  4  5
 3  4  5  6
 4  5  6  7
julia> A[2:3, 2:3]
2×2 Matrix{Int64}:
 3  4
 4  5
julia> A[2:3, 2:4]
2×3 Matrix{Int64}:
 3  4  5
 4  5  6
julia> A[2,:]
4-element Vector{Int64}:
 2
 3
 4
 5