ガウスの消去法

2021年4月8日 19:06

だいたいの理屈はここ

ガウスの消去法

       //前進消去
       public static (double[,] A, double[] b) ForwardElimination(double[,] A, double[] b)
       {
           int row_count = A.GetLength(0);//rowの個数
           int col_count = A.GetLength(1);//colの個数
           if (row_count != col_count) { throw new Exception(); }

           int pivot_count = row_count;
           //拡大係数行列
           var ex_matrix = AddCol(A, b);

           for (int pivot_index = 0; pivot_index < pivot_count - 1; pivot_index++)
           {
               //ピボット選択が必要
               if (ex_matrix[pivot_index, pivot_index] == 0)
               {
                   _pivoting(ex_matrix, pivot_index);
               }

               //ピボットから下の全ての行(pivot_count + 1)
               for (int row_index = pivot_index + 1; row_index < row_count; row_index++)
               {
                   var pivot_term = ex_matrix[pivot_index, pivot_index];
                   //各行の0でない先頭(この項をpivot行の先頭を用いて0にしていく)
                   var leading_coefficient = ex_matrix[row_index, pivot_index];
                   var ratio = -leading_coefficient / pivot_term;//ピボット行との比

                   var pivot_row = GetRow(ex_matrix, pivot_index);//ピボット行
                   var m_row = Multiple(pivot_row, ratio);//row_index行の先頭を0にできるようになったピボット行
                   var target_row = GetRow(ex_matrix, row_index);//先頭を0にすべき行
                   var added = Add(target_row, m_row);//0でない先頭が0に更新された
                   SetRow(ex_matrix, row_index, added);
               }
           }
           return (RemoveCol(ex_matrix, pivot_count), GetCol(ex_matrix, pivot_count));
       }
       //後退代入
       public static double[] BackSubstitution(double[,] A, double[] b)
       {
           double[] x = new double[b.Length];
           for (int pivot_index = x.Length - 1; pivot_index >= 0; pivot_index--)
           {
               var tb = b[pivot_index];
               //初回は実行されない
               for (int col_index = A.GetLength(1) - 1; col_index > pivot_index; col_index--)
               {
                   //pivot行のpivotより右のセルをbに統合していく
                   tb -= A[pivot_index, col_index] * x[col_index];
               }
               tb /= A[pivot_index, pivot_index];
               x[pivot_index] = tb;
           }
           return x;
       }
       
       //ガウス消去法
       public static double[] GaussianElimination(double[,] A, double[] b)
       {
           (var tA, var tb) = ForwardElimination(A, b);
           return BackSubstitution(tA, tb);
       }

Pivot

        public static int _get_pivoting_index(double[,] A, int pivot_index)
       {
           int row_count = A.GetLength(0);//rowの個数

           double abs_max = 0;
           int abs_max_row_index = pivot_index;
           for (int row_index = pivot_index; row_index < row_count; row_index++)
           {
               var temp = Math.Abs(A[row_index, pivot_index]);
               if (abs_max < temp)
               {
                   abs_max = temp;
                   abs_max_row_index = row_index;
               }
           }
           return abs_max_row_index;
       }

       //行交換実行関数
       //行交換したかしてないかは行列式の計算に関連する(したら*-1)
       //行t交換した=> return 1;
       //行交換してない=> return 0;
       public static int _pivoting(double[,] A, int pivot_index)
       {
           var abs_max_index = _get_pivoting_index(A, pivot_index);
           if (pivot_index == abs_max_index) { return 0; }

           if(ChangeRow(A, pivot_index, abs_max_index)==0)
           {
               //Changeしてない
               return 0;
           }
           else
           {
               //Changeした
               //pivot交換した
               return 1;
           }
       }

行と列のGet,Set

        #region 行と列の操作

       //一つの行列を左右に割る
       public static (double[,], double[,]) SplitEx(double[,] matrix, int split_col_count)
       {
           int row_count = matrix.GetLength(0);
           int col_count = matrix.GetLength(1);
           double[,] left = new double[row_count, split_col_count];
           double[,] right = new double[row_count, col_count - split_col_count];
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   if (split_col_count - 1 < col_index)
                   {
                       right[row_index, col_index - split_col_count] = matrix[row_index, col_index];
                   }
                   else
                   {
                       left[row_index, col_index] = matrix[row_index, col_index];
                   }
               }
           }
           return (left, right);
       }

       public static double[] GetRow(double[,] matrix, int row_index)
       {
           //int row_count = matrix.GetLength(0);//rowの個数
           int col_count = matrix.GetLength(1);//colの個数
           double[] ret = new double[col_count];
           //row固定、col走査
           for (int col_index = 0; col_index < col_count; col_index++)
           {
               ret[col_index] = matrix[row_index, col_index];
           }
           return ret;
       }
       
       public static double[] GetCol(double[,] matrix, int col_index)
       {
           int row_count = matrix.GetLength(0);//rowの個数
           //int col_count = matrix.GetLength(1);//colの個数
           double[] ret = new double[row_count];
           //col固定、row走査
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               ret[row_index] = matrix[row_index, col_index];
           }
           return ret;
       }
       public static void SetRow(double[,] matrix, int row_index, double[] new_row)
       {
           //int row_count = matrix.GetLength(0);//rowの個数
           int col_count = matrix.GetLength(1);//colの個数
           for (int col_index = 0; col_index < col_count; col_index++)
           {
               //入力配列が短い場合は0で埋める
               //必然的に、長い場合は切り落とされる
               if (new_row.Length - 1 < col_index)
               {
                   matrix[row_index, col_index] = 0;
               }
               matrix[row_index, col_index] = new_row[col_index];
           }
       }
       public static void SetCol(double[,] matrix, int col_index, double[] new_col)
       {
           int row_count = matrix.GetLength(0);//rowの個数
           //int col_count = matrix.GetLength(1);//colの個数
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               //入力配列が短い場合は0で埋める
               //必然的に、長い場合は切り落とされる
               if (new_col.Length - 1 < row_index)
               {
                   matrix[row_index, col_index] = 0;
               }
               matrix[row_index, col_index] = new_col[row_index];
           }
       }

       public static double[,] GetCols(double[,] matrix, int start_col_index)
       {
           int row_count = matrix.GetLength(0);//rowの個数
           int col_count = matrix.GetLength(1);//colの個数
           int ret_row_count = row_count;
           int ret_col_count = col_count - start_col_index;
           double[,] ret = new double[ret_row_count, ret_col_count];
           //ret側走査
           for (int row_index = 0; row_index < ret_row_count; row_index++)
           {
               //row固定、col走査
               for (int col_index = 0; col_index < ret_col_count; col_index++)
               {
                   ret[row_index, col_index] = matrix[row_index, start_col_index + col_index];
               }
           }
           return ret;
       }
       
       //指定indexから最後までとっちゃう
       public static double[,] GetRows(double[,] matrix, int start_row_index)
       {
           int row_count = matrix.GetLength(0);//rowの個数
           int col_count = matrix.GetLength(1);//colの個数
           int ret_row_count = row_count - start_row_index;
           int ret_col_count = col_count;
           double[,] ret = new double[ret_row_count, ret_col_count];
           //ret側走査
           for (int row_index = 0; row_index < ret_row_count; row_index++)
           {
               //row固定、col走査
               for (int col_index = 0; col_index < ret_col_count; col_index++)
               {
                   ret[row_index, col_index] = matrix[start_row_index + row_index, col_index];
               }
           }
           return ret;
       }

