lib/matrix/matrix.nit

   1 # This file is part of NIT ( http://www.nitlanguage.org ).
   2 #
   3 # Licensed under the Apache License, Version 2.0 (the "License");
   4 # you may not use this file except in compliance with the License.
   5 # You may obtain a copy of the License at
   6 #
   7 #     http://www.apache.org/licenses/LICENSE-2.0
   8 #
   9 # Unless required by applicable law or agreed to in writing, software
  10 # distributed under the License is distributed on an "AS IS" BASIS,
  11 # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  12 # See the License for the specific language governing permissions and
  13 # limitations under the License.
  14
  15 # Services for matrices of `Float` values
  16 module matrix
  17
  18 # A rectangular array of `Float`
  19 #
  20 # Require: `width > 0 and height > 0`
  21 class Matrix
  22         super Cloneable
  23
  24         # Number of columns
  25         var width: Int
  26
  27         # Number of rows
  28         var height: Int
  29
  30         # Items of this matrix, rows by rows
  31         private var items = new NativeDoubleArray(width*height) is lazy
  32
  33         # Create a matrix from nested sequences
  34         #
  35         # Require: all rows are of the same length
  36         #
  37         # ~~~
  38         # var matrix = new Matrix.from([[1.0, 2.0],
  39         #                               [3.0, 4.0]])
  40         # assert matrix.to_s == """
  41         # 1.0 2.0
  42         # 3.0 4.0"""
  43         # ~~~
  44         init from(items: SequenceRead[SequenceRead[Float]])
  45         do
  46                 if items.is_empty then
  47                         init(0, 0)
  48                         return
  49                 end
  50
  51                 init(items.first.length, items.length)
  52
  53                 for j in height.times do assert items[j].length == width
  54
  55                 for j in [0..height[ do
  56                         for i in [0..width[ do
  57                                 self[j, i] = items[j][i]
  58                         end
  59                 end
  60         end
  61
  62         # Get each row of this matrix in nested arrays
  63         #
  64         # ~~~
  65         # var items = [[1.0, 2.0],
  66         #              [3.0, 4.0]]
  67         # var matrix = new Matrix.from(items)
  68         # assert matrix.to_a == items
  69         # ~~~
  70         fun to_a: Array[Array[Float]]
  71         do
  72                 var a = new Array[Array[Float]]
  73                 for j in height.times do
  74                         var row = new Array[Float]
  75                         for i in width.times do
  76                                 row.add self[j, i]
  77                         end
  78                         a.add row
  79                 end
  80                 return a
  81         end
  82
  83         # Create a matrix from an `Array[Float]` composed of rows after rows
  84         #
  85         # Require: `width > 0 and height > 0`
  86         # Require: `array.length >= width*height`
  87         #
  88         # ~~~
  89         # var matrix = new Matrix.from_array(2, 2, [1.0, 2.0,
  90         #                                           3.0, 4.0])
  91         # assert matrix.to_s == """
  92         # 1.0 2.0
  93         # 3.0 4.0"""
  94         # ~~~
  95         init from_array(width, height: Int, array: SequenceRead[Float])
  96         do
  97                 assert width > 0
  98                 assert height > 0
  99                 assert array.length >= width*height
 100
 101                 init(width, height)
 102
 103                 for i in height.times do
 104                         for j in width.times do
 105                                 self[j, i] = array[i + j*width]
 106                         end
 107                 end
 108         end
 109
 110         # Create an identity matrix
 111         #
 112         # Require: `size >= 0`
 113         #
 114         # ~~~
 115         # var i = new Matrix.identity(3)
 116         # assert i.to_s == """
 117         # 1.0 0.0 0.0
 118         # 0.0 1.0 0.0
 119         # 0.0 0.0 1.0"""
 120         # ~~~
 121         new identity(size: Int)
 122         do
 123                 assert size >= 0
 124
 125                 var matrix = new Matrix(size, size)
 126                 for i in [0..size[ do
 127                         for j in [0..size[ do
 128                                 matrix[j, i] = if i == j then 1.0 else 0.0
 129                         end
 130                 end
 131                 return matrix
 132         end
 133
 134         # Create a new clone of this matrix
 135         redef fun clone
 136         do
 137                 var c = new Matrix(width, height)
 138                 for i in [0..width*height[ do c.items[i] = items[i]
 139                 return c
 140         end
 141
 142         # Get the value at column `y` and row `x`
 143         #
 144         # Require: `x >= 0 and x <= width and y >= 0 and y <= height`
 145         #
 146         # ~~~
 147         # var matrix = new Matrix.from([[0.0, 0.1],
 148         #                               [1.0, 1.1]])
 149         #
 150         # assert matrix[0, 0] == 0.0
 151         # assert matrix[0, 1] == 0.1
 152         # assert matrix[1, 0] == 1.0
 153         # assert matrix[1, 1] == 1.1
 154         # ~~~
 155         fun [](y, x: Int): Float
 156         do
 157                 assert x >= 0 and x < width
 158                 assert y >= 0 and y < height
 159
 160                 return items[x + y*width]
 161         end
 162
 163         # Set the `value` at row `y` and column `x`
 164         #
 165         # Require: `x >= 0 and x <= width and y >= 0 and y <= height`
 166         #
 167         # ~~~
 168         # var matrix = new Matrix.identity(2)
 169         #
 170         # matrix[0, 0] = 0.0
 171         # matrix[0, 1] = 0.1
 172         # matrix[1, 0] = 1.0
 173         # matrix[1, 1] = 1.1
 174         #
 175         # assert matrix.to_s == """
 176         # 0.0 0.1
 177         # 1.0 1.1"""
 178         # ~~~
 179         fun []=(y, x: Int, value: Float)
 180         do
 181                 assert x >= 0 and x < width
 182                 assert y >= 0 and y < height
 183
 184                 items[x + y*width] = value
 185         end
 186
 187         # Matrix product (×)
 188         #
 189         # Require: `self.width == other.height`
 190         #
 191         # ~~~
 192         # var m = new Matrix.from([[3.0, 4.0],
 193         #                          [5.0, 6.0]])
 194         # var i = new Matrix.identity(2)
 195         #
 196         # assert m * i == m
 197         # assert (m * m).to_s == """
 198         # 29.0 36.0
 199         # 45.0 56.0"""
 200         #
 201         # var a = new Matrix.from([[1.0, 2.0, 3.0],
 202         #                          [4.0, 5.0, 6.0]])
 203         # var b = new Matrix.from([[1.0],
 204         #                          [2.0],
 205         #                          [3.0]])
 206         # var c = a * b
 207         # assert c.to_s == """
 208         # 14.0
 209         # 32.0"""
 210         # ~~~
 211         fun *(other: Matrix): Matrix
 212         do
 213                 assert self.width == other.height
 214
 215                 var out = new Matrix(other.width, self.height)
 216                 out.items.mul(items, other.items, self.width, self.height, other.width)
 217                 return out
 218         end
 219
 220         # Get the transpose of this matrix
 221         #
 222         # ~~~
 223         # var matrix = new Matrix.from([[1.0, 2.0, 3.0],
 224         #                               [4.0, 5.0, 6.0]])
 225         # assert matrix.transposed.to_a == [[1.0, 4.0],
 226         #                                   [2.0, 5.0],
 227         #                                   [3.0, 6.0]]
 228         #
 229         # var i = new Matrix.identity(3)
 230         # assert i.transposed == i
 231         # ~~~
 232         fun transposed: Matrix
 233         do
 234                 var out = new Matrix(height, width)
 235                 for k, v in self do out[k.x, k.y] = v
 236                 return out
 237         end
 238
 239         # Iterate over the values in this matrix
 240         fun iterator: MapIterator[MatrixCoordinate, Float] do return new MatrixIndexIterator(self)
 241
 242         redef fun to_s
 243         do
 244                 var s = new FlatBuffer
 245                 for y in [0..height[ do
 246                         for x in [0..width[ do
 247                                 s.append items[y*width+x].to_s
 248                                 if x < width-1 then s.add ' '
 249                         end
 250                         if y < height-1 then s.add '\n'
 251                 end
 252                 return s.to_s
 253         end
 254
 255         redef fun ==(other) do return other isa Matrix and
 256                 width == other.width and height == other.height and
 257                 items.equal_items(items, width*height)
 258
 259         redef fun hash do return items.hash_items(width*height)
 260 end
 261
 262 private class MatrixIndexIterator
 263         super MapIterator[MatrixCoordinate, Float]
 264
 265         var matrix: Matrix
 266
 267         redef var key = new MatrixCoordinate(0, 0)
 268
 269         redef fun is_ok do return key.y < matrix.height
 270
 271         redef fun item
 272         do
 273                 assert is_ok
 274                 return matrix[key.y, key.x]
 275         end
 276
 277         redef fun next
 278         do
 279                 assert is_ok
 280                 var key = key
 281                 if key.x == matrix.width - 1 then
 282                         key.x = 0
 283                         key.y += 1
 284                 else
 285                         key.x += 1
 286                 end
 287         end
 288 end
 289
 290 # Position key when iterating over the values of a matrix
 291 class MatrixCoordinate
 292         # Index of the current column
 293         var x: Int
 294
 295         # Index of the current row
 296         var y: Int
 297
 298         redef fun to_s do return "({x},{y})"
 299 end
 300
 301 # Specialized native structure to store matrix items and avoid boxing cost
 302 private extern class NativeDoubleArray `{ double* `}
 303
 304         new(size: Int) do
 305                 var sizeof_double = 8
 306                 var buf = new CString(sizeof_double*size)
 307                 return new NativeDoubleArray.in_buffer(buf)
 308         end
 309
 310         new in_buffer(buffer: CString) `{ return (double*)buffer; `}
 311
 312         fun [](i: Int): Float `{ return self[i]; `}
 313
 314         fun []=(i: Int, value: Float) `{ self[i] = value; `}
 315
 316         fun equal_items(other: NativeDoubleArray, len: Int): Bool `{
 317                 int i;
 318                 for (i = 0; i < len; i ++)
 319                         if (self[i] != other[i])
 320                                 return 0;
 321                 return 1;
 322         `}
 323
 324         fun hash_items(len: Int): Int `{
 325                 // Adapted from `SequenceRead::hash`
 326                 long r = 17+len;
 327                 int i;
 328                 for (i = 0; i < len; i ++)
 329                         r = r * 3 / 2 + (long)(i*1024.0);
 330                 return r;
 331         `}
 332
 333         fun mul(a, b: NativeDoubleArray, a_width, a_height, b_width: Int) `{
 334                 int i, j, k;
 335                 for (j = 0; j < a_height; j ++)
 336                         for (i = 0; i < b_width; i ++) {
 337                                 float sum = 0.0;
 338                                 for (k = 0; k < a_width; k ++) sum += a[j*a_width + k] * b[k*b_width + i];
 339                                 self[j*b_width + i] = sum;
 340                         }
 341         `}
 342 end