add a trivial case for eigen values and vectors

src/function/matrix/eigs.js

@@ -2,15 +2,17 @@ import { factory } from '../../utils/factory.js'
 import { format } from '../../utils/string.js'
 import { createComplexEigs } from './eigs/complexEigs.js'
 import { createRealSymmetric } from './eigs/realSymetric.js'
+import { createTrivial } from './eigs/trivialCase.js'
 import { typeOf, isNumber, isBigNumber, isComplex, isFraction } from '../../utils/is.js'
 
 const name = 'eigs'
 
 // The absolute state of math.js's dependency system:
 const dependencies = ['config', 'typed', 'matrix', 'addScalar', 'equal', 'subtract', 'abs', 'atan', 'cos', 'sin', 'multiplyScalar', 'divideScalar', 'inv', 'bignumber', 'multiply', 'add', 'larger', 'column', 'flatten', 'number', 'complex', 'sqrt', 'diag', 'qr', 'usolve', 'usolveAll', 'im', 're', 'smaller', 'matrixFromColumns', 'dot']
 export const createEigs = /* #__PURE__ */ factory(name, dependencies, ({ config, typed, matrix, addScalar, subtract, equal, abs, atan, cos, sin, multiplyScalar, divideScalar, inv, bignumber, multiply, add, larger, column, flatten, number, complex, sqrt, diag, qr, usolve, usolveAll, im, re, smaller, matrixFromColumns, dot }) => {
-  const doRealSymetric = createRealSymmetric({ config, addScalar, subtract, column, flatten, equal, abs, atan, cos, sin, multiplyScalar, inv, bignumber, complex, multiply, add })
+  const doRealSymetric = createRealSymmetric({ config, addScalar, subtract, column, flatten, equal, abs, atan, cos, sin, multiplyScalar, inv, bignumber, complex, multiply, add, divideScalar, sqrt })
   const doComplexEigs = createComplexEigs({ config, addScalar, subtract, multiply, multiplyScalar, flatten, divideScalar, sqrt, abs, bignumber, diag, qr, inv, usolve, usolveAll, equal, complex, larger, smaller, matrixFromColumns, dot })
+  const doTrivial = createTrivial({ addScalar, subtract, multiplyScalar, sqrt, divideScalar })
 
   /**
    * Compute eigenvalues and eigenvectors of a matrix. The eigenvalues are sorted by their absolute value, ascending.
@@ -92,11 +94,17 @@ export const createEigs = /* #__PURE__ */ factory(name, dependencies, ({ config,
 
       if (isSymmetric(arr, N, prec)) {
         const type = coerceTypes(mat, arr, N)
+        if (N === 2) {
+          return doTrivial(arr)
+        }
         return doRealSymetric(arr, N, prec, type)
       }
     }
 
     const type = coerceTypes(mat, arr, N)
+    if (N === 2) {
+      return doTrivial(arr)
+    }
     return doComplexEigs(arr, N, prec, type)
   }
 

src/function/matrix/eigs/realSymetric.js

@@ -1,6 +1,6 @@
 import { clone } from '../../../utils/object.js'
 
-export function createRealSymmetric ({ config, addScalar, subtract, abs, atan, cos, sin, multiplyScalar, inv, bignumber, multiply, add }) {
+export function createRealSymmetric ({ config, addScalar, subtract, abs, atan, cos, sin, multiplyScalar, inv, bignumber, multiply, add, sqrt, divideScalar }) {
   /**
    * @param {number[] | BigNumber[]} arr
    * @param {number} N

src/function/matrix/eigs/trivialCase.js

@@ -0,0 +1,36 @@
+export function createTrivial ({ addScalar, subtract, multiplyScalar, sqrt, divideScalar }) {
+  function main (arr) {
+    const a = arr[0][0]
+    const b = arr[0][1]
+    const c = arr[1][0]
+    const d = arr[1][1]
+    const values = eigenvalues2x2(a, b, c, d)
+
+    if (values[0] === values[1]) {
+      return {
+        values: [values[0], values[0]],
+        vectors: [[0, 1], [0, 1]]
+      }
+    }
+
+    const eigenvector1 = [divideScalar(c, subtract(values[0], a)), 1]
+    const eigenvector2 = [divideScalar(c, subtract(values[1], a)), 1]
+
+    return {
+      values,
+      vectors: [eigenvector1, eigenvector2]
+    }
+  }
+
+  function eigenvalues2x2 (a, b, c, d) {
+    // λ± = ½ trA ± ½ √( tr²A - 4 detA )
+    const trA = addScalar(a, d)
+    const detA = subtract(multiplyScalar(a, d), multiplyScalar(b, c))
+    const x = multiplyScalar(trA, 0.5)
+    const y = multiplyScalar(sqrt(subtract(multiplyScalar(trA, trA), multiplyScalar(4, detA))), 0.5)
+
+    return [subtract(x, y), addScalar(x, y)]
+  }
+
+  return main
+}

test/unit-tests/function/matrix/eigs.test.js

@@ -63,6 +63,18 @@ describe('eigs', function () {
     approx.deepEqual(eigs(
       [[5, 2.3], [2.3, 1]]).values, [-0.04795013082563382, 6.047950130825635]
     )
+    approx.deepEqual(
+      eigs([[1, 2], [4, 3]]).vectors,
+      [[-2, 1], [1, 1]]
+    )
+    approx.deepEqual(
+      eigs([[5, 2.3], [2.3, 1]]).vectors,
+      [[-0.45563, 1], [2.19476, 1]]
+    )
+    approx.deepEqual(
+      eigs([[2.0, 1.0], [0.0, 2.0]]).vectors,
+      [[0, 1], [0, 1]]
+    )
   })
 
   it('calculates eigenvalues for 2x2 matrix with complex entries', function () {
@@ -201,8 +213,8 @@ describe('eigs', function () {
     const eig = eigs(B)
 
     assert.strictEqual(eig.values[0].toString(),
-      '-0.9999999999999999999999999999999999999999999999999999999999999999')
+      '-1')
     assert.strictEqual(eig.values[1].toString(),
-      '0.9999999999999999999999999999999999999999999999999999999999999999')
+      '1')
   })
 })