Add

       public static double[,] AddRow(double[,] matrix, double[] new_row)
       {
           int row_count = matrix.GetLength(0);//rowの個数
           int col_count = matrix.GetLength(1);//colの個数
           double[,] ret = new double[row_count+1, col_count];
           for (int row_index = 0; row_index < row_count+1; row_index++)
           {
               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   if(row_index==row_count)
                   {
                       ret[row_index, col_index] = new_row[col_index];
                   }
                   else
                   {
                       ret[row_index, col_index] = matrix[row_index, col_index];
                   }
               }
           }
           return ret;
       }
       
       public static double[,] AddRow(double[,] matrix, double[,] new_rows)
       {
           int source_row_count = matrix.GetLength(0);
           int total_row_count = matrix.GetLength(0) + new_rows.GetLength(0);//rowの個数
           int total_col_count = matrix.GetLength(1);//colの個数
           double[,] ret = new double[total_row_count, total_col_count];
           for (int row_index = 0; row_index < total_row_count; row_index++)
           {
               for (int col_index = 0; col_index < total_col_count; col_index++)
               {
                   if(row_index< source_row_count)
                   {
                       ret[row_index, col_index] = matrix[row_index, col_index];
                   }
                   else
                   {
                       ret[row_index, col_index] = new_rows[row_index - source_row_count, col_index];
                   }
               }
           }
           return ret;
       }
       
       public static double[,] AddCol(double[,] matrix, double[] new_col)
       {
           int row_count = matrix.GetLength(0);//rowの個数
           int col_count = matrix.GetLength(1);//colの個数
           double[,] ret = new double[row_count, col_count+1];
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               for (int col_index = 0; col_index < col_count+1; col_index++)
               {
                   if (col_index == col_count)
                   {
                       ret[row_index, col_index] = new_col[row_index];
                   }
                   else
                   {
                       ret[row_index, col_index] = matrix[row_index, col_index];
                   }
               }
           }
           return ret;
       }
       
       public static double[,] AddCol(double[,] matrix, double[,] new_cols)
       {
           int source_col_count = matrix.GetLength(1);//colの個数
           int row_count = matrix.GetLength(0);//rowの個数
           int total_col_count = matrix.GetLength(1)+new_cols.GetLength(1);//colの個数
           double[,] ret = new double[row_count, total_col_count];
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               for (int col_index = 0; col_index < total_col_count; col_index++)
               {
                   if (col_index < source_col_count)
                   {
                       ret[row_index, col_index] = matrix[row_index, col_index];
                   }
                   else
                   {
                       ret[row_index, col_index] = new_cols[row_index, col_index - source_col_count];
                   }                  
               }
           }
           return ret;
       }

Remove

       public static double[,] RemoveRow(double[,] matrix, int remove_row_index)
       {
           int row_count = matrix.GetLength(0);
           int col_count = matrix.GetLength(1);
           double[,] ret = new double[row_count - 1, col_count];
           int skip_count = 0;
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               //削除対象の行
               if (row_index == remove_row_index)
               {
                   skip_count += 1;
                   continue;
               }

               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   ret[row_index - skip_count, col_index] = matrix[row_index, col_index];
               }
           }
           return ret;
       }
       
       public static double[,] RemoveCol(double[,] matrix, int remove_col_index)
       {
           int row_count = matrix.GetLength(0);
           int col_count = matrix.GetLength(1);
           double[,] ret = new double[row_count, col_count-1];
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               int skip_count = 0;
               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   //削除対象の列
                   if (col_index == remove_col_index)
                   {
                       skip_count += 1;
                       continue;
                   }
                   ret[row_index, col_index - skip_count] = matrix[row_index , col_index];
               }
           }
           return ret;
       }

Insert

       public static double[,] InsertRow(double[,] matrix, int insert_row_index, double[] new_row)
       {
           double[,] ret = new double[matrix.GetLength(0) + 1, matrix.GetLength(1)];
           int row_count = ret.GetLength(0);
           int col_count = ret.GetLength(1);
           int skip_count = 0;
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   //追加対象の行
                   if (row_index == insert_row_index)
                   {
                       ret[row_index, col_index] = new_row[col_index];
                       if(col_index==col_count-1)
                       {
                           skip_count += 1;
                       }
                       continue;
                   }
                   else
                   {
                       ret[row_index, col_index] = matrix[row_index - skip_count, col_index];
                   }               
               }
           }
           return ret;
       }
       
       public static double[,] InsertCol(double[,] matrix, int insert_col_index, double[] new_col)
       {
           double[,] ret = new double[matrix.GetLength(0), matrix.GetLength(1)+1];
           int row_count = ret.GetLength(0);
           int col_count = ret.GetLength(1);         
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               int skip_count = 0;
               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   //追加対象の列
                   if (col_index == insert_col_index)
                   {
                       ret[row_index, col_index] = new_col[row_index];
                       skip_count += 1;
                   }
                   else
                   {
                       ret[row_index, col_index] = matrix[row_index, col_index - skip_count];
                   }
               }
           }
           return ret;
       }

