simple-squiggle

A restricted subset of Squiggle
Log | Files | Refs | README

det.js (4122B)

      1 "use strict";
      2 
      3 Object.defineProperty(exports, "__esModule", {
      4   value: true
      5 });
      6 exports.createDet = void 0;
      7 
      8 var _is = require("../../utils/is.js");
      9 
     10 var _object = require("../../utils/object.js");
     11 
     12 var _string = require("../../utils/string.js");
     13 
     14 var _factory = require("../../utils/factory.js");
     15 
     16 var name = 'det';
     17 var dependencies = ['typed', 'matrix', 'subtract', 'multiply', 'unaryMinus', 'lup'];
     18 var createDet = /* #__PURE__ */(0, _factory.factory)(name, dependencies, function (_ref) {
     19   var typed = _ref.typed,
     20       matrix = _ref.matrix,
     21       subtract = _ref.subtract,
     22       multiply = _ref.multiply,
     23       unaryMinus = _ref.unaryMinus,
     24       lup = _ref.lup;
     25 
     26   /**
     27    * Calculate the determinant of a matrix.
     28    *
     29    * Syntax:
     30    *
     31    *    math.det(x)
     32    *
     33    * Examples:
     34    *
     35    *    math.det([[1, 2], [3, 4]]) // returns -2
     36    *
     37    *    const A = [
     38    *      [-2, 2, 3],
     39    *      [-1, 1, 3],
     40    *      [2, 0, -1]
     41    *    ]
     42    *    math.det(A) // returns 6
     43    *
     44    * See also:
     45    *
     46    *    inv
     47    *
     48    * @param {Array | Matrix} x  A matrix
     49    * @return {number} The determinant of `x`
     50    */
     51   return typed(name, {
     52     any: function any(x) {
     53       return (0, _object.clone)(x);
     54     },
     55     'Array | Matrix': function det(x) {
     56       var size;
     57 
     58       if ((0, _is.isMatrix)(x)) {
     59         size = x.size();
     60       } else if (Array.isArray(x)) {
     61         x = matrix(x);
     62         size = x.size();
     63       } else {
     64         // a scalar
     65         size = [];
     66       }
     67 
     68       switch (size.length) {
     69         case 0:
     70           // scalar
     71           return (0, _object.clone)(x);
     72 
     73         case 1:
     74           // vector
     75           if (size[0] === 1) {
     76             return (0, _object.clone)(x.valueOf()[0]);
     77           } else {
     78             throw new RangeError('Matrix must be square ' + '(size: ' + (0, _string.format)(size) + ')');
     79           }
     80 
     81         case 2:
     82           {
     83             // two dimensional array
     84             var rows = size[0];
     85             var cols = size[1];
     86 
     87             if (rows === cols) {
     88               return _det(x.clone().valueOf(), rows, cols);
     89             } else {
     90               throw new RangeError('Matrix must be square ' + '(size: ' + (0, _string.format)(size) + ')');
     91             }
     92           }
     93 
     94         default:
     95           // multi dimensional array
     96           throw new RangeError('Matrix must be two dimensional ' + '(size: ' + (0, _string.format)(size) + ')');
     97       }
     98     }
     99   });
    100   /**
    101    * Calculate the determinant of a matrix
    102    * @param {Array[]} matrix  A square, two dimensional matrix
    103    * @param {number} rows     Number of rows of the matrix (zero-based)
    104    * @param {number} cols     Number of columns of the matrix (zero-based)
    105    * @returns {number} det
    106    * @private
    107    */
    108 
    109   function _det(matrix, rows, cols) {
    110     if (rows === 1) {
    111       // this is a 1 x 1 matrix
    112       return (0, _object.clone)(matrix[0][0]);
    113     } else if (rows === 2) {
    114       // this is a 2 x 2 matrix
    115       // the determinant of [a11,a12;a21,a22] is det = a11*a22-a21*a12
    116       return subtract(multiply(matrix[0][0], matrix[1][1]), multiply(matrix[1][0], matrix[0][1]));
    117     } else {
    118       // Compute the LU decomposition
    119       var decomp = lup(matrix); // The determinant is the product of the diagonal entries of U (and those of L, but they are all 1)
    120 
    121       var det = decomp.U[0][0];
    122 
    123       for (var _i = 1; _i < rows; _i++) {
    124         det = multiply(det, decomp.U[_i][_i]);
    125       } // The determinant will be multiplied by 1 or -1 depending on the parity of the permutation matrix.
    126       // This can be determined by counting the cycles. This is roughly a linear time algorithm.
    127 
    128 
    129       var evenCycles = 0;
    130       var i = 0;
    131       var visited = [];
    132 
    133       while (true) {
    134         while (visited[i]) {
    135           i++;
    136         }
    137 
    138         if (i >= rows) break;
    139         var j = i;
    140         var cycleLen = 0;
    141 
    142         while (!visited[decomp.p[j]]) {
    143           visited[decomp.p[j]] = true;
    144           j = decomp.p[j];
    145           cycleLen++;
    146         }
    147 
    148         if (cycleLen % 2 === 0) {
    149           evenCycles++;
    150         }
    151       }
    152 
    153       return evenCycles % 2 === 0 ? det : unaryMinus(det);
    154     }
    155   }
    156 });
    157 exports.createDet = createDet;