lusolve.js 3.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113
  1. "use strict";
  2. Object.defineProperty(exports, "__esModule", {
  3. value: true
  4. });
  5. exports.createLusolve = void 0;
  6. var _is = require("../../../utils/is.js");
  7. var _factory = require("../../../utils/factory.js");
  8. var _solveValidation = require("./utils/solveValidation.js");
  9. var _csIpvec = require("../sparse/csIpvec.js");
  10. var name = 'lusolve';
  11. var dependencies = ['typed', 'matrix', 'lup', 'slu', 'usolve', 'lsolve', 'DenseMatrix'];
  12. var createLusolve = /* #__PURE__ */(0, _factory.factory)(name, dependencies, function (_ref) {
  13. var typed = _ref.typed,
  14. matrix = _ref.matrix,
  15. lup = _ref.lup,
  16. slu = _ref.slu,
  17. usolve = _ref.usolve,
  18. lsolve = _ref.lsolve,
  19. DenseMatrix = _ref.DenseMatrix;
  20. var solveValidation = (0, _solveValidation.createSolveValidation)({
  21. DenseMatrix: DenseMatrix
  22. });
  23. /**
  24. * Solves the linear system `A * x = b` where `A` is an [n x n] matrix and `b` is a [n] column vector.
  25. *
  26. * Syntax:
  27. *
  28. * math.lusolve(A, b) // returns column vector with the solution to the linear system A * x = b
  29. * math.lusolve(lup, b) // returns column vector with the solution to the linear system A * x = b, lup = math.lup(A)
  30. *
  31. * Examples:
  32. *
  33. * const m = [[1, 0, 0, 0], [0, 2, 0, 0], [0, 0, 3, 0], [0, 0, 0, 4]]
  34. *
  35. * const x = math.lusolve(m, [-1, -1, -1, -1]) // x = [[-1], [-0.5], [-1/3], [-0.25]]
  36. *
  37. * const f = math.lup(m)
  38. * const x1 = math.lusolve(f, [-1, -1, -1, -1]) // x1 = [[-1], [-0.5], [-1/3], [-0.25]]
  39. * const x2 = math.lusolve(f, [1, 2, 1, -1]) // x2 = [[1], [1], [1/3], [-0.25]]
  40. *
  41. * const a = [[-2, 3], [2, 1]]
  42. * const b = [11, 9]
  43. * const x = math.lusolve(a, b) // [[2], [5]]
  44. *
  45. * See also:
  46. *
  47. * lup, slu, lsolve, usolve
  48. *
  49. * @param {Matrix | Array | Object} A Invertible Matrix or the Matrix LU decomposition
  50. * @param {Matrix | Array} b Column Vector
  51. * @param {number} [order] The Symbolic Ordering and Analysis order, see slu for details. Matrix must be a SparseMatrix
  52. * @param {Number} [threshold] Partial pivoting threshold (1 for partial pivoting), see slu for details. Matrix must be a SparseMatrix.
  53. *
  54. * @return {DenseMatrix | Array} Column vector with the solution to the linear system A * x = b
  55. */
  56. return typed(name, {
  57. 'Array, Array | Matrix': function ArrayArrayMatrix(a, b) {
  58. a = matrix(a);
  59. var d = lup(a);
  60. var x = _lusolve(d.L, d.U, d.p, null, b);
  61. return x.valueOf();
  62. },
  63. 'DenseMatrix, Array | Matrix': function DenseMatrixArrayMatrix(a, b) {
  64. var d = lup(a);
  65. return _lusolve(d.L, d.U, d.p, null, b);
  66. },
  67. 'SparseMatrix, Array | Matrix': function SparseMatrixArrayMatrix(a, b) {
  68. var d = lup(a);
  69. return _lusolve(d.L, d.U, d.p, null, b);
  70. },
  71. 'SparseMatrix, Array | Matrix, number, number': function SparseMatrixArrayMatrixNumberNumber(a, b, order, threshold) {
  72. var d = slu(a, order, threshold);
  73. return _lusolve(d.L, d.U, d.p, d.q, b);
  74. },
  75. 'Object, Array | Matrix': function ObjectArrayMatrix(d, b) {
  76. return _lusolve(d.L, d.U, d.p, d.q, b);
  77. }
  78. });
  79. function _toMatrix(a) {
  80. if ((0, _is.isMatrix)(a)) {
  81. return a;
  82. }
  83. if ((0, _is.isArray)(a)) {
  84. return matrix(a);
  85. }
  86. throw new TypeError('Invalid Matrix LU decomposition');
  87. }
  88. function _lusolve(l, u, p, q, b) {
  89. // verify decomposition
  90. l = _toMatrix(l);
  91. u = _toMatrix(u);
  92. // apply row permutations if needed (b is a DenseMatrix)
  93. if (p) {
  94. b = solveValidation(l, b, true);
  95. b._data = (0, _csIpvec.csIpvec)(p, b._data);
  96. }
  97. // use forward substitution to resolve L * y = b
  98. var y = lsolve(l, b);
  99. // use backward substitution to resolve U * x = y
  100. var x = usolve(u, y);
  101. // apply column permutations if needed (x is a DenseMatrix)
  102. if (q) {
  103. x._data = (0, _csIpvec.csIpvec)(q, x._data);
  104. }
  105. return x;
  106. }
  107. });
  108. exports.createLusolve = createLusolve;