Change,Paste,Replace





       #region Change

       //Changeした=> return 1;
       //Changeしてない=> return 0;
       public static int ChangeRow(double[,] matrix, int row_index1, int row_index2)
       {
           if (row_index1 == row_index2) { return 0; }
           var temp_row1 = GetRow(matrix, row_index1);
           var temp_row2 = GetRow(matrix, row_index2);
           SetRow(matrix, row_index1, temp_row2);
           SetRow(matrix, row_index2, temp_row1);
           return 1;
       }

       public static void _changeRow(double[,] matrix, int row_index1, int row_index2)
       {
           int col_count = matrix.GetLength(1);
           for (int col_index = 0; col_index < col_count; col_index++)
           {
               double temp = matrix[row_index1, col_index];
               matrix[row_index1, col_index] = matrix[row_index2, col_index];
               matrix[row_index2, col_index] = temp;
           }
       }

       public static int ChangeCol(double[,] matrix, int col_index1, int col_index2)
       {
           if (col_index1 == col_index2) { return 0; }
           var temp_col1 = GetCol(matrix, col_index1);
           var temp_col2 = GetCol(matrix, col_index2);
           SetCol(matrix, col_index1, temp_col2);
           SetCol(matrix, col_index2, temp_col1);
           return 1;
       }

       public static void _changeCol(double[,] matrix, int col_index1, int col_index2)
       {
           //colを固定してrowを走査する
           int row_count = matrix.GetLength(0);
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               double temp = matrix[row_index, col_index1];
               matrix[row_index, col_index1] = matrix[row_index, col_index2];
               matrix[row_index, col_index2] = temp;
           }
       }

       public static void _changeComponent(double[] vector, int index1, int index2)
       {
           double temp = vector[index1];
           vector[index1] = vector[index2];
           vector[index2] = temp;
       } 

       #endregion

       //newして戻り値返すタイプ
       //昨日的にはSetと変わらない
       public static double[,] PasteRow(double[,] matrix, int row_index, double[] new_row)
       {
           int row_count = matrix.GetLength(0);//rowの個数
           int col_count = matrix.GetLength(1);//colの個数
           double[,] ret = new double[row_count, col_count];
           for (int col_index = 0; col_index < col_count; col_index++)
           {
               //入力配列が短い場合は0で埋める
               //必然的に、長い場合は切り落とされる
               if (new_row.Length - 1 < col_index)
               {
                   ret[row_index, col_index] = 0;
               }
               ret[row_index, col_index] = new_row[col_index];
           }
           return ret;
       }
       
       public static double[,] PasteCol(double[,] matrix, int col_index, double[] new_col)
       {
           int row_count = matrix.GetLength(0);//rowの個数
           int col_count = matrix.GetLength(1);//colの個数
           double[,] ret = new double[row_count, col_count];
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               //入力配列が短い場合は0で埋める
               //必然的に、長い場合は切り落とされる
               if (new_col.Length - 1 < row_index)
               {
                   ret[row_index, col_index] = 0;
               }
               ret[row_index, col_index] = new_col[row_index];
           }
           return ret;
       }

       #region Reprace

       public static void Replace(double[] vector, params int[] indexes)
       {
           for (int col_index = 0; col_index < vector.Length; col_index++)
           {
               double temp = vector[col_index];
               vector[col_index] = vector[indexes[col_index]];
               vector[indexes[col_index]] = temp;
           }
       }
       
       //記録しておいた入れ換えを使用
       public static void ReplaceRow(double[,] matrix, params int[] indexes)
       {
           int row_count = matrix.GetLength(0);
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               _changeRow(matrix, row_index, indexes[row_index]);
           }
       }
       
       //記録しておいた入れ換えを使用
       public static void ReplaceCol(double[,] matrix, params int[] indexes)
       {
           int col_count = matrix.GetLength(1);
           for (int col_index = 0; col_index < col_count; col_index++)
           {
               _changeCol(matrix, col_index, indexes[col_index]);
           }
       }

       //今の所使わなくて良い
       public static void ReplaceBackRow(double[,] matrix, params int[] indexes)
       {
           int row_count = matrix.GetLength(0);
           for (int row_index = row_count - 1; row_index >= 0; row_index--)
           {
               _changeRow(matrix, row_index, indexes[row_index]);
           }
       }
       
       public static void ReplaceBackCol(double[,] matrix, params int[] indexes)
       {
           int col_count = matrix.GetLength(1);
           for (int col_index = col_count - 1; col_index >= 0; col_index--)
           {
               _changeCol(matrix, col_index, indexes[col_index]);
           }
       }

       #endregion

       #endregion

行列の積

       #region 行列の積

       public static double[,] MatrixMultiple_VH(double[] vertical, double[] horizontal)
       {
           int row_count = vertical.Length;
           int col_count = horizontal.Length;
           double[,] ret = new double[row_count, col_count];
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   ret[row_index, col_index] = vertical[row_index] * horizontal[col_index];
               }
           }
           return ret;
       }

       public static double[] MatrixMultiple(double[,] m, double[] v)
       {
           int row_count = m.GetLength(0);//row
           int col_count = m.GetLength(1);//col
           double[] ret = new double[col_count];
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               //右に向かって走査
               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   ret[row_index] += m[row_index, col_index] * v[col_index];
               }
           }
           return ret;
       }

       public static double[] MatrixMultiple(double[] h, double[,] m)
       {
           int row_count = m.GetLength(0);//row
           int col_count = m.GetLength(1);//col
           double[] ret = new double[col_count];
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               //右に向かって走査
               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   ret[row_index] += h[col_index] * m[row_index, col_index];
               }
           }
           return ret;
       }

       public static double[,] MatrixMultiple(double[,] m1, double[,] m2)
       {
           int row_count = m1.GetLength(0);//row
           int col_count = m1.GetLength(1);//col
           double[,] ret = new double[row_count, col_count];
           for (int row_index = 0; row_index < row_count; row_index++)
           {
               //右に向かって走査
               for (int col_index = 0; col_index < col_count; col_index++)
               {
                   double s = 0;
                   for (int k = 0; k < col_count; k++)
                   {
                       s += m1[row_index, k] * m2[k, col_index];
                   }
                   ret[row_index, col_index] = s;
               }
           }
           return ret;
       }

       #endregion

