123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161 |
- import { factory } from '../../../utils/factory.js';
- import { createSolveValidation } from './utils/solveValidation.js';
- var name = 'usolve';
- var dependencies = ['typed', 'matrix', 'divideScalar', 'multiplyScalar', 'subtract', 'equalScalar', 'DenseMatrix'];
- export var createUsolve = /* #__PURE__ */factory(name, dependencies, _ref => {
- var {
- typed,
- matrix,
- divideScalar,
- multiplyScalar,
- subtract,
- equalScalar,
- DenseMatrix
- } = _ref;
- var solveValidation = createSolveValidation({
- DenseMatrix
- });
- /**
- * Finds one solution of a linear equation system by backward substitution. Matrix must be an upper triangular matrix. Throws an error if there's no solution.
- *
- * `U * x = b`
- *
- * Syntax:
- *
- * math.usolve(U, b)
- *
- * Examples:
- *
- * const a = [[-2, 3], [2, 1]]
- * const b = [11, 9]
- * const x = usolve(a, b) // [[8], [9]]
- *
- * See also:
- *
- * usolveAll, lup, slu, usolve, lusolve
- *
- * @param {Matrix, Array} U A N x N matrix or array (U)
- * @param {Matrix, Array} b A column vector with the b values
- *
- * @return {DenseMatrix | Array} A column vector with the linear system solution (x)
- */
- return typed(name, {
- 'SparseMatrix, Array | Matrix': function SparseMatrixArrayMatrix(m, b) {
- return _sparseBackwardSubstitution(m, b);
- },
- 'DenseMatrix, Array | Matrix': function DenseMatrixArrayMatrix(m, b) {
- return _denseBackwardSubstitution(m, b);
- },
- 'Array, Array | Matrix': function ArrayArrayMatrix(a, b) {
- var m = matrix(a);
- var r = _denseBackwardSubstitution(m, b);
- return r.valueOf();
- }
- });
- function _denseBackwardSubstitution(m, b) {
- // make b into a column vector
- b = solveValidation(m, b, true);
- var bdata = b._data;
- var rows = m._size[0];
- var columns = m._size[1];
- // result
- var x = [];
- var mdata = m._data;
- // loop columns backwards
- for (var j = columns - 1; j >= 0; j--) {
- // b[j]
- var bj = bdata[j][0] || 0;
- // x[j]
- var xj = void 0;
- if (!equalScalar(bj, 0)) {
- // value at [j, j]
- var vjj = mdata[j][j];
- if (equalScalar(vjj, 0)) {
- // system cannot be solved
- throw new Error('Linear system cannot be solved since matrix is singular');
- }
- xj = divideScalar(bj, vjj);
- // loop rows
- for (var i = j - 1; i >= 0; i--) {
- // update copy of b
- bdata[i] = [subtract(bdata[i][0] || 0, multiplyScalar(xj, mdata[i][j]))];
- }
- } else {
- // zero value at j
- xj = 0;
- }
- // update x
- x[j] = [xj];
- }
- return new DenseMatrix({
- data: x,
- size: [rows, 1]
- });
- }
- function _sparseBackwardSubstitution(m, b) {
- // make b into a column vector
- b = solveValidation(m, b, true);
- var bdata = b._data;
- var rows = m._size[0];
- var columns = m._size[1];
- var values = m._values;
- var index = m._index;
- var ptr = m._ptr;
- // result
- var x = [];
- // loop columns backwards
- for (var j = columns - 1; j >= 0; j--) {
- var bj = bdata[j][0] || 0;
- if (!equalScalar(bj, 0)) {
- // non-degenerate row, find solution
- var vjj = 0;
- // upper triangular matrix values & index (column j)
- var jValues = [];
- var jIndices = [];
- // first & last indeces in column
- var firstIndex = ptr[j];
- var lastIndex = ptr[j + 1];
- // values in column, find value at [j, j], loop backwards
- for (var k = lastIndex - 1; k >= firstIndex; k--) {
- var i = index[k];
- // check row (rows are not sorted!)
- if (i === j) {
- vjj = values[k];
- } else if (i < j) {
- // store upper triangular
- jValues.push(values[k]);
- jIndices.push(i);
- }
- }
- // at this point we must have a value in vjj
- if (equalScalar(vjj, 0)) {
- throw new Error('Linear system cannot be solved since matrix is singular');
- }
- var xj = divideScalar(bj, vjj);
- for (var _k = 0, _lastIndex = jIndices.length; _k < _lastIndex; _k++) {
- var _i = jIndices[_k];
- bdata[_i] = [subtract(bdata[_i][0], multiplyScalar(xj, jValues[_k]))];
- }
- x[j] = [xj];
- } else {
- // degenerate row, we can choose any value
- x[j] = [0];
- }
- }
- return new DenseMatrix({
- data: x,
- size: [rows, 1]
- });
- }
- });
|