gauss_jordan_matrix_inversion.sf 999 B

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859
  1. #!/usr/bin/ruby
  2. #
  3. ## https://rosettacode.org/wiki/Gauss-Jordan_matrix_inversion
  4. #
  5. func rref (M) {
  6. var (j, rows, cols) = (0, M.len, M[0].len)
  7. for r in (^rows) {
  8. j < cols || return M
  9. var i = r
  10. while (!M[i][j]) {
  11. ++i == rows || next
  12. i = r
  13. ++j == cols && return M
  14. }
  15. M[i, r] = M[r, i] if (r != i)
  16. M[r] = (M[r] »/» M[r][j])
  17. for n in (^rows) {
  18. next if (n == r)
  19. M[n] = (M[n] »-« (M[r] »*» M[n][j]))
  20. }
  21. ++j
  22. }
  23. return M
  24. }
  25. func gauss_jordan_invert (M) {
  26. var I = M.len.of {|i|
  27. M.len.of {|j|
  28. i == j ? 1 : 0
  29. }
  30. }
  31. var A = gather {
  32. ^M -> each {|i| take(M[i] + I[i]) }
  33. }
  34. rref(A).map { .last(M.len) }
  35. }
  36. var n = irand(3, 7)
  37. var A = n.of {
  38. n.of { irand(-10, 10) }
  39. }
  40. say gauss_jordan_invert(A).map {
  41. .map { "%6s" % .as_rat }.join(" ")
  42. }.join("\n")
  43. assert_eq(A.inv, gauss_jordan_invert(A))