Elementwise


       public static double[,] ElementWise_Add(double[,] m1, double[,] m2)
       {
           return _calcuElementWise(m1, m2, _add);
       }
       
       public static double[,] ElementWise_Sub(double[,] m1, double[,] m2)
       {
           return _calcuElementWise(m1, m2, _sub);
       }
       
       public static double[,] ElementWise_Multiple(double[,] m1, double[,] m2)
       {
           return _calcuElementWise(m1, m2, _multiple);
       }
       
       public static double[,] ElementWise_Divide(double[,] m1, double[,] m2)
       {
           return _calcuElementWise(m1, m2, _divide);
       }

       
       //配列の定数倍
       public static double[] _calcuElementWise(double[] v1, double c, Func<double, double, double> operation)
       {
           double[] ret = new double[v1.Length];

           for (int i = 0; i < v1.Length; i++)
           {
               ret[i] = operation(v1[i], c);
           }
           return ret;
       }
       
       //配列×配列
       public static double[] _calcuElementWise(double[] v1, double[] v2, Func<double, double, double> operation)
       {
           if (v1.Length != v2.Length) { return null; }

           double[] ret = new double[v1.Length];

           for (int i = 0; i < v1.Length; i++)
           {
               ret[i] = operation(v1[i], v2[i]);
           }
           return ret;
       }
       
       //行列の定数倍
       public static double[,] _calcuElementWise(double[,] m1, double c, Func<double, double, double> operation)
       {
           double[,] ret = new double[m1.GetLength(0), m1.GetLength(1)];

           for (int m = 0; m < m1.GetLength(0); m++)
           {
               //右に向かって走査
               for (int n = 0; n < m1.GetLength(1); n++)
               {
                   ret[m, n] = operation(m1[m, n], c);
               }
           }
           return ret;
       }

       //行列×縦ベクトル
       public static double[] _calcuElementWise(double[,] m, double[] v, Func<double, double, double> operation)
       {
           double[] ret = new double[v.Length];

           for (int row_index = 0; row_index < m.GetLength(0); row_index++)
           {
               //右に向かって走査
               for (int col_index = 0; col_index < m.GetLength(1); col_index++)
               {
                   ret[row_index] = operation(m[row_index, col_index], v[col_index]);
               }
           }
           return ret;
       }

       //行列×行列
       public static double[,] _calcuElementWise(double[,] m1, double[,] m2, Func<double, double, double> operation)
       {
           if (m1.GetLength(0) != m2.GetLength(0)) { return null; }
           if (m1.GetLength(1) != m2.GetLength(1)) { return null; }

           double[,] ret = new double[m1.GetLength(0), m1.GetLength(1)];

           //一行ずつ取り出す[同型なのでどちらの行列を基準にしても良い]
           for (int m = 0; m < m1.GetLength(0); m++)
           {
               //右に向かって走査
               for (int n = 0; n < m1.GetLength(1); n++)
               {
                   ret[m, n] = operation(m1[m, n], m2[m, n]);
               }
           }
           return ret;
       }

配列の定数倍(Elementwise)

       //配列の定数倍
       public static double[] Add(double[] v, double c)
       {
           return _calcuElementWise(v, c, _add);
       }
       public static double[] Sub(double[] v, double c)
       {
           return _calcuElementWise(v, c, _sub);
       }
       //これはElementWiseであって行列の積ではない
       public static double[] Multiple(double c, double[] v)
       {
           return _calcuElementWise(v, c, _multiple);
       }
       public static double[] Multiple(double[] v, double c)
       {
           return _calcuElementWise(v, c, _multiple);
       }
       public static double[] Divide(double[] v, double c)
       {
           return _calcuElementWise(v, c, _divide);
       }
       public static double[] Mod(double[] v, double c)
       {
           return _calcuElementWise(v, c, _mod);
       }

配列×配列(Elementwise)

       //配列×配列
       public static double[] Add(double[] v1, double[] v2)
       {
           return _calcuElementWise(v1, v2, _add);
       }
       public static double[] Sub(double[] v1, double[] v2)
       {
           return _calcuElementWise(v1, v2, _sub);
       }
       //これはElementWiseであって行列の積ではない
       public static double[] Multiple(double[] v1, double[] v2)
       {
           return _calcuElementWise(v1, v2, _multiple);
       }
       public static double[] Divide(double[] v1, double[] v2)
       {
           return _calcuElementWise(v1, v2, _divide);
       }
       public static double[] Mod(double[] v1, double[] v2)
       {
           return _calcuElementWise(v1, v2, _mod);
       }

行列の定数倍(Elementwise)

       //行列の定数倍
       public static double[,] Add(double[,] m1, double c)
       {
           return _calcuElementWise(m1, c, _add);
       }
       public static double[,] Sub(double[,] m1, double c)
       {
           return _calcuElementWise(m1, c, _sub);
       }
       //これはElementWiseであって行列の積ではない
       public static double[,] Multiple(double c, double[,] m1)
       {
           return _calcuElementWise(m1, c, _multiple);
       }
       public static double[,] Multiple(double[,] m1, double c)
       {
           return _calcuElementWise(m1, c, _multiple);
       }
       public static double[,] Divide(double[,] m1, double c)
       {
           return _calcuElementWise(m1, c, _divide);
       }
       public static double[,] Mod(double[,] m1, double c)
       {
           return _calcuElementWise(m1, c, _mod);
       }

横ベクトル×行列(Elementwise)

       //横ベクトル×行列
       public static double[] Add(double[] v, double[,] m)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _add);
       }
       public static double[] Sub(double[] v, double[,] m)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _sub);
       }
       //これはElementWiseであって行列の積ではない
       public static double[] Multiple(double[] v, double[,] m)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _multiple);
       }
       public static double[] Divide(double[] v, double[,] m)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _divide);
       }
       public static double[] Mod(double[] v, double[,] m)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _mod);
       }

行列×縦ベクトル(Elementwise)

       //行列×縦ベクトル
       public static double[] Add(double[,] m, double[] v)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _add);
       }
       public static double[] Sub(double[,] m, double[] v)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _sub);
       }
       //これはElementWiseであって行列の積ではない
       public static double[] Multiple(double[,] m, double[] v)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _multiple);
       }
       public static double[] Divide(double[,] m, double[] v)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _divide);
       }
       public static double[] Mod(double[,] m, double[] v)
       {
           //vは縦ベクトル
           return _calcuElementWise(m, v, _mod);
       }

行列×行列(Elementwise)

       //行列×行列
       public static double[,] Add(double[,] m1, double[,] m2)
       {
           return _calcuElementWise(m1, m2, _add);
       }
       public static double[,] Sub(double[,] m1, double[,] m2)
       {
           return _calcuElementWise(m1, m2, _sub);
       }
       //これはElementWiseであって行列の積ではない
       public static double[,] Multiple(double[,] m1, double[,] m2)
       {
           return _calcuElementWise(m1, m2, _multiple);
       }
       public static double[,] Divide(double[,] m1, double[,] m2)
       {
           return _calcuElementWise(m1, m2, _divide);
       }
       public static double[,] Mod(double[,] m1, double[,] m2)
       {
           return _calcuElementWise(m1, m2, _mod);
       }

