diff --git a/QuadTree.Tests/Tests.Matrix.fs b/QuadTree.Tests/Tests.Matrix.fs index 90b678d..f5fcc84 100644 --- a/QuadTree.Tests/Tests.Matrix.fs +++ b/QuadTree.Tests/Tests.Matrix.fs @@ -226,3 +226,131 @@ let ``Simple Matrix.map2. Square where number of cols and rows are not power of let eq = actual = expected Assert.True(eq) + +[] +let ``Conversion identity`` () = + let id = toCoordinateList << fromCoordinateList + + let nrows = 10UL + let ncols = 12UL + + let data = + [ 0UL, 3UL, 10 + 3UL, 3UL, 33 + 9UL, 2UL, 5 + 3UL, 11UL, 1 ] + |> List.sort + + let coordinates = CoordinateList(nrows, ncols, data) + + let expected = coordinates + let actual = id coordinates + + Assert.Equal(expected, actual) + +[] +let ``Simple addition`` () = + let nrows = 10UL + let ncols = 12UL + + let d1 = + [ 0UL, 3UL, 4 + 9UL, 2UL, 5 + 3UL, 11UL, 2 ] + + let d2 = + [ 0UL, 3UL, 6 + 3UL, 3UL, 33 + 3UL, 11UL, -1 ] + + let expected = + let expectedList = + [ 0UL, 3UL, 10 + 3UL, 3UL, 33 + 9UL, 2UL, 5 + 3UL, 11UL, 1 ] + |> List.sort + + CoordinateList(nrows, ncols, expectedList) + + let actual = + let c1 = CoordinateList(nrows, ncols, d1) + let c2 = CoordinateList(nrows, ncols, d2) + let m1 = fromCoordinateList c1 + let m2 = fromCoordinateList c2 + + let addition o1 o2 = + match o1, o2 with + | Some x, Some y -> Some(x + y) + | Some x, None + | None, Some x -> Some x + | None, None -> None + + let result = + match map2 m1 m2 addition with + | Result.Success x -> x + | _ -> failwith "Unreachable" + + toCoordinateList result + + Assert.Equal(expected, actual) + +[] +let ``Condensation of empty`` () = + let clist = CoordinateList(2UL, 3UL, []) + + let actual = fromCoordinateList clist + + // 2 * 3 = 5 + // 4 * 4 None and Dummy + // NN N D + // NN N D + // DDDD + // DDDD + let tree = + qtree.Node( + qtree.Leaf <| UserValue None, + qtree.Node(qtree.Leaf <| UserValue None, qtree.Leaf Dummy, qtree.Leaf <| UserValue None, qtree.Leaf Dummy), + qtree.Leaf Dummy, + qtree.Leaf Dummy + ) + + let expected = + SparseMatrix(2UL, 3UL, 0UL, Storage(4UL, 4UL, tree)) + + Assert.Equal(expected.storage.data, actual.storage.data) + +[] +let ``Condensation of sparse`` () = + let clist = + CoordinateList(4UL, 3UL, [ 0UL, 2UL, 2; 3UL, 2UL, 4 ]) + + let actual = fromCoordinateList clist + + // NN2D + // NNND + // NNND + // NN4D + + let tree = + qtree.Node( + qtree.Leaf <| UserValue None, + qtree.Node( + qtree.Leaf << UserValue <| Some 2, + qtree.Leaf Dummy, + qtree.Leaf <| UserValue None, + qtree.Leaf Dummy + ), + qtree.Leaf <| UserValue None, + qtree.Node( + qtree.Leaf <| UserValue None, + qtree.Leaf Dummy, + qtree.Leaf << UserValue <| Some 4, + qtree.Leaf Dummy + ) + ) + + let expected = + SparseMatrix(4UL, 3UL, 0UL, Storage(4UL, 4UL, tree)) + + Assert.Equal(expected.storage.data, actual.storage.data) diff --git a/QuadTree.Tests/Tests.Vector.fs b/QuadTree.Tests/Tests.Vector.fs index 1f3f32f..6aa6d57 100644 --- a/QuadTree.Tests/Tests.Vector.fs +++ b/QuadTree.Tests/Tests.Vector.fs @@ -57,7 +57,6 @@ let ``Simple Vector.map2. Length is power of two.`` () = Assert.True(eq) - [] let ``Simple Vector.map2. Length is not power of two.`` () = let v1 = @@ -100,3 +99,71 @@ let ``Simple Vector.map2. Length is not power of two.`` () = let eq = actual = expected Assert.True(eq) + +[] +let ``Conversion identity`` () = + let id = toCoordinateList << fromCoordinateList + + let dataLength = 10UL + + let data = + [ 0UL, 3; 3UL, -1; 7UL, 2; 8UL, 2; 9UL, 2 ] + + let coordinates = CoordinateList(dataLength, data) + + let expected = coordinates + let actual = id coordinates + + Assert.Equal(expected, actual) + +[] +let ``Simple addition`` () = + let dataLength = 10UL + + let d1 = [ 0UL, 2; 9UL, 1 ] + let d2 = [ 0UL, 3; 8UL, 1 ] + + let expected = + let expectedList = [ 0UL, 5; 8UL, 1; 9UL, 1 ] + CoordinateList(dataLength, expectedList) + + let actual = + let c1 = CoordinateList(dataLength, d1) + let c2 = CoordinateList(dataLength, d2) + let v1 = fromCoordinateList c1 + let v2 = fromCoordinateList c2 + + let addition o1 o2 = + match o1, o2 with + | Some x, Some y -> Some(x + y) + | Some x, None + | None, Some x -> Some x + | None, None -> None + + let result = + match map2 v1 v2 addition with + | Result.Success x -> x + | _ -> failwith "Unreachable" + + toCoordinateList result + + Assert.Equal(expected, actual) + +[] +let ``Condensation of empty`` () = + let clist = CoordinateList(10UL, []) + + let actual = fromCoordinateList clist + + // 16 elements total None and Dummy: NNNNNNNN | NN DD | DDDD + let tree = + + btree.Node( + btree.Leaf <| UserValue None, + btree.Node(btree.Node(btree.Leaf <| UserValue None, btree.Leaf Dummy), btree.Leaf Dummy) + ) + + let expected = + SparseVector(clist.length, 0UL, Storage(16UL, tree)) + + Assert.Equal(expected, actual) diff --git a/QuadTree/Matrix.fs b/QuadTree/Matrix.fs index 4133f60..057b0c8 100644 --- a/QuadTree/Matrix.fs +++ b/QuadTree/Matrix.fs @@ -57,6 +57,98 @@ let mkNode x1 x2 x3 x4 = | Leaf(v1), Leaf(v2), Leaf(v3), Leaf(v4) when v1 = v2 && v2 = v3 && v3 = v4 -> Leaf(v1) | _ -> Node(x1, x2, x3, x4) +[] +type rowindex + +[] +type colindex + +type COOEntry<'value> = uint64 * uint64 * 'value + +[] +type CoordinateList<'value> = + val nrows: uint64 + val ncols: uint64 + val list: COOEntry<'value> list + + new(_nrows, _ncols, _list) = + { nrows = _nrows + ncols = _ncols + list = _list } + +let private getQuadrantCoords (pr, pc) halfSize = + (pr, pc), // NORTH WEST + (pr, pc + halfSize * 1UL), // NORTH EAST + (pr + halfSize * 1UL, pc), // SOUTH WEST + (pr + halfSize * 1UL, pc + halfSize * 1UL) // SOUTH EAST + +let fromCoordinateList (coo: CoordinateList<'a>) = + let nvals = (uint64 <| List.length coo.list) * 1UL + let nrows = coo.nrows + let ncols = coo.ncols + + // the resulting matrix is always square + let storageSize = getNearestUpperPowerOfTwo (max (uint64 nrows) (uint64 ncols)) + + let isEntryInQuadrant (pr, pc) size (entry: COOEntry<'a>) = + let (i, j, _) = entry + + i >= pr + && j >= pc + && i < pr + size * 1UL + && j < pc + size * 1UL + + let rec traverse coordinates (pr, pc) size = + match coordinates with + | [] when (uint64 pr) + size < uint64 nrows && (uint64 pc) + size < uint64 ncols -> Leaf <| UserValue None + | [] when uint64 pr >= uint64 nrows || uint64 pc >= uint64 ncols -> Leaf Dummy + | (i, j, value) :: _ when pr = i && pc = j && size = 1UL -> Leaf << UserValue <| Some value + | _ -> + let halfSize = size / 2UL + let nwp, nep, swp, sep = getQuadrantCoords (pr, pc) halfSize + let nwCoo = coordinates |> List.filter (isEntryInQuadrant nwp halfSize) + let neCoo = coordinates |> List.filter (isEntryInQuadrant nep halfSize) + let swCoo = coordinates |> List.filter (isEntryInQuadrant swp halfSize) + let seCoo = coordinates |> List.filter (isEntryInQuadrant sep halfSize) + + mkNode + (traverse nwCoo nwp halfSize) + (traverse neCoo nep halfSize) + (traverse swCoo swp halfSize) + (traverse seCoo sep halfSize) + + let tree = traverse coo.list (0UL, 0UL) storageSize + + SparseMatrix(nrows, ncols, nvals, Storage(storageSize * 1UL, storageSize * 1UL, tree)) + +let toCoordinateList (matrix: SparseMatrix<'a>) = + let nrows = matrix.nrows + let ncols = matrix.ncols + + let rec traverse tree (pr, pc) size = + match tree with + | Leaf Dummy + | Leaf(UserValue None) -> [] + | Leaf(UserValue(Some value)) -> + [ for i in uint64 pr .. (uint64 pr) + size - 1UL do + for j in uint64 pc .. (uint64 pc) + size - 1UL -> (i * 1UL, j * 1UL, value) ] + | Node(nw, ne, sw, se) -> + let halfSize = size / 2UL + let nwp, nep, swp, sep = getQuadrantCoords (pr, pc) halfSize + + traverse nw nwp halfSize + @ traverse ne nep halfSize + @ traverse sw swp halfSize + @ traverse se sep halfSize + + let coo = + traverse + matrix.storage.data + (0UL, 0UL) + (max (uint64 matrix.storage.hSize) (uint64 matrix.storage.vSize)) + + CoordinateList(nrows, ncols, coo) + let map2 (matrix1: SparseMatrix<_>) (matrix2: SparseMatrix<_>) f = let rec inner (vSize: uint64) (hSize: uint64) matrix1 matrix2 = let _do x1 x2 x3 x4 y1 y2 y3 y4 = diff --git a/QuadTree/Vector.fs b/QuadTree/Vector.fs index 35c945b..dd0d8a1 100644 --- a/QuadTree/Vector.fs +++ b/QuadTree/Vector.fs @@ -38,11 +38,68 @@ let foldValues state f tree = | Leaf *) + let mkNode t1 t2 = match (t1, t2) with | Leaf(v1), Leaf(v2) when v1 = v2 -> Leaf(v1) | _ -> Node(t1, t2) +[] +type index + +[] +type CoordinateList<'value> = + val length: uint64 + val data: (uint64 * 'value) list + new(_length, _data) = { length = _length; data = _data } + +let fromCoordinateList (lst: CoordinateList<'a>) : SparseVector<'a> = + let length = lst.length + let nvals = (uint64 <| List.length lst.data) * 1UL + let storageSize = (getNearestUpperPowerOfTwo <| uint64 length) * 1UL + + let rec traverse coordinates pointer size = + match coordinates with + | [] when uint64 (pointer + size) < uint64 (length) -> Leaf <| UserValue None, [] + | [] when uint64 pointer >= uint64 length -> Leaf Dummy, [] + | (idx, _) :: _ when idx > pointer + size -> Leaf <| UserValue None, coordinates + | (idx, value) :: xs when idx = pointer && size = 1UL -> Leaf << UserValue <| Some value, xs + | _ -> + let halfSize = size / 2UL + + let left, lCoordinates = traverse coordinates pointer halfSize + let right, rCoordinates = traverse lCoordinates (pointer + halfSize) halfSize + + mkNode left right, rCoordinates + + let sortedCoordinates = List.sort lst.data + + let tree, _ = + traverse sortedCoordinates 0UL ((uint64 storageSize) * 1UL) + + SparseVector(length, nvals, Storage(storageSize, tree)) + +let toCoordinateList (vector: SparseVector<'a>) = + let length = vector.length + + let rec traverse tree accum (pointer: uint64) (size: uint64) = + match tree with + | Leaf Dummy + | Leaf(UserValue(None)) -> accum + | Leaf(UserValue(Some value)) -> + accum + @ [ for idx in 0UL .. uint64 (size - 1UL) -> (pointer + idx * 1UL, value) ] + | Node(left, right) -> + let halfSize = size / 2UL + let lAccum = traverse left accum pointer halfSize + let rAccum = traverse right lAccum (pointer + halfSize) halfSize + rAccum + + let lst = + traverse vector.storage.data [] 0UL ((uint64 vector.storage.size) * 1UL) + + CoordinateList(length, lst) + let map2 (vector1: SparseVector<'a>) (vector2: SparseVector<'b>) f = let len1 = vector1.length