1)abbrev domain IMATRIX IndexedMatrix 2++ Author: Grabmeier, Gschnitzer, Williamson 3++ Date Created: 1987 4++ Basic Operations: 5++ Related Domains: Matrix, RectangularMatrix, SquareMatrix, 6++ StorageEfficientMatrixOperations 7++ Also See: 8++ AMS Classifications: 9++ Keywords: matrix, linear algebra 10++ Examples: 11++ References: 12++ Description: 13++ An \spad{IndexedMatrix} is a matrix where the minimal row and column 14++ indices are parameters of the type. The domains Row and Col 15++ are both IndexedVectors. 16++ The index of the 'first' row may be obtained by calling the 17++ function \spadfun{minRowIndex}. The index of the 'first' column may 18++ be obtained by calling the function \spadfun{minColIndex}. The index of 19++ the first element of a 'Row' is the same as the index of the 20++ first column in a matrix and vice versa. 21IndexedMatrix(R, mnRow, mnCol) : Exports == Implementation where 22 R : AbelianMonoid 23 mnRow, mnCol : Integer 24 Row ==> IndexedVector(R, mnCol) 25 Col ==> IndexedVector(R, mnRow) 26 MATLIN ==> MatrixLinearAlgebraFunctions(R, Row, Col, %) 27 28 Exports ==> MatrixCategory(R, Row, Col) 29 30 Implementation ==> 31 InnerIndexedTwoDimensionalArray(R, mnRow, mnCol, Row, Col) add 32 33 Qelt2 ==> QAREF2O$Lisp 34 Qsetelt2 ==> QSETAREF2O$Lisp 35 36 swapRows!(x, i1, i2) == 37 (i1 < minRowIndex(x)) or (i1 > maxRowIndex(x)) or _ 38 (i2 < minRowIndex(x)) or (i2 > maxRowIndex(x)) => 39 error "swapRows!: index out of range" 40 i1 = i2 => x 41 ro := mnRow 42 co := mnCol 43 for j in co..maxColIndex(x) repeat 44 t1 : R := Qelt2(x, i1, j, ro, co) 45 t2 : R := Qelt2(x, i2, j, ro, co) 46 Qsetelt2(x, i1, j, t2, ro, co) 47 Qsetelt2(x, i2, j, t1, ro, co) 48 x 49 50 if R has CommutativeRing then 51 52 determinant x == determinant(x)$MATLIN 53 minordet x == minordet(x)$MATLIN 54 55 if R has EuclideanDomain then 56 57 rowEchelon x == rowEchelon(x)$MATLIN 58 59 if R has IntegralDomain then 60 61 rank x == rank(x)$MATLIN 62 nullity x == nullity(x)$MATLIN 63 nullSpace x == nullSpace(x)$MATLIN 64 65 if R has Field then 66 67 inverse x == inverse(x)$MATLIN 68 69)abbrev domain MATRIX Matrix 70++ Author: Grabmeier, Gschnitzer, Williamson 71++ Date Created: 1987 72++ Basic Operations: 73++ Related Domains: IndexedMatrix, RectangularMatrix, SquareMatrix 74++ Also See: 75++ AMS Classifications: 76++ Keywords: matrix, linear algebra 77++ Examples: 78++ References: 79++ Description: 80++ \spadtype{Matrix} is a matrix domain where 1-based indexing is used 81++ for both rows and columns. 82Matrix(R) : Exports == Implementation where 83 R : AbelianMonoid 84 Row ==> Vector R 85 Col ==> Vector R 86 mnRow ==> 1 87 mnCol ==> 1 88 MATLIN ==> MatrixLinearAlgebraFunctions(R, Row, Col, %) 89 90 Exports ==> MatrixCategory(R, Row, Col) with 91 diagonalMatrix : Vector R -> % 92 ++ \spad{diagonalMatrix(v)} returns a diagonal matrix where the elements 93 ++ of v appear on the diagonal. 94 95 if R has ConvertibleTo InputForm then ConvertibleTo InputForm 96 97 if R has IntegralDomain then 98 invertIfCan : % -> Union(%,"failed") 99 ++ \spad{invertIfCan(m)} returns the inverse of the matrix m. 100 ++ If the matrix is not invertible, "failed" is returned. 101 ++ Error: if the matrix is not square. 102-- matrix: Vector Vector R -> % 103-- ++ \spad{matrix(v)} converts the vector of vectors v to a matrix, where 104-- ++ the vector of vectors is viewed as a vector of the rows of the 105-- ++ matrix 106-- diagonalMatrix: Vector % -> % 107-- ++ \spad{diagonalMatrix([m1, ..., mk])} creates a block diagonal matrix 108-- ++ M with block matrices {\em m1}, ..., {\em mk} down the diagonal, 109-- ++ with 0 block matrices elsewhere. 110-- vectorOfVectors: % -> Vector Vector R 111-- ++ \spad{vectorOfVectors(m)} returns the rows of the matrix m as a 112-- ++ vector of vectors 113 114 Implementation ==> 115 InnerIndexedTwoDimensionalArray(R, mnRow, mnCol, Row, Col) add 116 minr ==> minRowIndex 117 maxr ==> maxRowIndex 118 minc ==> minColIndex 119 maxc ==> maxColIndex 120 mini ==> minIndex 121 maxi ==> maxIndex 122 Qelt2 ==> QAREF2O$Lisp 123 Qsetelt2 ==> QSETAREF2O$Lisp 124 125 import from List(List(R)) 126 127 minRowIndex x == mnRow 128 minColIndex x == mnCol 129 130 -- qelt, qsetelt!, swapRows! and copy are logically unnecessary, 131 -- but good for performance 132 133 qelt(m, i, j) == Qelt2(m, i, j, mnRow@Integer, mnCol@Integer) 134 qsetelt!(m, i, j, r) == Qsetelt2(m, i, j, r, mnRow@Integer, mnCol@Integer) 135 136 swapRows!(x, i1, i2) == 137 (i1 < minRowIndex(x)) or (i1 > maxRowIndex(x)) or _ 138 (i2 < minRowIndex(x)) or (i2 > maxRowIndex(x)) => 139 error "swapRows!: index out of range" 140 i1 = i2 => x 141 for j in minColIndex(x)..maxColIndex(x) repeat 142 t1 : R := Qelt2(x, i1, j, mnRow@Integer, mnCol@Integer) 143 t2 : R := Qelt2(x, i2, j, mnRow@Integer, mnCol@Integer) 144 Qsetelt2(x, i1, j, t2, mnRow@Integer, mnCol@Integer) 145 Qsetelt2(x, i2, j, t1, mnRow@Integer, mnCol@Integer) 146 x 147 148 copy m == 149 ans : % := MAKE_MATRIX(nrows m, ncols m)$Lisp 150 for i in minRowIndex(m)..maxRowIndex(m) repeat 151 for j in minColIndex(m)..maxColIndex(m) repeat 152 qsetelt!(ans, i, j, qelt(m, i, j)) 153 ans 154 155 if R has CommutativeRing then 156 157 determinant x == determinant(x)$MATLIN 158 minordet x == minordet(x)$MATLIN 159 160 if R has EuclideanDomain then 161 162 rowEchelon x == rowEchelon(x)$MATLIN 163 164 if R has IntegralDomain then 165 166 rank x == rank(x)$MATLIN 167 nullity x == nullity(x)$MATLIN 168 nullSpace x == nullSpace(x)$MATLIN 169 170 if R has Field then 171 172 inverse x == inverse(x)$MATLIN 173 174 if R has IntegralDomain then 175 176 invertIfCan(x) == invertIfCan(x)$MATLIN 177 178-- matrix(v: Vector Vector R) == 179-- (rows := # v) = 0 => new(0, 0, 0) 180-- -- error check: this is a top level function 181-- cols := # v.mini(v) 182-- for k in (mini(v) + 1)..maxi(v) repeat 183-- cols ~= # v.k => error "matrix: rows of different lengths" 184-- ans := new(rows, cols, 0) 185-- for i in minr(ans)..maxr(ans) for k in mini(v)..maxi(v) repeat 186-- vv := v.k 187-- for j in minc(ans)..maxc(ans) for l in mini(vv)..maxi(vv) repeat 188-- ans(i, j) := vv.l 189-- ans 190 191 diagonalMatrix(v : Vector R) == 192 n := #v; ans := zero(n, n) 193 for i in minr(ans)..maxr(ans) for j in minc(ans)..maxc(ans) _ 194 for k in mini(v)..maxi(v) repeat qsetelt!(ans, i, j, qelt(v, k)) 195 ans 196 197-- diagonalMatrix(vec: Vector %) == 198-- rows : NonNegativeInteger := 0 199-- cols : NonNegativeInteger := 0 200-- for r in mini(vec)..maxi(vec) repeat 201-- mat := vec.r 202-- rows := rows + nrows mat; cols := cols + ncols mat 203-- ans := zero(rows, cols) 204-- loR := minr ans; loC := minc ans 205-- for r in mini(vec)..maxi(vec) repeat 206-- mat := vec.r 207-- hiR := loR + nrows(mat) - 1; hiC := loC + nrows(mat) - 1 208-- for i in loR..hiR for k in minr(mat)..maxr(mat) repeat 209-- for j in loC..hiC for l in minc(mat)..maxc(mat) repeat 210-- ans(i, j) := mat(k, l) 211-- loR := hiR + 1; loC := hiC + 1 212-- ans 213 214-- vectorOfVectors x == 215-- vv : Vector Vector R := new(nrows x, 0) 216-- cols := ncols x 217-- for k in mini(vv)..maxi(vv) repeat 218-- vv.k := new(cols, 0) 219-- for i in minr(x)..maxr(x) for k in mini(vv)..maxi(vv) repeat 220-- v := vv.k 221-- for j in minc(x)..maxc(x) for l in mini(v)..maxi(v) repeat 222-- v.l := x(i, j) 223-- vv 224 225 if R has ConvertibleTo InputForm then 226 convert(x : %) : InputForm == 227 convert [convert('matrix)@InputForm, 228 convert listOfLists x]$List(InputForm) 229 230)abbrev domain RMATRIX RectangularMatrix 231++ Author: Grabmeier, Gschnitzer, Williamson 232++ Date Created: 1987 233++ Basic Operations: 234++ Related Domains: IndexedMatrix, Matrix, SquareMatrix 235++ Also See: 236++ AMS Classifications: 237++ Keywords: matrix, linear algebra 238++ Examples: 239++ References: 240++ Description: 241++ \spadtype{RectangularMatrix} is a matrix domain where the number of rows 242++ and the number of columns are parameters of the domain. 243RectangularMatrix(m, n, R) : Exports == Implementation where 244 m, n : NonNegativeInteger 245 R : Join(SemiRng, AbelianMonoid) 246 Row ==> DirectProduct(n, R) 247 Col ==> DirectProduct(m, R) 248 Exports ==> Join(RectangularMatrixCategory(m, n, R, Row, Col), _ 249 CoercibleTo Matrix R) with 250 251 if R has ConvertibleTo InputForm then ConvertibleTo InputForm 252 253 rectangularMatrix : Matrix R -> % 254 ++ \spad{rectangularMatrix(m)} converts a matrix of type \spadtype{Matrix} 255 ++ to a matrix of type \spad{RectangularMatrix}. 256 coerce : % -> Matrix R 257 ++ \spad{coerce(m)} converts a matrix of type \spadtype{RectangularMatrix} 258 ++ to a matrix of type \spad{Matrix}. 259 260 Implementation ==> Matrix R add 261 minr ==> minRowIndex 262 maxr ==> maxRowIndex 263 minc ==> minColIndex 264 maxc ==> maxColIndex 265 mini ==> minIndex 266 maxi ==> maxIndex 267 268 ZERO := new(m, n, 0)$Matrix(R) pretend % 269 0 == ZERO 270 271 coerce(x : %) : OutputForm == coerce(x pretend Matrix R)$Matrix(R) 272 273 matrix(l : List List R) == 274 -- error check: this is a top level function 275 #l ~= m => error "matrix: wrong number of rows" 276 for ll in l repeat 277 #ll ~= n => error "matrix: wrong number of columns" 278 ans : Matrix R := new(m, n, 0) 279 for i in minr(ans)..maxr(ans) for ll in l repeat 280 for j in minc(ans)..maxc(ans) for r in ll repeat 281 qsetelt!(ans, i, j, r) 282 ans pretend % 283 284 row(x, i) == directProduct row(x pretend Matrix(R), i) 285 column(x, j) == directProduct column(x pretend Matrix(R), j) 286 287 coerce(x : %) : Matrix(R) == copy(x pretend Matrix(R)) 288 289 rectangularMatrix x == 290 (nrows(x) ~= m) or (ncols(x) ~= n) => 291 error "rectangularMatrix: matrix of bad dimensions" 292 copy(x) pretend % 293 294 if R has EuclideanDomain then 295 296 rowEchelon x == rowEchelon(x pretend Matrix(R)) pretend % 297 298 columnSpace x == 299 [directProduct c for c in columnSpace(x pretend Matrix R)] 300 301 if R has IntegralDomain then 302 303 rank x == rank(x pretend Matrix(R)) 304 nullity x == nullity(x pretend Matrix(R)) 305 nullSpace x == 306 [directProduct c for c in nullSpace(x pretend Matrix(R))] 307 308 if R has ConvertibleTo InputForm then 309 convert(x : %) : InputForm == 310 convert [convert('rectangularMatrix)@InputForm, 311 convert(x::Matrix(R))]$List(InputForm) 312 313)abbrev domain SQMATRIX SquareMatrix 314++ Author: Grabmeier, Gschnitzer, Williamson 315++ Date Created: 1987 316++ Basic Operations: 317++ Related Domains: IndexedMatrix, Matrix, RectangularMatrix 318++ Also See: 319++ AMS Classifications: 320++ Keywords: matrix, linear algebra 321++ Examples: 322++ References: 323++ Description: 324++ \spadtype{SquareMatrix} is a matrix domain of square matrices, where the 325++ number of rows (= number of columns) is a parameter of the type. 326SquareMatrix(ndim, R) : Exports == Implementation where 327 ndim : NonNegativeInteger 328 R : Join(SemiRng, AbelianMonoid) 329 Row ==> DirectProduct(ndim, R) 330 Col ==> DirectProduct(ndim, R) 331 MATLIN ==> MatrixLinearAlgebraFunctions(R, Row, Col, %) 332 333 Exports ==> Join(SquareMatrixCategory(ndim, R, Row, Col), _ 334 CoercibleTo Matrix R) with 335 336 transpose : % -> % 337 ++ \spad{transpose(m)} returns the transpose of the matrix m. 338 squareMatrix : Matrix R -> % 339 ++ \spad{squareMatrix(m)} converts a matrix of type \spadtype{Matrix} 340 ++ to a matrix of type \spadtype{SquareMatrix}. 341 coerce : % -> Matrix R 342 ++ \spad{coerce(m)} converts a matrix of type \spadtype{SquareMatrix} 343 ++ to a matrix of type \spadtype{Matrix}. 344-- symdecomp : % -> Record(sym: %, antisym: %) 345-- ++ \spad{symdecomp(m)} decomposes the matrix m as a sum of a symmetric 346-- ++ matrix \spad{m1} and an antisymmetric matrix \spad{m2}. The object 347-- ++ returned is the Record \spad{[m1, m2]} 348-- if R has CommutativeStar then 349-- minorsVect: -> Vector(Union(R,"uncomputed")) --range: 1..2^n-1 350-- ++ \spad{minorsVect(m)} returns a vector of the minors of the matrix m 351 if R has CommutativeStar and R has unitsKnown then unitsKnown 352 ++ the invertible matrices are simply the matrices whose determinants 353 ++ are units in the Ring R. 354 if R has ConvertibleTo InputForm then ConvertibleTo InputForm 355 356 Implementation ==> Matrix R add 357 minr ==> minRowIndex 358 maxr ==> maxRowIndex 359 minc ==> minColIndex 360 maxc ==> maxColIndex 361 mini ==> minIndex 362 maxi ==> maxIndex 363 364 ZERO := scalarMatrix 0 365 0 == ZERO 366 367 if R has Monoid then 368 369 ONE := scalarMatrix 1 370 1 == ONE 371 372 if R has Ring then 373 374 characteristic() == characteristic()$R 375 376 matrix(l : List List R) == 377 -- error check: this is a top level function 378 #l ~= ndim => error "matrix: wrong number of rows" 379 for ll in l repeat 380 #ll ~= ndim => error "matrix: wrong number of columns" 381 ans : Matrix R := new(ndim, ndim, 0) 382 for i in minr(ans)..maxr(ans) for ll in l repeat 383 for j in minc(ans)..maxc(ans) for r in ll repeat 384 qsetelt!(ans, i, j, r) 385 ans pretend % 386 387 row(x, i) == directProduct row(x pretend Matrix(R), i) 388 column(x, j) == directProduct column(x pretend Matrix(R), j) 389 coerce(x : %) : OutputForm == coerce(x pretend Matrix R)$Matrix(R) 390 391 scalarMatrix r == scalarMatrix(ndim, r)$Matrix(R) pretend % 392 393 diagonalMatrix l == 394 #l ~= ndim => 395 error "diagonalMatrix: wrong number of entries in list" 396 diagonalMatrix(l)$Matrix(R) pretend % 397 398 coerce(x : %) : Matrix(R) == copy(x pretend Matrix(R)) 399 400 squareMatrix x == 401 (nrows(x) ~= ndim) or (ncols(x) ~= ndim) => 402 error "squareMatrix: matrix of bad dimensions" 403 copy(x) pretend % 404 405 x : % * v : Col == 406 directProduct((x pretend Matrix(R)) * (v :: Vector(R))) 407 408 v : Row * x : % == 409 directProduct((v :: Vector(R)) * (x pretend Matrix(R))) 410 411 if R has CommutativeRing then 412 413 determinant x == determinant(x pretend Matrix(R)) 414 minordet x == minordet(x pretend Matrix(R)) 415 Pfaffian x == Pfaffian(x pretend Matrix(R)) 416 417 if R has EuclideanDomain then 418 419 rowEchelon x == rowEchelon(x pretend Matrix(R)) pretend % 420 421 columnSpace x == 422 [directProduct c for c in columnSpace(x pretend Matrix R)] 423 424 if R has IntegralDomain then 425 426 rank x == rank(x pretend Matrix(R)) 427 nullity x == nullity(x pretend Matrix(R)) 428 nullSpace x == 429 [directProduct c for c in nullSpace(x pretend Matrix(R))] 430 recip(x) == 431 (u := invertIfCan(x pretend Matrix(R))) case "failed" => "failed" 432 (u :: Matrix(R)) pretend % 433 434 if R has Field then 435 436 -- dimension() == (m * n) :: CardinalNumber 437 438 inverse x == 439 (u := inverse(x pretend Matrix(R))) case "failed" => "failed" 440 (u :: Matrix(R)) pretend % 441 442 x : % ^ n : Integer == 443 ((x pretend Matrix(R)) ^ n) pretend % 444 445 recip x == inverse x 446 447 if R has ConvertibleTo InputForm then 448 convert(x : %) : InputForm == 449 convert [convert('squareMatrix)@InputForm, 450 convert(x::Matrix(R))]$List(InputForm) 451 452 453--Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd. 454--All rights reserved. 455-- 456--Redistribution and use in source and binary forms, with or without 457--modification, are permitted provided that the following conditions are 458--met: 459-- 460-- - Redistributions of source code must retain the above copyright 461-- notice, this list of conditions and the following disclaimer. 462-- 463-- - Redistributions in binary form must reproduce the above copyright 464-- notice, this list of conditions and the following disclaimer in 465-- the documentation and/or other materials provided with the 466-- distribution. 467-- 468-- - Neither the name of The Numerical ALgorithms Group Ltd. nor the 469-- names of its contributors may be used to endorse or promote products 470-- derived from this software without specific prior written permission. 471-- 472--THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS 473--IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED 474--TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A 475--PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER 476--OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, 477--EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 478--PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 479--PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF 480--LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING 481--NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 482--SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 483