式


       #region 式

       public static double _add(double a, double b)
       {
           return a + b;
       }
       public static double _sub(double a, double b)
       {
           return a - b;
       }
       public static double _multiple(double a, double b)
       {
           return a * b;
       }
       public static double _divide(double a, double b)
       {
           if (b == 0)
           {
               return 0;
           }
           return a / b;
       }
       public static double _mod(double a, double b)
       {
           return a % b;
       }

       public static double _exp(double a)
       {
           return Math.Exp(a);
       }
       //床関数
       public static double _floor(double a)
       {
           return (int)a;
       }
       //天井関数
       public static double _ceil(double a)
       {
           return (int)a + 1;
       }
       //四捨五入
       public static double _round(double a)
       {
           return _floor(a + 0.5);
       }
       //aの小数部分
       public static double _fraction(double a)
       {
           return a - (int)a;
       }
       public static double _rfraction(double a)
       {
           return 1 - (a - (int)a);
       }

       #endregion

processing:ガウスの消去法



class ExMatrix
{
 double[][] A;
 double[] b;
 
 public ExMatrix(double[][] A_, double[] b_)
 {
   this.A = A_;
   this.b = b_;
 }
}

int _get_pivoting_index(double[][] A, int pivot_index)
{
   //pivot列の
   //pivot+1行以降で
   //絶対値最大の行のインデックスを返す

   int row_count = A.length;//rowの個数

   double abs_max = 0;
   int abs_max_row_index = pivot_index;

   //pivotエレメントは最後の一つまで走査される場合と
   //一つ手前で止まる場合がある
   //if (A.GetLength(0) - 1 < abs_max_row_index) { abs_max_row_index = pivot_index; }

   for (int row_index = pivot_index; row_index < row_count; row_index++)
   {
       double temp = Math.abs(A[row_index][pivot_index]);
       if (abs_max < temp)
       {
           abs_max = temp;
           abs_max_row_index = row_index;
       }
   }
   return abs_max_row_index;
}
       
//行交換実行関数
//行交換したかしてないかは行列式の計算に関連する(したら*-1)
//行t交換した=> return 1;
//行交換してない=> return 0;
int _pivoting(double[][] A, int pivot_index)
{
   int abs_max_index = _get_pivoting_index(A, pivot_index);
   if (pivot_index == abs_max_index) { return 0; }

   if(ChangeRow(A, pivot_index, abs_max_index)==0)
   {
       //Changeしてない
       return 0;
   }
   else
   {
       //Changeした
       //pivot交換した
       return 1;
   }
}        

//前進消去
ExMatrix ForwardElimination(double[][] A, double[] b)
{
   int row_count = A.length;//rowの個数
   int col_count = A[0].length;//colの個数
   //if (row_count != col_count) { throw new Exception(); }

   int pivot_count = row_count;
   //拡大係数行列
   double[][] ex_matrix = AddCol(A, b);

   for (int pivot_index = 0; pivot_index < pivot_count - 1; pivot_index++)
   {
       //ピボット選択が必要
       if (ex_matrix[pivot_index][pivot_index] == 0)
       {
           _pivoting(ex_matrix, pivot_index);
       }

       //ピボットから下の全ての行(pivot_count + 1)
       for (int row_index = pivot_index + 1; row_index < row_count; row_index++)
       {
           double pivot_term = ex_matrix[pivot_index][pivot_index];

           //各行の0でない先頭(この項をpivot行の先頭を用いて0にしていく)
           double leading_coefficient = ex_matrix[row_index][pivot_index];
           double ratio = -leading_coefficient / pivot_term;//ピボット行との比

           double[] pivot_row = GetRow(ex_matrix, pivot_index);//ピボット行
           double[] m_row = Multiple(pivot_row, ratio);//row_index行の先頭を0にできるようになったピボット行
           double[] target_row = GetRow(ex_matrix, row_index);//先頭を0にすべき行
           double[] added = Add(target_row, m_row);//0でない先頭が0に更新された

           SetRow(ex_matrix, row_index, added);
       }
   }
   
   //拡大係数行列を返す
   ExMatrix ret = new ExMatrix(RemoveCol(ex_matrix, pivot_count), GetCol(ex_matrix, pivot_count));    
   return ret;
}
//後退代入
double[] BackSubstitution(double[][] A, double[] b)
{
   double[] x = new double[b.length];

   for (int pivot_index = x.length - 1; pivot_index >= 0; pivot_index--)
   {
       double tb = b[pivot_index];
       //初回は実行されない
       for (int col_index = A[0].length - 1; col_index > pivot_index; col_index--)
       {
           //pivot行のpivotより右のセルをbに統合していく
           tb -= A[pivot_index][col_index] * x[col_index];
       }
       //Amn*Xmn=Bm(二週目以降はこのBmをcol-pivotループで構築する)
       //Xmn=Bm/Amn
       tb /= A[pivot_index][pivot_index];

       x[pivot_index] = tb;
   }
   return x;
}

//ガウス消去法
double[] GaussianElimination(double[][] A, double[] b)
{
   ExMatrix ex = ForwardElimination(A, b);
   return BackSubstitution(ex.A, ex.b);
}

processing:Get,Set,Add,Remove,Insert


//一つの行列を左右に割る
ArrayList<double[][]> SplitEx(double[][] matrix, int split_col_count)
{
 ArrayList<double[][]> ret = new ArrayList<double[][]>();
 int row_count = matrix.length;
 int col_count = matrix[0].length;
 double[][] left = new double[row_count][split_col_count];
 double[][] right = new double[row_count][col_count - split_col_count];
 for (int row_index = 0; row_index < row_count; row_index++)
 {
     for (int col_index = 0; col_index < col_count; col_index++)
     {
         if (split_col_count - 1 < col_index)
         {
             right[row_index][col_index - split_col_count] = matrix[row_index][col_index];
         }
         else
         {
             left[row_index][col_index] = matrix[row_index][col_index];
         }
     }
 }
 ret.add(left);
 ret.add(right);
 return ret;
}

double[] GetRow(double[][] matrix, int row_index)
{
   int col_count = matrix[0].length;//colの個数
   double[] ret = new double[col_count];
   //row固定、col走査
   for (int col_index = 0; col_index < col_count; col_index++)
   {
       ret[col_index] = matrix[row_index][col_index];
   }
   return ret;
}

double[] GetCol(double[][] matrix, int col_index)
{
   int row_count = matrix.length;//rowの個数
   double[] ret = new double[row_count];
   //col固定、row走査
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       ret[row_index] = matrix[row_index][col_index];
   }
   return ret;
}

