lower matrix ordering needs to change

Author: jacobbogers

jest.config.cjs

@@ -1,6 +1,7 @@
 const testRegex = [
     'lib/l3/complex/csyrk/__test__/test.ts',
-    'lib/l3/complex/csyrk/__test__/wasm-test.ts'
+    'lib/l3/complex/csyrk/__test__/wasm-test.ts',
+    'lib/utils/__test__/matrix-triangular.test.ts'
 ];
 
 const collectCoverageFrom = [
@@ -26,7 +27,7 @@ module.exports = {
     //testMatch: ['**/__tests__/**/*.[t]s?(x)', '**/?(*.)+(spec|test).[t]s?(x)'],
     testRegex,
     transform: {
-        "\\.test\\.ts$" :["ts-jest", {
+        "\\.test\\.ts$": ["ts-jest", {
             compiler: 'typescript',
             tsconfig: 'tsconfig.json',
             diagnostics: {

lib/utils/__test__/matrix-triangular.test.ts

@@ -0,0 +1,41 @@
+import { transformUpperPackedIdxToColumnRow, transformLowerPackedIdxToColumnRow } from '../transforms';
+import { back } from '../transforms';
+back;
+transformLowerPackedIdxToColumnRow;
+
+describe('coordinate transform testing', function () {
+    describe('upper packed', () => {
+
+        it('n=91*92/2 => col=91 row = 0w = 0', () => {
+            /**0 1 2 3 4  etc
+             * 0 1 3 6 10
+             * . 2 4 7 11
+             * . . 5 8 12
+             * . . . 9 13
+             * . . . . 14
+             */
+            const col = 91;
+            expect(transformUpperPackedIdxToColumnRow(col * (col + 1) / 2)).toEqual({ col, row: 0 });
+            expect(transformUpperPackedIdxToColumnRow(col * (col + 1) / 2 - 1)).toEqual({ col: col - 1, row: col - 1 });
+        });
+
+        describe('upper packed', () => {
+            it('N=5, for idx= 7, 10 14, 0, 3 and 9', () => {
+                /**
+                     * 0 1 .2 .3 .4
+                     * ----------
+                     * 0 . .. .. ..
+                     * 1 5 .. .. .. 
+                     * 2 6 .9 .. ..
+                     * 3 7 10 12 ..
+                     * 4 8 11 13 14
+                 */
+                const N = 5;
+                expect(transformLowerPackedIdxToColumnRow(7, N)).toEqual({ col: 1, row: 3 });
+                expect(transformLowerPackedIdxToColumnRow(10, N)).toEqual({ col: 2, row: 3 });
+                expect(transformLowerPackedIdxToColumnRow(14, N)).toEqual({ col: 4, row: 4 });
+                expect(transformLowerPackedIdxToColumnRow(0, N)).toEqual({ col: 0, row: 0 });
+            })
+        });
+    });
+});

lib/utils/isqrt.ts

@@ -0,0 +1,28 @@
+
+// from https://en.wikipedia.org/wiki/Methods_of_computing_square_roots#Binary_numeral_system_(base_2)
+export default function isqrt(n: number): number {
+    // Xₙ₊₁
+    let x = n;
+
+    // cₙ
+    let c = 0;
+
+    // dₙ which starts at the highest power of four <= n
+    let d = 1 << 30; // The second-to-top bit is set.
+    // Same as ((unsigned) INT32_MAX + 1) / 2.
+    while (d > n)
+        d >>= 2;
+
+    // for dₙ … d₀
+    while (d != 0) {
+        if (x >= c + d) {      // if Xₘ₊₁ ≥ Yₘ then aₘ = 2ᵐ
+            x -= c + d;        // Xₘ = Xₘ₊₁ - Yₘ
+            c = (c >> 1) + d;  // cₘ₋₁ = cₘ/2 + dₘ (aₘ is 2ᵐ)
+        }
+        else {
+            c >>= 1;           // cₘ₋₁ = cₘ/2      (aₘ is 0)
+        }
+        d >>= 2;               // dₘ₋₁ = dₘ/4
+    }
+    return c;                  // c₋₁
+}

lib/utils/transforms.ts

@@ -0,0 +1,64 @@
+import isqrt from './isqrt';
+
+/*
+  0 4 7 9  
+  1 5 8 .
+  2 6 . .
+  3 . . .
+*/
+
+export function transformUpperPackedIdxToColumnRow(idx: number): { col: number, row: number } {
+    const col = (isqrt((idx << 3) + 1) - 1) >> 1;
+    /**
+        * 0 1 2 3 4
+        * ----------
+        * 0 1 3 6 10
+        * . 2 4 7 11
+        * . . 5 8 12
+        * . . . 9 13
+        * . . . . 14
+    */
+    const row = idx - (col * (col + 1) >> 1);
+    return { col, row };
+}
+
+/*
+  0 . . .  
+  1 4 . .
+  2 5 7 .
+  3 6 8 9
+*/
+
+export function transformLowerPackedIdxToColumnRow(idx: number, N: number): { col: number, row: number } {
+    const b = ((N << 1) + 1);
+    const bb = b * b;
+    const ac4 = idx << 3;
+    const D = bb - ac4;
+    const col = (b - Math.sqrt(D)) >> 1;
+    /**
+        * 0 1 .2 .3 .4
+        * ----------
+        * 0 . .. .. ..
+        * 1 5 .. .. .. 
+        * 2 6 .9 .. ..
+        * 3 7 10 12 ..
+        * 4 8 11 13 14
+    */
+    const base = -col * (col - b) >> 1;
+    const row = idx - base + col;
+    return { col, row };
+}
+
+
+
+export function back(N: number, r: number): number {
+    const b = ((N << 1) + 1);
+    const bf = (2 * N + 1);
+    const bb = bf * bf;
+    //const bb = ((N * N) << 2) + (N << 2) + 1
+    const ac4 = r << 3;
+    const D = bb - ac4;
+    //console.log(`isqrt:${isqrt(D)}, sqrt:${Math.sqrt(D)}, b-isqrt:${b - isqrt(D)}, b-sqrt:${b - Math.sqrt(D)}`);
+    const col1 = (b - Math.sqrt(D)) >> 1;
+    return col1;
+}