void SetRow(double[][] matrix, int row_index, double[] new_row)
{
   int col_count = matrix[0].length;//colの個数
   for (int col_index = 0; col_index < col_count; col_index++)
   {
       //入力配列が短い場合は0で埋める
       //必然的に、長い場合は切り落とされる
       if (new_row.length - 1 < col_index)
       {
           matrix[row_index][col_index] = 0;
       }
       matrix[row_index][col_index] = new_row[col_index];
   }
}

void SetCol(double[][] matrix, int col_index, double[] new_col)
{
   int row_count = matrix.length;//rowの個数
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       //入力配列が短い場合は0で埋める
       //必然的に、長い場合は切り落とされる
       if (new_col.length - 1 < row_index)
       {
           matrix[row_index][col_index] = 0;
       }
       matrix[row_index][col_index] = new_col[row_index];
   }
}

double[][] GetCols(double[][] matrix, int start_col_index)
{
   int row_count = matrix.length;//rowの個数
   int col_count = matrix[0].length;//colの個数
   int ret_row_count = row_count;
   int ret_col_count = col_count - start_col_index;
   double[][] ret = new double[ret_row_count][ret_col_count];
   //ret側走査
   for (int row_index = 0; row_index < ret_row_count; row_index++)
   {
       //row固定、col走査
       for (int col_index = 0; col_index < ret_col_count; col_index++)
       {
           ret[row_index][col_index] = matrix[row_index][start_col_index + col_index];
       }
   }
   return ret;
}

//指定indexから最後までとっちゃう
double[][] GetRows(double[][] matrix, int start_row_index)
{
   int row_count = matrix.length;//rowの個数
   int col_count = matrix[0].length;//colの個数
   int ret_row_count = row_count - start_row_index;
   int ret_col_count = col_count;
   double[][] ret = new double[ret_row_count][ret_col_count];
   //ret側走査
   for (int row_index = 0; row_index < ret_row_count; row_index++)
   {
       //row固定、col走査
       for (int col_index = 0; col_index < ret_col_count; col_index++)
       {
           ret[row_index][col_index] = matrix[start_row_index + row_index][col_index];
       }
   }
   return ret;
}

double[][] AddRow(double[][] matrix, double[] new_row)
{
   int row_count = matrix.length;//rowの個数
   int col_count = matrix[0].length;//colの個数
   double[][] ret = new double[row_count+1][col_count];
   for (int row_index = 0; row_index < row_count+1; row_index++)
   {
       for (int col_index = 0; col_index < col_count; col_index++)
       {
           if(row_index==row_count)
           {
               ret[row_index][col_index] = new_row[col_index];
           }
           else
           {
               ret[row_index][col_index] = matrix[row_index][col_index];
           }
       }
   }
   return ret;
}

double[][] AddRow(double[][] matrix, double[][] new_rows)
{
   int source_row_count = matrix.length;
   int total_row_count = matrix.length + new_rows.length;//rowの個数
   int total_col_count = matrix[0].length;//colの個数
   double[][] ret = new double[total_row_count][total_col_count];
   for (int row_index = 0; row_index < total_row_count; row_index++)
   {
       for (int col_index = 0; col_index < total_col_count; col_index++)
       {
           if(row_index< source_row_count)
           {
               ret[row_index][col_index] = matrix[row_index][col_index];
           }
           else
           {
               ret[row_index][col_index] = new_rows[row_index - source_row_count][col_index];
           }
       }
   }
   return ret;
}
double[][] AddCol(double[][] matrix, double[] new_col)
{
   int row_count = matrix.length;//rowの個数
   int col_count = matrix[0].length;//colの個数
   double[][] ret = new double[row_count][col_count+1];
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       for (int col_index = 0; col_index < col_count+1; col_index++)
       {
           if (col_index == col_count)
           {
               ret[row_index][col_index] = new_col[row_index];
           }
           else
           {
               ret[row_index][col_index] = matrix[row_index][col_index];
           }
       }
   }
   return ret;
}
double[][] AddCol(double[][] matrix, double[][] new_cols)
{
   int source_col_count = matrix[0].length;//colの個数
   int row_count = matrix.length;//rowの個数
   int total_col_count = matrix[0].length+new_cols[0].length;//colの個数
   double[][] ret = new double[row_count][total_col_count];
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       for (int col_index = 0; col_index < total_col_count; col_index++)
       {
           if (col_index < source_col_count)
           {
               ret[row_index][col_index] = matrix[row_index][col_index];
           }
           else
           {
               ret[row_index][col_index] = new_cols[row_index][col_index - source_col_count];
           }                  
       }
   }
   return ret;
}

double[][] RemoveRow(double[][] matrix, int remove_row_index)
{
   int row_count = matrix.length;
   int col_count = matrix[0].length;
   double[][] ret = new double[row_count - 1][col_count];
   int skip_count = 0;
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       //削除対象の行
       if (row_index == remove_row_index)
       {
           skip_count += 1;
           continue;
       }
       for (int col_index = 0; col_index < col_count; col_index++)
       {
           ret[row_index - skip_count][col_index] = matrix[row_index][col_index];
       }
   }
   return ret;
}

double[][] RemoveCol(double[][] matrix, int remove_col_index)
{
   int row_count = matrix.length;
   int col_count = matrix[0].length;
   double[][] ret = new double[row_count][col_count-1];
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       int skip_count = 0;
       for (int col_index = 0; col_index < col_count; col_index++)
       {
           //削除対象の列
           if (col_index == remove_col_index)
           {
               skip_count += 1;
               continue;
           }
           ret[row_index][col_index - skip_count] = matrix[row_index][col_index];
       }
   }
   return ret;
}

double[][] InsertRow(double[][] matrix, int insert_row_index, double[] new_row)
{
   double[][] ret = new double[matrix.length + 1][ matrix[0].length];
   int row_count = ret.length;
   int col_count = ret[0].length;
   int skip_count = 0;
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       for (int col_index = 0; col_index < col_count; col_index++)
       {
           //追加対象の行
           if (row_index == insert_row_index)
           {
               ret[row_index][col_index] = new_row[col_index];
               if(col_index==col_count-1)
               {
                   skip_count += 1;
               }
               continue;
           }
           else
           {
               ret[row_index][col_index] = matrix[row_index - skip_count][col_index];
           }               
       }
   }
   return ret;
}

double[][] InsertCol(double[][] matrix, int insert_col_index, double[] new_col)
{
   double[][] ret = new double[matrix.length][matrix[0].length+1];
   int row_count = ret.length;
   int col_count = ret[0].length;  
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       int skip_count = 0;
       for (int col_index = 0; col_index < col_count; col_index++)
       {
           //追加対象の列
           if (col_index == insert_col_index)
           {
               ret[row_index][col_index] = new_col[row_index];
               skip_count += 1;
           }
           else
           {
               ret[row_index][col_index] = matrix[row_index][col_index - skip_count];
           }
       }
   }
   return ret;
}

processing:Change,Paste,Replace

int ChangeRow(double[][] matrix, int row_index1, int row_index2)
{
   if (row_index1 == row_index2) { return 0; }
   double[] temp_row1 = GetRow(matrix, row_index1);
   double[] temp_row2 = GetRow(matrix, row_index2);
   SetRow(matrix, row_index1, temp_row2);
   SetRow(matrix, row_index2, temp_row1);
   return 1;
}

void _changeRow(double[][] matrix, int row_index1, int row_index2)
{
   int col_count = matrix[0].length;
   for (int col_index = 0; col_index < col_count; col_index++)
   {
       double temp = matrix[row_index1][col_index];
       matrix[row_index1][col_index] = matrix[row_index2][col_index];
       matrix[row_index2][col_index] = temp;
   }
}

int ChangeCol(double[][] matrix, int col_index1, int col_index2)
{
   if (col_index1 == col_index2) { return 0; }
   double[] temp_col1 = GetCol(matrix, col_index1);
   double[] temp_col2 = GetCol(matrix, col_index2);
   SetCol(matrix, col_index1, temp_col2);
   SetCol(matrix, col_index2, temp_col1);
   return 1;
}

void _changeCol(double[][] matrix, int col_index1, int col_index2)
{
   //colを固定してrowを走査する
   int row_count = matrix.length;
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       double temp = matrix[row_index][col_index1];
       matrix[row_index][col_index1] = matrix[row_index][col_index2];
       matrix[row_index][col_index2] = temp;
   }
}

void _changeComponent(double[] vector, int index1, int index2)
{
   double temp = vector[index1];
   vector[index1] = vector[index2];
   vector[index2] = temp;
} 

//newして戻り値返すタイプ
//昨日的にはSetと変わらない
double[][] PasteRow(double[][] matrix, int row_index, double[] new_row)
{
   int row_count = matrix.length;//rowの個数
   int col_count = matrix[0].length;//colの個数
   double[][] ret = new double[row_count][col_count];
   for (int col_index = 0; col_index < col_count; col_index++)
   {
       //入力配列が短い場合は0で埋める
       //必然的に、長い場合は切り落とされる
       if (new_row.length - 1 < col_index)
       {
           ret[row_index][col_index] = 0;
       }
       ret[row_index][col_index] = new_row[col_index];
   }
   return ret;
}

double[][] PasteCol(double[][] matrix, int col_index, double[] new_col)
{
   int row_count = matrix.length;//rowの個数
   int col_count = matrix[0].length;//colの個数
   double[][] ret = new double[row_count][col_count];
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       //入力配列が短い場合は0で埋める
       //必然的に、長い場合は切り落とされる
       if (new_col.length - 1 < row_index)
       {
           ret[row_index][col_index] = 0;
       }
       ret[row_index][col_index] = new_col[row_index];
   }
   return ret;
}

void Replace(double[] vector, int... indexes)
{
   for (int col_index = 0; col_index < vector.length; col_index++)
   {
       double temp = vector[col_index];
       vector[col_index] = vector[indexes[col_index]];
       vector[indexes[col_index]] = temp;
   }
}

//記録しておいた入れ換えを使用
void ReplaceRow(double[][] matrix, int... indexes)
{
   int row_count = matrix.length;
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       _changeRow(matrix, row_index, indexes[row_index]);
   }
}

//記録しておいた入れ換えを使用
void ReplaceCol(double[][] matrix, int... indexes)
{
   int col_count = matrix[0].length;
   for (int col_index = 0; col_index < col_count; col_index++)
   {
       _changeCol(matrix, col_index, indexes[col_index]);
   }
}

//今の所使わなくて良い
void ReplaceBackRow(double[][] matrix, int... indexes)
{
   int row_count = matrix.length;
   for (int row_index = row_count - 1; row_index >= 0; row_index--)
   {
       _changeRow(matrix, row_index, indexes[row_index]);
   }
}

void ReplaceBackCol(double[][] matrix, int... indexes)
{
   int col_count = matrix[0].length;
   for (int col_index = col_count - 1; col_index >= 0; col_index--)
   {
       _changeCol(matrix, col_index, indexes[col_index]);
   }
}

processing:行列の積



double[][] MatrixMultiple_VH(double[] vertical, double[] horizontal)
{
   int row_count = vertical.length;
   int col_count = horizontal.length;
   double[][] ret = new double[row_count][col_count];
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       for (int col_index = 0; col_index < col_count; col_index++)
       {
           ret[row_index][col_index] = vertical[row_index] * horizontal[col_index];
       }
   }
   return ret;
}

double[] MatrixMultiple(double[][] m, double[] v)
{
   int row_count = m.length;//row
   int col_count = m[0].length;//col
   double[] ret = new double[col_count];
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       //右に向かって走査
       for (int col_index = 0; col_index < col_count; col_index++)
       {
           ret[row_index] += m[row_index][col_index] * v[col_index];
       }
   }
   return ret;
}

double[] MatrixMultiple(double[] h, double[][] m)
{
   int row_count = m.length;//row
   int col_count = m[0].length;//col
   double[] ret = new double[col_count];
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       //右に向かって走査
       for (int col_index = 0; col_index < col_count; col_index++)
       {
           ret[row_index] += h[col_index] * m[row_index][col_index];
       }
   }
   return ret;
}

double[][] MatrixMultiple(double[][] m1, double[][] m2)
{
   int row_count = m1.length;//row
   int col_count = m1[0].length;//col
   double[][] ret = new double[row_count][col_count];
   for (int row_index = 0; row_index < row_count; row_index++)
   {
       //右に向かって走査
       for (int col_index = 0; col_index < col_count; col_index++)
       {
           double s = 0;
           for (int k = 0; k < col_count; k++)
           {
               s += m1[row_index][k] * m2[k][col_index];
           }

           ret[row_index][col_index] = s;
       }
   }
   return ret;
}

processing:Elementwise

processingはstatic使うとワケ分からんなるので使用せず



//配列の定数倍
double[] _calcuElementWise(double[] v1, double c, IFomula2 operation)
{
   double[] ret = new double[v1.length];

   for (int i = 0; i < v1.length; i++)
   {
       ret[i] = operation.exe(v1[i], c);
   }
   return ret;
}

//配列×配列
double[] _calcuElementWise(double[] v1, double[] v2, IFomula2 operation)
{
   if (v1.length != v2.length) { return null; }

   double[] ret = new double[v1.length];

   for (int i = 0; i < v1.length; i++)
   {
       ret[i] = operation.exe(v1[i], v2[i]);
   }
   return ret;
}

//行列の定数倍
double[][] _calcuElementWise(double[][] m1, double c, IFomula2 operation)
{
   double[][] ret = new double[m1.length][m1[0].length];

   for (int m = 0; m < m1.length; m++)
   {
       //右に向かって走査
       for (int n = 0; n < m1[0].length; n++)
       {
           ret[m][n] = operation.exe(m1[m][n], c);
       }
   }
   return ret;
}

//行列×縦ベクトル
double[] _calcuElementWise(double[][] m, double[] v, IFomula2 operation)
{
   double[] ret = new double[v.length];

   for (int row_index = 0; row_index < m.length; row_index++)
   {
       //右に向かって走査
       for (int col_index = 0; col_index < m[0].length; col_index++)
       {
           ret[row_index] = operation.exe(m[row_index][col_index], v[col_index]);
       }
   }
   return ret;
}

//行列×行列
double[][] _calcuElementWise(double[][] m1, double[][] m2, IFomula2 operation)
{
   if (m1.length != m2.length) { return null; }
   if (m1[0].length != m2[0].length) { return null; }

   double[][] ret = new double[m1.length][m1[0].length];

   //一行ずつ取り出す[同型なのでどちらの行列を基準にしても良い]
   for (int m = 0; m < m1.length; m++)
   {
       //右に向かって走査
       for (int n = 0; n < m1[0].length; n++)
       {
           ret[m][n] = operation.exe(m1[m][n], m2[m][n]);
       }
   }
   return ret;
}

//配列の定数倍
double[] Add(double[] v, double c){return _calcuElementWise(v, c, new Add());}
double[] Sub(double[] v, double c){return _calcuElementWise(v, c, new Sub());}
//これはElementWiseであって行列の積ではない
double[] Multiple(double c, double[] v){return _calcuElementWise(v, c, new Mult());}
double[] Multiple(double[] v, double c){return _calcuElementWise(v, c, new Mult());}
double[] Divide(double[] v, double c){return _calcuElementWise(v, c, new Div());}
double[] Mod(double[] v, double c){return _calcuElementWise(v, c, new Mod());}

//配列×配列
double[] Add(double[] v1, double[] v2){return _calcuElementWise(v1, v2, new Add());}
double[] Sub(double[] v1, double[] v2){return _calcuElementWise(v1, v2, new Sub());}
//これはElementWiseであって行列の積ではない
double[] Multiple(double[] v1, double[] v2){return _calcuElementWise(v1, v2, new Mult());}
double[] Divide(double[] v1, double[] v2){return _calcuElementWise(v1, v2, new Div());}
double[] Mod(double[] v1, double[] v2){return _calcuElementWise(v1, v2, new Mod());}

//行列の定数倍
double[][] Add(double[][] m1, double c){return _calcuElementWise(m1, c, new Add());}
double[][] Sub(double[][] m1, double c){return _calcuElementWise(m1, c, new Sub());}
//これはElementWiseであって行列の積ではない
double[][] Multiple(double c, double[][] m1){return _calcuElementWise(m1, c, new Mult());}
double[][] Multiple(double[][] m1, double c){return _calcuElementWise(m1, c, new Mult());}
double[][] Divide(double[][] m1, double c){return _calcuElementWise(m1, c, new Div());}
double[][] Mod(double[][] m1, double c){return _calcuElementWise(m1, c, new Mod());}

//横ベクトル×行列
double[] Add(double[] v, double[][] m){return _calcuElementWise(m, v, new Add());}
double[] Sub(double[] v, double[][] m){return _calcuElementWise(m, v, new Sub());}
//これはElementWiseであって行列の積ではない
double[] Multiple(double[] v, double[][] m){return _calcuElementWise(m, v, new Mult());}
double[] Divide(double[] v, double[][] m){return _calcuElementWise(m, v, new Div());}
double[] Mod(double[] v, double[][] m){return _calcuElementWise(m, v, new Mod());}
//行列×縦ベクトル
double[] Add(double[][] m, double[] v){return _calcuElementWise(m, v, new Add());}
double[] Sub(double[][] m, double[] v){return _calcuElementWise(m, v, new Sub());}
//これはElementWiseであって行列の積ではない
double[] Multiple(double[][] m, double[] v){return _calcuElementWise(m, v, new Mult());}
double[] Divide(double[][] m, double[] v){return _calcuElementWise(m, v, new Div());}
double[] Mod(double[][] m, double[] v){return _calcuElementWise(m, v, new Mod());}

//行列×行列
double[][] Add(double[][] m1, double[][] m2){return _calcuElementWise(m1, m2, new Add());}
double[][] Sub(double[][] m1, double[][] m2){return _calcuElementWise(m1, m2, new Sub());}
//これはElementWiseであって行列の積ではない
double[][] Multiple(double[][] m1, double[][] m2){return _calcuElementWise(m1, m2, new Mult());}
double[][] Divide(double[][] m1, double[][] m2){return _calcuElementWise(m1, m2, new Div());}
double[][] Mod(double[][] m1, double[][] m2){return _calcuElementWise(m1, m2, new Mod());}               

interface IFomula2{double exe(double a, double b);}
class Add implements IFomula2{double exe(double a, double b){return a+b;}}
class Sub implements IFomula2{double exe(double a, double b){return a-b;}}
class Mult implements IFomula2{double exe(double a, double b){return a*b;}}
class Div implements IFomula2{double exe(double a, double b){return a/b;}}
class Mod implements IFomula2{double exe(double a, double b){return a%b;}}
interface IFomula1{double exe(double a);}
class Exp implements IFomula1{double exe(double a){return Math.exp(a);}}
class Floor implements IFomula1{double exe(double a){return (int)a;}}
class Ceil implements IFomula1{double exe(double a){return (int)a+1;}}
class Round implements IFomula1
{
 double exe(double a)
 {
   Floor floor = new Floor();
   return floor.exe(a+0.5);
 }
}
class Fraction implements IFomula1{double exe(double a){return a-(int)a;}}

参考文献

行列あるいは線形代数なら

ガウスの消去法、あるいは反復法等、連立方程式を数値的に解くならば

ソースなり疑似コードから見るなら

最適化方面から

この記事が気に入ったらサポートをしてみませんか？