Mercurial > repos > public > sbplib_julia

changeset 1736:863385aae454 feature/grids/curvilinear

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Merge default
author Jonatan Werpers <jonatan@werpers.com>
date Tue, 10 Sep 2024 21:59:10 +0200
parents 36986b75bf98 (current diff) f215ac2a5c66 (diff)
children a4e1721a7109
files Manifest.toml Project.toml src/Grids/Grids.jl src/Grids/equidistant_grid.jl src/Grids/grid.jl src/Grids/mapped_grid.jl src/Grids/tensor_grid.jl test/Grids/mapped_grid_test.jl test/Grids/tensor_grid_test.jl
diffstat 67 files changed, 611 insertions(+), 433 deletions(-) [+]
line wrap: on
line diff
--- a/Manifest.toml	Tue Sep 10 21:35:55 2024 +0200
+++ b/Manifest.toml	Tue Sep 10 21:59:10 2024 +0200
@@ -1,6 +1,6 @@
 # This file is machine-generated - editing it directly is not advised
 
-julia_version = "1.10.4"
+julia_version = "1.10.5"
 manifest_format = "2.0"
 project_hash = "32fac879810099260f177c27318d3f26de4a00cc"
 
@@ -16,9 +16,9 @@
 
 [[deps.ArrayInterface]]
 deps = ["Adapt", "LinearAlgebra"]
-git-tree-sha1 = "f54c23a5d304fb87110de62bace7777d59088c34"
+git-tree-sha1 = "3640d077b6dafd64ceb8fd5c1ec76f7ca53bcf76"
 uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9"
-version = "7.15.0"
+version = "7.16.0"
 
     [deps.ArrayInterface.extensions]
     ArrayInterfaceBandedMatricesExt = "BandedMatrices"
@@ -184,4 +184,4 @@
 [[deps.libblastrampoline_jll]]
 deps = ["Artifacts", "Libdl"]
 uuid = "8e850b90-86db-534c-a0d3-1478176c7d93"
-version = "5.8.0+1"
+version = "5.11.0+0"
--- a/Project.toml	Tue Sep 10 21:35:55 2024 +0200
+++ b/Project.toml	Tue Sep 10 21:59:10 2024 +0200
@@ -1,7 +1,7 @@
-name = "Sbplib"
+name = "Diffinitive"
 uuid = "5a373a26-915f-4769-bcab-bf03835de17b"
-authors = ["Jonatan Werpers <jonatan@werpers.com>", "Vidar Stiernström <vidar.stiernstrom@it.uu.se>, and contributors"]
-version = "0.1.1"
+authors = ["Jonatan Werpers <jonatan@werpers.com>", "Vidar Stiernström <vidar.stiernstrom@gmail.com>, and contributors"]
+version = "0.1.2"
 
 [deps]
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
@@ -11,9 +11,14 @@
 
 [weakdeps]
 Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
+SparseArrayKit = "a9a3c162-d163-4c15-8926-b8794fbefed2"
+SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+Tokens = "040c2ec2-8d69-4aca-bf03-7d3a7092f2f6"
 
 [extensions]
-SbplibMakieExt = "Makie"
+DiffinitiveMakieExt = "Makie"
+DiffinitiveSparseArrayKitExt = ["SparseArrayKit", "Tokens"]
+DiffinitiveSparseArraysExt = ["SparseArrays", "Tokens"]
 
 [compat]
 julia = "1.5"
--- a/README.md	Tue Sep 10 21:35:55 2024 +0200
+++ b/README.md	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,10 @@
-# Sbplib
+# Diffinitive
 
 ## Running tests
 To run all tests simply run
 ```
 (@v1.5) pkg> activate .
-(Sbplib) pkg> test
+(Diffinitive) pkg> test
 ```
 
 If you want to run tests from a specific file in `test/`, you can do
@@ -59,8 +59,8 @@
 
 ```julia
 using PkgBenchmark
-import Sbplib
-r = benchmarkpkg(Sbplib)
+import Diffinitive
+r = benchmarkpkg(Diffinitive)
 
 export_markdown(stdout, r)
 ```
--- a/TODO.md	Tue Sep 10 21:35:55 2024 +0200
+++ b/TODO.md	Tue Sep 10 21:59:10 2024 +0200
@@ -33,6 +33,7 @@
  - [ ] Kolla att vi har @inline på rätt ställen
  - [ ] Profilera
 
+ - [ ] Keep a lookout for allowing dependencies of package extensions (https://github.com/JuliaLang/Pkg.jl/issues/3641) This should be used to simplify the matrix extensions so that you don't have to load Tokens which is only used internally to the extension
 
 ### Grids
 
--- a/benchmark/Manifest.toml	Tue Sep 10 21:35:55 2024 +0200
+++ b/benchmark/Manifest.toml	Tue Sep 10 21:59:10 2024 +0200
@@ -1,8 +1,8 @@
 # This file is machine-generated - editing it directly is not advised
 
-julia_version = "1.10.4"
+julia_version = "1.10.5"
 manifest_format = "2.0"
-project_hash = "25bba7b4a00465d5a2b00b589eb10e3301c31f2a"
+project_hash = "ecfc3e12aca5be17a874aba6134ff821abf61540"
 
 [[deps.AbstractTrees]]
 git-tree-sha1 = "2d9c9a55f9c93e8887ad391fbae72f8ef55e1177"
@@ -25,9 +25,9 @@
 
 [[deps.ArrayInterface]]
 deps = ["Adapt", "LinearAlgebra"]
-git-tree-sha1 = "f54c23a5d304fb87110de62bace7777d59088c34"
+git-tree-sha1 = "3640d077b6dafd64ceb8fd5c1ec76f7ca53bcf76"
 uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9"
-version = "7.15.0"
+version = "7.16.0"
 
     [deps.ArrayInterface.extensions]
     ArrayInterfaceBandedMatricesExt = "BandedMatrices"
@@ -99,6 +99,23 @@
 deps = ["Printf"]
 uuid = "ade2ca70-3891-5945-98fb-dc099432e06a"
 
+[[deps.Diffinitive]]
+deps = ["StaticArrays", "TOML", "TiledIteration"]
+path = ".."
+uuid = "5a373a26-915f-4769-bcab-bf03835de17b"
+version = "0.1.1"
+
+    [deps.Diffinitive.extensions]
+    DiffinitiveMakieExt = "Makie"
+    DiffinitiveSparseArrayKitExt = ["SparseArrayKit", "Tokens"]
+    DiffinitiveSparseArraysExt = ["SparseArrays", "Tokens"]
+
+    [deps.Diffinitive.weakdeps]
+    Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
+    SparseArrayKit = "a9a3c162-d163-4c15-8926-b8794fbefed2"
+    SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+    Tokens = "040c2ec2-8d69-4aca-bf03-7d3a7092f2f6"
+
 [[deps.Downloads]]
 deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"]
 uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6"
@@ -273,18 +290,6 @@
 uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce"
 version = "0.7.0"
 
-[[deps.Sbplib]]
-deps = ["StaticArrays", "TOML", "TiledIteration"]
-path = ".."
-uuid = "5a373a26-915f-4769-bcab-bf03835de17b"
-version = "0.1.1"
-
-    [deps.Sbplib.extensions]
-    SbplibMakieExt = "Makie"
-
-    [deps.Sbplib.weakdeps]
-    Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
-
 [[deps.Serialization]]
 uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b"
 
@@ -391,7 +396,7 @@
 [[deps.libblastrampoline_jll]]
 deps = ["Artifacts", "Libdl"]
 uuid = "8e850b90-86db-534c-a0d3-1478176c7d93"
-version = "5.8.0+1"
+version = "5.11.0+0"
 
 [[deps.nghttp2_jll]]
 deps = ["Artifacts", "Libdl"]
--- a/benchmark/Project.toml	Tue Sep 10 21:35:55 2024 +0200
+++ b/benchmark/Project.toml	Tue Sep 10 21:59:10 2024 +0200
@@ -1,5 +1,5 @@
 [deps]
 BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
+Diffinitive = "5a373a26-915f-4769-bcab-bf03835de17b"
 Mustache = "ffc61752-8dc7-55ee-8c37-f3e9cdd09e70"
-PkgBenchmark = "32113eaa-f34f-5b0d-bd6c-c81e245fc73d"
-Sbplib = "5a373a26-915f-4769-bcab-bf03835de17b"
+PkgBenchmark = "32113eaa-f34f-5b0d-bd6c-c81e245fc73d"
\ No newline at end of file
--- a/benchmark/benchmark_laplace.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/benchmark/benchmark_laplace.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,7 +1,7 @@
-using Sbplib
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.RegionIndices
+using Diffinitive
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.RegionIndices
 using BenchmarkTools
 
 # TODO: Move the below benchmarks into the benchmark suite
--- a/benchmark/benchmark_utils.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/benchmark/benchmark_utils.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -3,11 +3,11 @@
 import Mustache
 import Dates
 
-import Sbplib
+import Diffinitive
 
-const sbplib_root = splitpath(pathof(Sbplib))[1:end-2] |> joinpath
-const results_dir = mkpath(joinpath(sbplib_root, "benchmark/results"))
-const template_path = joinpath(sbplib_root, "benchmark/result.tmpl")
+const diffinitive_root = splitpath(pathof(Diffinitive))[1:end-2] |> joinpath
+const results_dir = mkpath(joinpath(diffinitive_root, "benchmark/results"))
+const template_path = joinpath(diffinitive_root, "benchmark/result.tmpl")
 
 """
     mainmain(;rev=nothing, target=nothing, baseline=nothing , kwargs...)
@@ -49,7 +49,7 @@
 `PkgBenchmark.BenchmarkResult`
 """
 function run_benchmark(;kwargs...)
-    r = PkgBenchmark.benchmarkpkg(Sbplib; kwargs...)
+    r = PkgBenchmark.benchmarkpkg(Diffinitive; kwargs...)
 
     rev = hg_rev() # Should be changed to hg_id() when the html can handle it.
 
@@ -158,17 +158,17 @@
 
 
 function hg_id()
-    cmd = Cmd(`hg id`, dir=sbplib_root)
+    cmd = Cmd(`hg id`, dir=diffinitive_root)
     return readchomp(addenv(cmd, "HGPLAIN"=>""))
 end
 
 function hg_rev()
-    cmd = Cmd(`hg id -i`, dir=sbplib_root)
+    cmd = Cmd(`hg id -i`, dir=diffinitive_root)
     return readchomp(addenv(cmd, "HGPLAIN"=>""))
 end
 
 function hg_update(rev)
-    cmd = Cmd(`hg update --check -r $rev`, dir=sbplib_root)
+    cmd = Cmd(`hg update --check -r $rev`, dir=diffinitive_root)
     run(addenv(cmd, "HGPLAIN"=>""))
 
     return nothing
@@ -182,9 +182,9 @@
 """
 function hg_commit(msg; secret=false)
     if secret
-        cmd = Cmd(`hg commit --verbose --secret --message $msg`, dir=sbplib_root)
+        cmd = Cmd(`hg commit --verbose --secret --message $msg`, dir=diffinitive_root)
     else
-        cmd = Cmd(`hg commit --verbose          --message $msg`, dir=sbplib_root)
+        cmd = Cmd(`hg commit --verbose          --message $msg`, dir=diffinitive_root)
     end
 
     out = readchomp(addenv(cmd, "HGPLAIN"=>""))
@@ -200,9 +200,9 @@
 """
 function hg_strip(rev; keep=false)
     if keep
-        cmd = Cmd(`hg --config extensions.strip= strip --keep -r $rev`, dir=sbplib_root)
+        cmd = Cmd(`hg --config extensions.strip= strip --keep -r $rev`, dir=diffinitive_root)
     else
-        cmd = Cmd(`hg --config extensions.strip= strip        -r $rev`, dir=sbplib_root)
+        cmd = Cmd(`hg --config extensions.strip= strip        -r $rev`, dir=diffinitive_root)
     end
 
     run(addenv(cmd, "HGPLAIN"=>""))
@@ -216,7 +216,7 @@
 Return true if the repositopry has uncommited changes.
 """
 function hg_is_dirty()
-    cmd = Cmd(`hg identify --id`, dir=sbplib_root)
+    cmd = Cmd(`hg identify --id`, dir=diffinitive_root)
     out = readchomp(addenv(cmd, "HGPLAIN"=>""))
 
     return endswith(out, "+")
--- a/benchmark/benchmarks.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/benchmark/benchmarks.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,9 @@
 using BenchmarkTools
 
-using Sbplib
-using Sbplib.Grids
-using Sbplib.SbpOperators
-using Sbplib.RegionIndices
-using Sbplib.LazyTensors
+using Diffinitive
+using Diffinitive.Grids
+using Diffinitive.SbpOperators
+using Diffinitive.LazyTensors
 
 using LinearAlgebra
 
@@ -168,20 +167,20 @@
 Dxx = second_derivative(g2, stencil_set, 1)
 Dyy = second_derivative(g2, stencil_set, 2)
 
-e₁ₗ = boundary_restriction(g2, stencil_set, CartesianBoundary{1,Lower}())
-e₁ᵤ = boundary_restriction(g2, stencil_set, CartesianBoundary{1,Upper}())
-e₂ₗ = boundary_restriction(g2, stencil_set, CartesianBoundary{2,Lower}())
-e₂ᵤ = boundary_restriction(g2, stencil_set, CartesianBoundary{2,Upper}())
+e₁ₗ = boundary_restriction(g2, stencil_set, CartesianBoundary{1,LowerBoundary}())
+e₁ᵤ = boundary_restriction(g2, stencil_set, CartesianBoundary{1,UpperBoundary}())
+e₂ₗ = boundary_restriction(g2, stencil_set, CartesianBoundary{2,LowerBoundary}())
+e₂ᵤ = boundary_restriction(g2, stencil_set, CartesianBoundary{2,UpperBoundary}())
 
-d₁ₗ = normal_derivative(g2, stencil_set, CartesianBoundary{1,Lower}())
-d₁ᵤ = normal_derivative(g2, stencil_set, CartesianBoundary{1,Upper}())
-d₂ₗ = normal_derivative(g2, stencil_set, CartesianBoundary{2,Lower}())
-d₂ᵤ = normal_derivative(g2, stencil_set, CartesianBoundary{2,Upper}())
+d₁ₗ = normal_derivative(g2, stencil_set, CartesianBoundary{1,LowerBoundary}())
+d₁ᵤ = normal_derivative(g2, stencil_set, CartesianBoundary{1,UpperBoundary}())
+d₂ₗ = normal_derivative(g2, stencil_set, CartesianBoundary{2,LowerBoundary}())
+d₂ᵤ = normal_derivative(g2, stencil_set, CartesianBoundary{2,UpperBoundary}())
 
-H₁ₗ = inner_product(boundary_grid(g2, CartesianBoundary{1,Lower}()), stencil_set)
-H₁ᵤ = inner_product(boundary_grid(g2, CartesianBoundary{1,Upper}()), stencil_set)
-H₂ₗ = inner_product(boundary_grid(g2, CartesianBoundary{2,Lower}()), stencil_set)
-H₂ᵤ = inner_product(boundary_grid(g2, CartesianBoundary{2,Upper}()), stencil_set)
+H₁ₗ = inner_product(boundary_grid(g2, CartesianBoundary{1,LowerBoundary}()), stencil_set)
+H₁ᵤ = inner_product(boundary_grid(g2, CartesianBoundary{1,UpperBoundary}()), stencil_set)
+H₂ₗ = inner_product(boundary_grid(g2, CartesianBoundary{2,LowerBoundary}()), stencil_set)
+H₂ᵤ = inner_product(boundary_grid(g2, CartesianBoundary{2,UpperBoundary}()), stencil_set)
 
 SUITE["boundary_terms"]["pre_composition"] = @benchmarkable $u2 .= $(H⁻¹∘e₁ₗ'∘H₁ₗ∘d₁ₗ)*$v2
 SUITE["boundary_terms"]["composition"]     = @benchmarkable $u2 .= ($H⁻¹∘$e₁ₗ'∘$H₁ₗ∘$d₁ₗ)*$v2
--- a/docs/Manifest.toml	Tue Sep 10 21:35:55 2024 +0200
+++ b/docs/Manifest.toml	Tue Sep 10 21:59:10 2024 +0200
@@ -1,8 +1,8 @@
 # This file is machine-generated - editing it directly is not advised
 
-julia_version = "1.10.4"
+julia_version = "1.10.5"
 manifest_format = "2.0"
-project_hash = "4f0756199bb5f6739a5f4697152617efc4e0705c"
+project_hash = "c04450bb5c379e77d137cc05c4c0ab58eb1bfae9"
 
 [[deps.ANSIColoredPrinters]]
 git-tree-sha1 = "574baf8110975760d391c710b6341da1afa48d8c"
@@ -30,9 +30,9 @@
 
 [[deps.ArrayInterface]]
 deps = ["Adapt", "LinearAlgebra"]
-git-tree-sha1 = "f54c23a5d304fb87110de62bace7777d59088c34"
+git-tree-sha1 = "3640d077b6dafd64ceb8fd5c1ec76f7ca53bcf76"
 uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9"
-version = "7.15.0"
+version = "7.16.0"
 
     [deps.ArrayInterface.extensions]
     ArrayInterfaceBandedMatricesExt = "BandedMatrices"
@@ -94,6 +94,23 @@
 deps = ["Printf"]
 uuid = "ade2ca70-3891-5945-98fb-dc099432e06a"
 
+[[deps.Diffinitive]]
+deps = ["StaticArrays", "TOML", "TiledIteration"]
+path = ".."
+uuid = "5a373a26-915f-4769-bcab-bf03835de17b"
+version = "0.1.1"
+
+    [deps.Diffinitive.extensions]
+    DiffinitiveMakieExt = "Makie"
+    DiffinitiveSparseArrayKitExt = ["SparseArrayKit", "Tokens"]
+    DiffinitiveSparseArraysExt = ["SparseArrays", "Tokens"]
+
+    [deps.Diffinitive.weakdeps]
+    Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
+    SparseArrayKit = "a9a3c162-d163-4c15-8926-b8794fbefed2"
+    SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+    Tokens = "040c2ec2-8d69-4aca-bf03-7d3a7092f2f6"
+
 [[deps.DocStringExtensions]]
 deps = ["LibGit2"]
 git-tree-sha1 = "2fb1e02f2b635d0845df5d7c167fec4dd739b00d"
@@ -102,9 +119,9 @@
 
 [[deps.Documenter]]
 deps = ["ANSIColoredPrinters", "AbstractTrees", "Base64", "CodecZlib", "Dates", "DocStringExtensions", "Downloads", "Git", "IOCapture", "InteractiveUtils", "JSON", "LibGit2", "Logging", "Markdown", "MarkdownAST", "Pkg", "PrecompileTools", "REPL", "RegistryInstances", "SHA", "TOML", "Test", "Unicode"]
-git-tree-sha1 = "9d29b99b6b2b6bc2382a4c8dbec6eb694f389853"
+git-tree-sha1 = "5a1ee886566f2fa9318df1273d8b778b9d42712d"
 uuid = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
-version = "1.6.0"
+version = "1.7.0"
 
 [[deps.Downloads]]
 deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"]
@@ -149,9 +166,9 @@
 
 [[deps.JLLWrappers]]
 deps = ["Artifacts", "Preferences"]
-git-tree-sha1 = "7e5d6779a1e09a36db2a7b6cff50942a0a7d0fca"
+git-tree-sha1 = "f389674c99bfcde17dc57454011aa44d5a260a40"
 uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210"
-version = "1.5.0"
+version = "1.6.0"
 
 [[deps.JSON]]
 deps = ["Dates", "Mmap", "Parsers", "Unicode"]
@@ -246,9 +263,9 @@
 
 [[deps.OpenSSL_jll]]
 deps = ["Artifacts", "JLLWrappers", "Libdl"]
-git-tree-sha1 = "a028ee3cb5641cccc4c24e90c36b0a4f7707bdf5"
+git-tree-sha1 = "1b35263570443fdd9e76c76b7062116e2f374ab8"
 uuid = "458c3c95-2e84-50aa-8efc-19380b2a3a95"
-version = "3.0.14+0"
+version = "3.0.15+0"
 
 [[deps.PCRE2_jll]]
 deps = ["Artifacts", "Libdl"]
@@ -306,18 +323,6 @@
 uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce"
 version = "0.7.0"
 
-[[deps.Sbplib]]
-deps = ["StaticArrays", "TOML", "TiledIteration"]
-path = ".."
-uuid = "5a373a26-915f-4769-bcab-bf03835de17b"
-version = "0.1.1"
-
-    [deps.Sbplib.extensions]
-    SbplibMakieExt = "Makie"
-
-    [deps.Sbplib.weakdeps]
-    Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
-
 [[deps.Serialization]]
 uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b"
 
@@ -400,7 +405,7 @@
 [[deps.libblastrampoline_jll]]
 deps = ["Artifacts", "Libdl"]
 uuid = "8e850b90-86db-534c-a0d3-1478176c7d93"
-version = "5.8.0+1"
+version = "5.11.0+0"
 
 [[deps.nghttp2_jll]]
 deps = ["Artifacts", "Libdl"]
--- a/docs/Project.toml	Tue Sep 10 21:35:55 2024 +0200
+++ b/docs/Project.toml	Tue Sep 10 21:59:10 2024 +0200
@@ -1,3 +1,3 @@
 [deps]
+Diffinitive = "5a373a26-915f-4769-bcab-bf03835de17b"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
-Sbplib = "5a373a26-915f-4769-bcab-bf03835de17b"
--- a/docs/make.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/docs/make.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,14 +1,14 @@
 using Documenter
-using Sbplib
+using Diffinitive
 
-using Sbplib.DiffOps
-using Sbplib.Grids
-using Sbplib.LazyTensors
-using Sbplib.RegionIndices
-using Sbplib.SbpOperators
-using Sbplib.StaticDicts
+using Diffinitive.DiffOps
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
+using Diffinitive.RegionIndices
+using Diffinitive.SbpOperators
+using Diffinitive.StaticDicts
 
-sitename = "Sbplib.jl"
+sitename = "Diffinitive.jl"
 
 remotes = nothing
 edit_link = nothing
@@ -31,6 +31,7 @@
     "Home" => "index.md",
     "operator_file_format.md",
     "grids_and_grid_functions.md",
+    "matrix_and_tensor_representations.md",
     "Submodules" => [
         "submodules/grids.md",
         "submodules/diff_ops.md",
--- a/docs/src/index.md	Tue Sep 10 21:35:55 2024 +0200
+++ b/docs/src/index.md	Tue Sep 10 21:59:10 2024 +0200
@@ -1,4 +1,4 @@
-# Sbplib.jl
+# Diffinitive.jl
 
 ```@contents
 Depth = 1
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/docs/src/matrix_and_tensor_representations.md	Tue Sep 10 21:59:10 2024 +0200
@@ -0,0 +1,5 @@
+# Matrix and tensor representations
+
+Sparse matrix and sparse tensor representations of lazy tensors can be constructed by loading [Tokens.jl](http://) and one of SparseArrays.jl or [SparseArrayKit.jl](http://). Through package extensions the following methods `sparse(::LazyTensor)` and `SparseArray(::LazyTensor)` are provided.
+
+<!-- TODO figure out how to add the docstrings here --/>
--- a/docs/src/operator_file_format.md	Tue Sep 10 21:35:55 2024 +0200
+++ b/docs/src/operator_file_format.md	Tue Sep 10 21:59:10 2024 +0200
@@ -1,6 +1,6 @@
 # Operator file format
 
-The intention is that Sbplib.jl should be a general and extensible framework
+The intention is that Diffinitive.jl should be a general and extensible framework
 for working with finite difference methods. It therefore includes a set of
 tools for storing and sharing operator definitions as well as a set of widely
 used operators.
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/ext/DiffinitiveMakieExt.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -0,0 +1,88 @@
+module DiffinitiveMakieExt
+
+using Diffinitive.Grids
+using Makie
+using StaticArrays
+
+
+function verticies_and_faces_and_values(g::Grid{<:Any,2}, gf::AbstractArray{<:Any, 2})
+    ps = map(Tuple, g)[:]
+    values = gf[:]
+    faces = Vector{NTuple{3,Int}}()
+
+    n,m = size(g)
+    Li = LinearIndices((1:n, 1:m))
+    for i ∈ 1:n-1, j = 1:m-1
+
+        # Add point in the middle of the patch to preserve symmetries
+        push!(ps, Tuple((g[i,j] + g[i+1,j] + g[i+1,j+1] + g[i,j+1])/4))
+        push!(values, (gf[i,j] + gf[i+1,j] + gf[i+1,j+1] + gf[i,j+1])/4)
+
+        push!(faces, (Li[i,j],     Li[i+1,j],   length(ps)))
+        push!(faces, (Li[i+1,j],   Li[i+1,j+1], length(ps)))
+        push!(faces, (Li[i+1,j+1], Li[i,j+1],   length(ps)))
+        push!(faces, (Li[i,j+1],   Li[i,j],     length(ps)))
+    end
+
+    verticies = permutedims(reinterpret(reshape,eltype(eltype(ps)), ps))
+    faces = permutedims(reinterpret(reshape,Int, faces))
+
+    return verticies, faces, values
+end
+
+
+## Grids
+
+Makie.convert_arguments(::Type{<:Scatter}, g::Grid) = (reshape(map(Point,g),:),) # (map(Point,collect(g)[:]),)
+function Makie.convert_arguments(::Type{<:Lines}, g::Grid{<:Any,2})
+    M = collect(g)
+
+    function cat_with_NaN(a,b)
+        vcat(a,[@SVector[NaN,NaN]],b)
+    end
+
+    xlines = reduce(cat_with_NaN, eachrow(M))
+    ylines = reduce(cat_with_NaN, eachcol(M))
+
+    return (cat_with_NaN(xlines,ylines),)
+end
+
+Makie.plot!(plot::Plot(Grid{<:Any,2})) = lines!(plot, plot.attributes, plot[1])
+
+
+## Grid functions
+
+### 1D
+function Makie.convert_arguments(::Type{<:Lines}, g::Grid{<:Any,1}, gf::AbstractArray{<:Any, 1})
+    (collect(g), gf)
+end
+
+function Makie.convert_arguments(::Type{<:Scatter}, g::Grid{<:Any,1}, gf::AbstractArray{<:Any, 1})
+    (collect(g), gf)
+end
+
+Makie.plot!(plot::Plot(Grid{<:Any,1}, AbstractArray{<:Any,1})) = lines!(plot, plot.attributes, plot[1], plot[2])
+
+### 2D
+function Makie.convert_arguments(::Type{<:Surface}, g::Grid{<:Any,2}, gf::AbstractArray{<:Any, 2})
+    (getindex.(g,1), getindex.(g,2), gf)
+end
+
+function Makie.plot!(plot::Plot(Grid{<:Any,2},AbstractArray{<:Any, 2}))
+    r = @lift verticies_and_faces_and_values($(plot[1]), $(plot[2]))
+    v,f,c = (@lift $r[1]), (@lift $r[2]), (@lift $r[3])
+    mesh!(plot, plot.attributes, v, f;
+        color=c,
+        shading = NoShading,
+    )
+end
+# TBD: Can we define `mesh` instead of the above function and then forward plot! to that?
+
+function Makie.convert_arguments(::Type{<:Scatter}, g::Grid{<:Any,2}, gf::AbstractArray{<:Any, 2})
+    ps = map(g,gf) do (x,y), z
+        @SVector[x,y,z]
+    end
+    (reshape(ps,:),)
+end
+
+end
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/ext/DiffinitiveSparseArrayKitExt.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -0,0 +1,22 @@
+module DiffinitiveSparseArrayKitExt
+
+using Diffinitive
+using Diffinitive.LazyTensors
+
+using SparseArrayKit
+using Tokens
+
+"""
+    SparseArray(t::LazyTensor)
+
+The sparse tensor representation of `t` with range dimensions to the left and
+domain dimensions to the right. If `L` is a `LazyTensor` with range and
+domain dimension 2 and `v` a 2-tensor, then `A = SparseArray(t)` is
+constructed so that `∑ₖ∑ₗA[i,j,k,l]*v[k,l] == L*v`.
+"""
+function SparseArrayKit.SparseArray(t::LazyTensor)
+    v = ArrayToken(:v, domain_size(t)...)
+    return Tokens._to_tensor(t*v, range_size(t), domain_size(t))
+end
+
+end
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/ext/DiffinitiveSparseArraysExt.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -0,0 +1,25 @@
+module DiffinitiveSparseArraysExt
+
+using Diffinitive
+using Diffinitive.LazyTensors
+
+using SparseArrays
+using Tokens
+
+"""
+    sparse(t::LazyTensor)
+
+The sparse matrix representation of `t`.
+
+If `L` is a `LazyTensor` and `v` a tensor, then `A = sparse(L)` is constructed
+so that `A*reshape(v,:) == reshape(L*v,:)`.
+"""
+function SparseArrays.sparse(t::LazyTensor)
+    v = ArrayToken(:v, prod(domain_size(t)))
+
+    v̄ = reshape(v,domain_size(t)...)
+    tv = reshape(t*v̄, :)
+    return Tokens._to_matrix(tv, prod(range_size(t)), prod(domain_size(t)))
+end
+
+end
--- a/ext/SbplibMakieExt.jl	Tue Sep 10 21:35:55 2024 +0200
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,88 +0,0 @@
-module SbplibMakieExt
-
-using Sbplib.Grids
-using Makie
-using StaticArrays
-
-
-function verticies_and_faces_and_values(g::Grid{<:Any,2}, gf::AbstractArray{<:Any, 2})
-    ps = map(Tuple, g)[:]
-    values = gf[:]
-    faces = Vector{NTuple{3,Int}}()
-
-    n,m = size(g)
-    Li = LinearIndices((1:n, 1:m))
-    for i ∈ 1:n-1, j = 1:m-1
-
-        # Add point in the middle of the patch to preserve symmetries
-        push!(ps, Tuple((g[i,j] + g[i+1,j] + g[i+1,j+1] + g[i,j+1])/4))
-        push!(values, (gf[i,j] + gf[i+1,j] + gf[i+1,j+1] + gf[i,j+1])/4)
-
-        push!(faces, (Li[i,j],     Li[i+1,j],   length(ps)))
-        push!(faces, (Li[i+1,j],   Li[i+1,j+1], length(ps)))
-        push!(faces, (Li[i+1,j+1], Li[i,j+1],   length(ps)))
-        push!(faces, (Li[i,j+1],   Li[i,j],     length(ps)))
-    end
-
-    verticies = permutedims(reinterpret(reshape,eltype(eltype(ps)), ps))
-    faces = permutedims(reinterpret(reshape,Int, faces))
-
-    return verticies, faces, values
-end
-
-
-## Grids
-
-Makie.convert_arguments(::Type{<:Scatter}, g::Grid) = (reshape(map(Point,g),:),) # (map(Point,collect(g)[:]),)
-function Makie.convert_arguments(::Type{<:Lines}, g::Grid{<:Any,2})
-    M = collect(g)
-
-    function cat_with_NaN(a,b)
-        vcat(a,[@SVector[NaN,NaN]],b)
-    end
-
-    xlines = reduce(cat_with_NaN, eachrow(M))
-    ylines = reduce(cat_with_NaN, eachcol(M))
-
-    return (cat_with_NaN(xlines,ylines),)
-end
-
-Makie.plot!(plot::Plot(Grid{<:Any,2})) = lines!(plot, plot.attributes, plot[1])
-
-
-## Grid functions
-
-### 1D
-function Makie.convert_arguments(::Type{<:Lines}, g::Grid{<:Any,1}, gf::AbstractArray{<:Any, 1})
-    (collect(g), gf)
-end
-
-function Makie.convert_arguments(::Type{<:Scatter}, g::Grid{<:Any,1}, gf::AbstractArray{<:Any, 1})
-    (collect(g), gf)
-end
-
-Makie.plot!(plot::Plot(Grid{<:Any,1}, AbstractArray{<:Any,1})) = lines!(plot, plot.attributes, plot[1], plot[2])
-
-### 2D
-function Makie.convert_arguments(::Type{<:Surface}, g::Grid{<:Any,2}, gf::AbstractArray{<:Any, 2})
-    (getindex.(g,1), getindex.(g,2), gf)
-end
-
-function Makie.plot!(plot::Plot(Grid{<:Any,2},AbstractArray{<:Any, 2}))
-    r = @lift verticies_and_faces_and_values($(plot[1]), $(plot[2]))
-    v,f,c = (@lift $r[1]), (@lift $r[2]), (@lift $r[3])
-    mesh!(plot, plot.attributes, v, f;
-        color=c,
-        shading = NoShading,
-    )
-end
-# TBD: Can we define `mesh` instead of the above function and then forward plot! to that?
-
-function Makie.convert_arguments(::Type{<:Scatter}, g::Grid{<:Any,2}, gf::AbstractArray{<:Any, 2})
-    ps = map(g,gf) do (x,y), z
-        @SVector[x,y,z]
-    end
-    (reshape(ps,:),)
-end
-
-end
--- a/sbp.jl	Tue Sep 10 21:35:55 2024 +0200
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,9 +0,0 @@
-module sbp
-
-using Sbplib.Grids
-using Sbplib.RegionIndices
-using Sbplib.SbpOperators
-using Sbplib.DiffOps
-
-include("TimeStepper.jl")
-end  # module
--- a/sbpPlot.jl	Tue Sep 10 21:35:55 2024 +0200
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,12 +0,0 @@
-include("sbp.jl")
-using Makie
-import .sbp.Grid
-function plotgridfunction(grid::sbp.Grid.EquidistantGrid, gridfunction::AbstractArray)
-    if sbp.Grid.dimension(grid) == 1
-        plot(sbp.Grid.pointsalongdim(grid,1), gridfunction)
-    elseif sbp.Grid.dimension(grid) == 2
-        scene = surface(sbp.Grid.pointsalongdim(grid,1),sbp.Grid.pointsalongdim(grid,2), gridfunction)
-    else
-        error(string("Plot not implemented for dimension ", string(dimension(grid))))
-    end
-end
--- a/src/DiffOps/DiffOps.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/src/DiffOps/DiffOps.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,9 +1,9 @@
 module DiffOps
 
-using Sbplib.RegionIndices
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
+using Diffinitive.RegionIndices
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
 
 """
     DiffOp
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/Diffinitive.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -0,0 +1,16 @@
+module Diffinitive
+
+include("StaticDicts/StaticDicts.jl")
+include("RegionIndices/RegionIndices.jl")
+include("LazyTensors/LazyTensors.jl")
+include("Grids/Grids.jl")
+include("SbpOperators/SbpOperators.jl")
+include("DiffOps/DiffOps.jl")
+
+export RegionIndices
+export LazyTensors
+export Grids
+export SbpOperators
+export DiffOps
+
+end
--- a/src/Grids/Grids.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/src/Grids/Grids.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,8 +1,7 @@
 # TODO: Double check that the interfaces for indexing and iterating are fully implemented and tested for all grids.
 module Grids
 
-using Sbplib.RegionIndices
-using Sbplib.LazyTensors
+using Diffinitive.LazyTensors
 using StaticArrays
 using LinearAlgebra
 
@@ -26,6 +25,8 @@
 export BoundaryIdentifier
 export TensorGridBoundary
 export CartesianBoundary
+export LowerBoundary
+export UpperBoundary
 
 export TensorGrid
 export ZeroDimGrid
--- a/src/Grids/equidistant_grid.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/src/Grids/equidistant_grid.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -49,12 +49,30 @@
 
 min_spacing(g::EquidistantGrid) = spacing(g)
 
+"""
+    LowerBoundary <: BoundaryIdentifier
 
-boundary_identifiers(::EquidistantGrid) = (Lower(), Upper())
-boundary_grid(g::EquidistantGrid, id::Lower) = ZeroDimGrid(g[begin])
-boundary_grid(g::EquidistantGrid, id::Upper) = ZeroDimGrid(g[end])
-boundary_indices(g::EquidistantGrid, id::Lower) = (1,)
-boundary_indices(g::EquidistantGrid, id::Upper) = (length(g),)
+Boundary identifier for the the lower (left) boundary of a one-dimensional grid.
+
+See also: [`BoundaryIdentifier`](@ref)
+"""
+struct LowerBoundary <: BoundaryIdentifier end
+
+"""
+    UpperBoundary <: BoundaryIdentifier
+
+Boundary identifier for the the upper (right)  boundary of a one-dimensional grid.
+
+See also: [`BoundaryIdentifier`](@ref)
+"""
+struct UpperBoundary <: BoundaryIdentifier end
+
+
+boundary_identifiers(::EquidistantGrid) = (LowerBoundary(), UpperBoundary())
+boundary_grid(g::EquidistantGrid, id::LowerBoundary) = ZeroDimGrid(g[begin])
+boundary_grid(g::EquidistantGrid, id::UpperBoundary) = ZeroDimGrid(g[end])
+boundary_indices(g::EquidistantGrid, id::LowerBoundary) = (firstindex(g),)
+boundary_indices(g::EquidistantGrid, id::UpperBoundary) = (lastindex(g),)
 
 """
     refine(g::EquidistantGrid, r::Int)
--- a/src/Grids/grid.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/src/Grids/grid.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -109,6 +109,13 @@
 function coarsen end
 
 """
+   BoundaryIdentifier
+
+An identifier for a boundary of a grid.
+"""
+abstract type BoundaryIdentifier end
+
+"""
     boundary_identifiers(g::Grid)
 
 Identifiers for all the boundaries of `g`.
--- a/src/Grids/mapped_grid.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/src/Grids/mapped_grid.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -140,11 +140,11 @@
 end
 
 function _boundary_sign(T, boundary)
-    if boundary_id(boundary) == Upper()
+    if boundary_id(boundary) == UpperBoundary()
         return one(T)
-    elseif boundary_id(boundary) == Lower()
+    elseif boundary_id(boundary) == LowerBoundary()
         return -one(T)
     else
-        throw(ArgumentError("The boundary identifier must be either `Lower()` or `Upper()`"))
+        throw(ArgumentError("The boundary identifier must be either `LowerBoundary()` or `UpperBoundary()`"))
     end
 end
--- a/src/Grids/tensor_grid.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/src/Grids/tensor_grid.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -117,7 +117,7 @@
 end
 
 """
-   grid_and_local_dim_index(nds, d)
+    grid_and_local_dim_index(nds, d)
 
 Given a tuple of number of dimensions `nds`, and a global dimension index `d`,
 calculate which grid index, and local dimension, `d` corresponds to.
--- a/src/SbpOperators/SbpOperators.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/src/SbpOperators/SbpOperators.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -40,9 +40,9 @@
 export sat_tensors
 
 # Using
-using Sbplib.RegionIndices
-using Sbplib.LazyTensors
-using Sbplib.Grids
+using Diffinitive.RegionIndices
+using Diffinitive.LazyTensors
+using Diffinitive.Grids
 
 # Includes
 include("stencil.jl")
--- a/src/SbpOperators/boundaryops/boundary_operator.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/src/SbpOperators/boundaryops/boundary_operator.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,27 +1,27 @@
 """
-    BoundaryOperator{T,R,N} <: LazyTensor{T,0,1}
+    BoundaryOperator{T,B,N} <: LazyTensor{T,0,1}
 
 Implements the boundary operator `op` for 1D as a `LazyTensor`
 
 `op` is the restriction of a grid function to the boundary using some closure
-`Stencil{T,N}`. The boundary to restrict to is determined by `R`. `op'` is the
+`Stencil{T,N}`. The boundary to restrict to is determined by `B`. `op'` is the
 prolongation of a zero dimensional array to the whole grid using the same
 closure stencil.
 """
-struct BoundaryOperator{T,R<:Region,N} <: LazyTensor{T,0,1}
+struct BoundaryOperator{T,B<:BoundaryIdentifier,N} <: LazyTensor{T,0,1}
     stencil::Stencil{T,N}
     size::Int
 end
 
 """
-    BoundaryOperator(grid::EquidistantGrid, closure_stencil, region)
+    BoundaryOperator(grid::EquidistantGrid, closure_stencil, boundary)
 
 Constructs the BoundaryOperator with stencil `closure_stencil` for a
 `EquidistantGrid` `grid`, restricting to to the boundary specified by
-`region`.
+`boundary`.
 """
-function BoundaryOperator(grid::EquidistantGrid, closure_stencil::Stencil{T,N}, region::Region) where {T,N}
-    return BoundaryOperator{T,typeof(region),N}(closure_stencil,size(grid)[1])
+function BoundaryOperator(grid::EquidistantGrid, closure_stencil::Stencil{T,N}, boundary::BoundaryIdentifier) where {T,N}
+    return BoundaryOperator{T,typeof(boundary),N}(closure_stencil,size(grid)[1])
 end
 
 """
@@ -29,24 +29,24 @@
 
 The size of the closure stencil.
 """
-closure_size(::BoundaryOperator{T,R,N}) where {T,R,N} = N
+closure_size(::BoundaryOperator{T,B,N}) where {T,B,N} = N
 
 LazyTensors.range_size(op::BoundaryOperator) = ()
 LazyTensors.domain_size(op::BoundaryOperator) = (op.size,)
 
-function LazyTensors.apply(op::BoundaryOperator{<:Any,Lower}, v::AbstractVector)
+function LazyTensors.apply(op::BoundaryOperator{<:Any,LowerBoundary}, v::AbstractVector)
     apply_stencil(op.stencil,v,1)
 end
 
-function LazyTensors.apply(op::BoundaryOperator{<:Any,Upper}, v::AbstractVector)
+function LazyTensors.apply(op::BoundaryOperator{<:Any,UpperBoundary}, v::AbstractVector)
     apply_stencil_backwards(op.stencil,v,op.size)
 end
 
-function LazyTensors.apply_transpose(op::BoundaryOperator{<:Any,Lower}, v::AbstractArray{<:Any,0}, i::Index{Lower})
+function LazyTensors.apply_transpose(op::BoundaryOperator{<:Any,LowerBoundary}, v::AbstractArray{<:Any,0}, i::Index{Lower})
     return op.stencil[Int(i)-1]*v[]
 end
 
-function LazyTensors.apply_transpose(op::BoundaryOperator{<:Any,Upper}, v::AbstractArray{<:Any,0}, i::Index{Upper})
+function LazyTensors.apply_transpose(op::BoundaryOperator{<:Any,UpperBoundary}, v::AbstractArray{<:Any,0}, i::Index{Upper})
     return op.stencil[op.size[1] - Int(i)]*v[]
 end
 
--- a/src/SbpOperators/volumeops/laplace/laplace.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/src/SbpOperators/volumeops/laplace/laplace.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -60,7 +60,7 @@
 condition. `H_tuning` and `R_tuning` are used to specify the strength of the
 penalty.
 
-See also: [`sat`](@ref),[`DirichletCondition`](@ref), [`positivity_decomposition`](@ref).
+See also: [`sat`](@ref), [`DirichletCondition`](@ref), [`positivity_decomposition`](@ref).
 """
 function sat_tensors(Δ::Laplace, g::Grid, bc::DirichletCondition; H_tuning = 1., R_tuning = 1.)
     id = boundary(bc)
@@ -69,7 +69,7 @@
     Hᵧ = inner_product(boundary_grid(g, id), set)
     e = boundary_restriction(g, set, id)
     d = normal_derivative(g, set, id)
-    B = positivity_decomposition(Δ, g, bc; H_tuning, R_tuning)
+    B = positivity_decomposition(Δ, g, boundary(bc); H_tuning, R_tuning)
     penalty_tensor = H⁻¹∘(d' - B*e')∘Hᵧ
     return penalty_tensor, e
 end
@@ -94,7 +94,7 @@
 end
 
 """
-    positivity_decomposition(Δ::Laplace, g::Grid, bc::DirichletCondition; H_tuning, R_tuning)
+    positivity_decomposition(Δ::Laplace, g::Grid, b::BoundaryIdentifier; H_tuning, R_tuning)
 
 Constructs the scalar `B` such that `d' - 1/2*B*e'` is symmetric positive
 definite with respect to the boundary quadrature. Here `d` is the normal
@@ -102,28 +102,26 @@
 to form a symmetric and energy stable penalty for a Dirichlet condition. The
 parameters `H_tuning` and `R_tuning` are used to specify the strength of the
 penalty and must be greater than 1. For details we refer to
-https://doi.org/10.1016/j.jcp.2020.109294
+<https://doi.org/10.1016/j.jcp.2020.109294>
 """
-function positivity_decomposition(Δ::Laplace, g::Grid, bc::DirichletCondition; H_tuning, R_tuning)
+function positivity_decomposition(Δ::Laplace, g::Grid, b::BoundaryIdentifier; H_tuning, R_tuning)
     @assert(H_tuning ≥ 1.)
     @assert(R_tuning ≥ 1.)
-    Nτ_H, τ_R = positivity_limits(Δ,g,bc)
+    Nτ_H, τ_R = positivity_limits(Δ,g,b)
     return H_tuning*Nτ_H + R_tuning*τ_R
 end
 
-# TODO: We should consider implementing a proper BoundaryIdentifier for EquidistantGrid and then
-# change bc::BoundaryCondition to id::BoundaryIdentifier
-function positivity_limits(Δ::Laplace, g::EquidistantGrid, bc::DirichletCondition)
+function positivity_limits(Δ::Laplace, g::EquidistantGrid, b::BoundaryIdentifier)
     h = spacing(g)
     θ_H = parse_scalar(Δ.stencil_set["H"]["closure"][1])
     θ_R = parse_scalar(Δ.stencil_set["D2"]["positivity"]["theta_R"])
 
-    τ_H = 1/(h*θ_H)
-    τ_R = 1/(h*θ_R)
+    τ_H = one(eltype(Δ))/(h*θ_H)
+    τ_R = one(eltype(Δ))/(h*θ_R)
     return τ_H, τ_R
 end
 
-function positivity_limits(Δ::Laplace, g::TensorGrid, bc::DirichletCondition)
-    τ_H, τ_R = positivity_limits(Δ, g.grids[grid_id(boundary(bc))], bc)
+function positivity_limits(Δ::Laplace, g::TensorGrid, b::BoundaryIdentifier)
+    τ_H, τ_R = positivity_limits(Δ, g.grids[grid_id(b)], b)
     return τ_H*ndims(g), τ_R
 end
--- a/src/Sbplib.jl	Tue Sep 10 21:35:55 2024 +0200
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,16 +0,0 @@
-module Sbplib
-
-include("StaticDicts/StaticDicts.jl")
-include("RegionIndices/RegionIndices.jl")
-include("LazyTensors/LazyTensors.jl")
-include("Grids/Grids.jl")
-include("SbpOperators/SbpOperators.jl")
-include("DiffOps/DiffOps.jl")
-
-export RegionIndices
-export LazyTensors
-export Grids
-export SbpOperators
-export DiffOps
-
-end
--- a/test/DiffOps/DiffOps_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/DiffOps/DiffOps_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,9 +1,9 @@
 using Test
-using Sbplib.DiffOps
-using Sbplib.Grids
-using Sbplib.SbpOperators
-using Sbplib.RegionIndices
-using Sbplib.LazyTensors
+using Diffinitive.DiffOps
+using Diffinitive.Grids
+using Diffinitive.SbpOperators
+using Diffinitive.RegionIndices
+using Diffinitive.LazyTensors
 
 #
 # @testset "BoundaryValue" begin
--- a/test/Grids/equidistant_grid_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/Grids/equidistant_grid_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,7 +1,6 @@
-using Sbplib.Grids
+using Diffinitive.Grids
 using Test
-using Sbplib.RegionIndices
-using Sbplib.LazyTensors
+using Diffinitive.LazyTensors
 
 
 @testset "EquidistantGrid" begin
@@ -63,24 +62,24 @@
 
     @testset "boundary_identifiers" begin
         g = EquidistantGrid(0:0.1:10)
-        @test boundary_identifiers(g) == (Lower(), Upper())
+        @test boundary_identifiers(g) == (LowerBoundary(), UpperBoundary())
         @inferred boundary_identifiers(g)
     end
 
     @testset "boundary_grid" begin
         g = EquidistantGrid(0:0.1:1)
-        @test boundary_grid(g, Lower()) == ZeroDimGrid(0.0)
-        @test boundary_grid(g, Upper()) == ZeroDimGrid(1.0)
+        @test boundary_grid(g, LowerBoundary()) == ZeroDimGrid(0.0)
+        @test boundary_grid(g, UpperBoundary()) == ZeroDimGrid(1.0)
     end
 
     @testset "boundary_indices" begin
         g = EquidistantGrid(0:0.1:1)
-        @test boundary_indices(g, Lower()) == (1,)
-        @test boundary_indices(g, Upper()) == (11,)
+        @test boundary_indices(g, LowerBoundary()) == (1,)
+        @test boundary_indices(g, UpperBoundary()) == (11,)
 
         g = EquidistantGrid(2:0.1:10)
-        @test boundary_indices(g, Lower()) == (1,)
-        @test boundary_indices(g, Upper()) == (81,)
+        @test boundary_indices(g, LowerBoundary()) == (1,)
+        @test boundary_indices(g, UpperBoundary()) == (81,)
 
     end
 
--- a/test/Grids/grid_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/Grids/grid_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,6 +1,6 @@
 using Test
-using Sbplib.Grids
-using Sbplib.LazyTensors
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
 using StaticArrays
 
 @testset "Grid" begin
--- a/test/Grids/mapped_grid_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/Grids/mapped_grid_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,5 +1,5 @@
-using Sbplib.Grids
-using Sbplib.RegionIndices
+using Diffinitive.Grids
+using Diffinitive.RegionIndices
 using Test
 using StaticArrays
 using LinearAlgebra
@@ -173,9 +173,9 @@
         J = map(ξ̄ -> @SArray(fill(2., 2, 2)), lg)
         mg = MappedGrid(lg, x̄, J)
 
-        @test boundary_indices(mg, CartesianBoundary{1,Lower}()) == boundary_indices(lg,CartesianBoundary{1,Lower}())
-        @test boundary_indices(mg, CartesianBoundary{2,Lower}()) == boundary_indices(lg,CartesianBoundary{2,Lower}())
-        @test boundary_indices(mg, CartesianBoundary{1,Upper}()) == boundary_indices(lg,CartesianBoundary{1,Upper}())
+        @test boundary_indices(mg, CartesianBoundary{1,LowerBoundary}()) == boundary_indices(lg,CartesianBoundary{1,LowerBoundary}())
+        @test boundary_indices(mg, CartesianBoundary{2,LowerBoundary}()) == boundary_indices(lg,CartesianBoundary{2,LowerBoundary}())
+        @test boundary_indices(mg, CartesianBoundary{1,UpperBoundary}()) == boundary_indices(lg,CartesianBoundary{1,UpperBoundary}())
     end
 
     @testset "boundary_grid" begin
@@ -208,10 +208,10 @@
             end
         end
 
-        @testset test_boundary_grid(mg, TensorGridBoundary{1, Lower}(), J2)
-        @testset test_boundary_grid(mg, TensorGridBoundary{1, Upper}(), J2)
-        @testset test_boundary_grid(mg, TensorGridBoundary{2, Lower}(), J1)
-        @testset test_boundary_grid(mg, TensorGridBoundary{2, Upper}(), J1)
+        @testset test_boundary_grid(mg, TensorGridBoundary{1, LowerBoundary}(), J2)
+        @testset test_boundary_grid(mg, TensorGridBoundary{1, UpperBoundary}(), J2)
+        @testset test_boundary_grid(mg, TensorGridBoundary{2, LowerBoundary}(), J1)
+        @testset test_boundary_grid(mg, TensorGridBoundary{2, UpperBoundary}(), J1)
     end
 end
 
@@ -277,10 +277,10 @@
 @testset "normal" begin
     g = mapped_grid(_partially_curved_mapping()...,10, 11)
 
-    @test normal(g, CartesianBoundary{1,Lower}()) == fill(@SVector[-1,0], 11)
-    @test normal(g, CartesianBoundary{1,Upper}()) == fill(@SVector[1,0], 11)
-    @test normal(g, CartesianBoundary{2,Lower}()) == fill(@SVector[0,-1], 10)
-    @test normal(g, CartesianBoundary{2,Upper}()) ≈ map(boundary_grid(g,CartesianBoundary{2,Upper}())|>logicalgrid) do ξ̄
+    @test normal(g, CartesianBoundary{1,LowerBoundary}()) == fill(@SVector[-1,0], 11)
+    @test normal(g, CartesianBoundary{1,UpperBoundary}()) == fill(@SVector[1,0], 11)
+    @test normal(g, CartesianBoundary{2,LowerBoundary}()) == fill(@SVector[0,-1], 10)
+    @test normal(g, CartesianBoundary{2,UpperBoundary}()) ≈ map(boundary_grid(g,CartesianBoundary{2,UpperBoundary}())|>logicalgrid) do ξ̄
         α = 1-2ξ̄[1]
         @SVector[α,1]/√(α^2 + 1)
     end
@@ -288,28 +288,28 @@
     g = mapped_grid(_fully_curved_mapping()...,5,4)
 
     unit(v) = v/norm(v)
-    @testset let bId = CartesianBoundary{1,Lower}()
+    @testset let bId = CartesianBoundary{1,LowerBoundary}()
         lbg = boundary_grid(logicalgrid(g), bId)
         @test normal(g, bId) ≈ map(lbg) do (ξ, η)
             -unit(@SVector[1/2,  η/3-1/6])
         end
     end
 
-    @testset let bId = CartesianBoundary{1,Upper}()
+    @testset let bId = CartesianBoundary{1,UpperBoundary}()
         lbg = boundary_grid(logicalgrid(g), bId)
         @test normal(g, bId) ≈ map(lbg) do (ξ, η)
             unit(@SVector[7/2, 2η-1]/(5 + 3η + 2η^2))
         end
     end
 
-    @testset let bId = CartesianBoundary{2,Lower}()
+    @testset let bId = CartesianBoundary{2,LowerBoundary}()
         lbg = boundary_grid(logicalgrid(g), bId)
         @test normal(g, bId) ≈ map(lbg) do (ξ, η)
             -unit(@SVector[-2ξ, 2]/(6 + ξ^2 - 2ξ))
         end
     end
 
-    @testset let bId = CartesianBoundary{2,Upper}()
+    @testset let bId = CartesianBoundary{2,UpperBoundary}()
         lbg = boundary_grid(logicalgrid(g), bId)
         @test normal(g, bId) ≈ map(lbg) do (ξ, η)
             unit(@SVector[-3ξ, 2]/(6 + ξ^2 + 3ξ))
--- a/test/Grids/tensor_grid_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/Grids/tensor_grid_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,7 +1,6 @@
 using Test
-using Sbplib.Grids
+using Diffinitive.Grids
 using StaticArrays
-using Sbplib.RegionIndices
 
 @testset "TensorGrid" begin
     g₁ = EquidistantGrid(range(0,1,length=11))
@@ -170,13 +169,13 @@
     end
 
     @testset "boundary_identifiers" begin
-        @test boundary_identifiers(TensorGrid(g₁, g₂)) == map((n,id)->TensorGridBoundary{n,id}(), (1,1,2,2), (Lower,Upper,Lower,Upper))
-        @test boundary_identifiers(TensorGrid(g₁, g₄)) == (TensorGridBoundary{1,Lower}(),TensorGridBoundary{1,Upper}())
+        @test boundary_identifiers(TensorGrid(g₁, g₂)) == map((n,id)->TensorGridBoundary{n,id}(), (1,1,2,2), (LowerBoundary,UpperBoundary,LowerBoundary,UpperBoundary))
+        @test boundary_identifiers(TensorGrid(g₁, g₄)) == (TensorGridBoundary{1,LowerBoundary}(),TensorGridBoundary{1,UpperBoundary}())
     end
 
     @testset "boundary_grid" begin
-        @test boundary_grid(TensorGrid(g₁, g₂), TensorGridBoundary{1, Upper}()) == TensorGrid(ZeroDimGrid(g₁[end]), g₂)
-        @test boundary_grid(TensorGrid(g₁, g₄), TensorGridBoundary{1, Upper}()) == TensorGrid(ZeroDimGrid(g₁[end]), g₄)
+        @test boundary_grid(TensorGrid(g₁, g₂), TensorGridBoundary{1, UpperBoundary}()) == TensorGrid(ZeroDimGrid(g₁[end]), g₂)
+        @test boundary_grid(TensorGrid(g₁, g₄), TensorGridBoundary{1, UpperBoundary}()) == TensorGrid(ZeroDimGrid(g₁[end]), g₄)
     end
 
     @testset "boundary_indices" begin
@@ -184,14 +183,14 @@
         g₂ = EquidistantGrid(range(2,3,length=6))
         g₄ = ZeroDimGrid(@SVector[1,2])
 
-        @test boundary_indices(TensorGrid(g₁, g₂), TensorGridBoundary{1, Lower}()) == (1,:)
-        @test boundary_indices(TensorGrid(g₁, g₂), TensorGridBoundary{1, Upper}()) == (11,:)
-        @test boundary_indices(TensorGrid(g₁, g₂), TensorGridBoundary{2, Lower}()) == (:,1)
-        @test boundary_indices(TensorGrid(g₁, g₂), TensorGridBoundary{2, Upper}()) == (:,6)
-        @test boundary_indices(TensorGrid(g₁, g₄), TensorGridBoundary{1, Lower}()) == (1,)
-        @test boundary_indices(TensorGrid(g₁, g₄), TensorGridBoundary{1, Upper}()) == (11,)
-        @test boundary_indices(TensorGrid(g₄,g₁), TensorGridBoundary{2, Lower}()) == (1,)
-        @test boundary_indices(TensorGrid(g₄,g₁), TensorGridBoundary{2, Upper}()) == (11,)
+        @test boundary_indices(TensorGrid(g₁, g₂), TensorGridBoundary{1, LowerBoundary}()) == (1,:)
+        @test boundary_indices(TensorGrid(g₁, g₂), TensorGridBoundary{1, UpperBoundary}()) == (11,:)
+        @test boundary_indices(TensorGrid(g₁, g₂), TensorGridBoundary{2, LowerBoundary}()) == (:,1)
+        @test boundary_indices(TensorGrid(g₁, g₂), TensorGridBoundary{2, UpperBoundary}()) == (:,6)
+        @test boundary_indices(TensorGrid(g₁, g₄), TensorGridBoundary{1, LowerBoundary}()) == (1,)
+        @test boundary_indices(TensorGrid(g₁, g₄), TensorGridBoundary{1, UpperBoundary}()) == (11,)
+        @test boundary_indices(TensorGrid(g₄,g₁), TensorGridBoundary{2, LowerBoundary}()) == (1,)
+        @test boundary_indices(TensorGrid(g₄,g₁), TensorGridBoundary{2, UpperBoundary}()) == (11,)
     end
 end
 
--- a/test/Grids/zero_dim_grid_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/Grids/zero_dim_grid_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,5 +1,5 @@
 using Test
-using Sbplib.Grids
+using Diffinitive.Grids
 using StaticArrays
 
 @testset "ZeroDimGrid" begin
--- a/test/LazyTensors/lazy_array_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/LazyTensors/lazy_array_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,6 +1,6 @@
 using Test
-using Sbplib.LazyTensors
-using Sbplib.RegionIndices
+using Diffinitive.LazyTensors
+using Diffinitive.RegionIndices
 
 
 @testset "LazyArray" begin
--- a/test/LazyTensors/lazy_tensor_operations_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/LazyTensors/lazy_tensor_operations_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,6 +1,6 @@
 using Test
-using Sbplib.LazyTensors
-using Sbplib.RegionIndices
+using Diffinitive.LazyTensors
+using Diffinitive.RegionIndices
 
 using Tullio
 
--- a/test/LazyTensors/lazy_tensor_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/LazyTensors/lazy_tensor_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,5 +1,5 @@
 using Test
-using Sbplib.LazyTensors
+using Diffinitive.LazyTensors
 
 @testset "Generic Mapping methods" begin
     struct DummyMapping{T,R,D} <: LazyTensor{T,R,D} end
--- a/test/LazyTensors/tensor_types_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/LazyTensors/tensor_types_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,5 +1,5 @@
 using Test
-using Sbplib.LazyTensors
+using Diffinitive.LazyTensors
 using BenchmarkTools
 
 @testset "IdentityTensor" begin
--- a/test/LazyTensors/tuple_manipulation_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/LazyTensors/tuple_manipulation_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,5 +1,5 @@
 using Test
-using Sbplib.LazyTensors
+using Diffinitive.LazyTensors
 
 @testset "split_index" begin
     @test LazyTensors.split_index(2,1,2,2, 1,2,3,4,5,6) == ((1,2,:,5,6),(3,4))
--- a/test/Manifest.toml	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/Manifest.toml	Tue Sep 10 21:59:10 2024 +0200
@@ -1,8 +1,8 @@
 # This file is machine-generated - editing it directly is not advised
 
-julia_version = "1.10.4"
+julia_version = "1.10.5"
 manifest_format = "2.0"
-project_hash = "f2b0634c12bbed93a17efc88d466604d5a07c465"
+project_hash = "9dddd5385164ee197d1b3f22302bc95701c1f5e5"
 
 [[deps.ArgTools]]
 uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f"
@@ -74,9 +74,9 @@
 
 [[deps.JLLWrappers]]
 deps = ["Artifacts", "Preferences"]
-git-tree-sha1 = "7e5d6779a1e09a36db2a7b6cff50942a0a7d0fca"
+git-tree-sha1 = "f389674c99bfcde17dc57454011aa44d5a260a40"
 uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210"
-version = "1.5.0"
+version = "1.6.0"
 
 [[deps.JSON]]
 deps = ["Dates", "Mmap", "Parsers", "Unicode"]
@@ -176,6 +176,12 @@
 uuid = "efe28fd5-8261-553b-a9e1-b2916fc3738e"
 version = "0.5.5+0"
 
+[[deps.PackageExtensionCompat]]
+git-tree-sha1 = "fb28e33b8a95c4cee25ce296c817d89cc2e53518"
+uuid = "65ce6f38-6b18-4e1d-a461-8949797d7930"
+version = "1.0.2"
+weakdeps = ["Requires", "TOML"]
+
 [[deps.Parsers]]
 deps = ["Dates", "PrecompileTools", "UUIDs"]
 git-tree-sha1 = "8489905bcdbcfac64d1daa51ca07c0d8f0283821"
@@ -231,6 +237,20 @@
 [[deps.Sockets]]
 uuid = "6462fe0b-24de-5631-8697-dd941f90decc"
 
+[[deps.SparseArrayKit]]
+deps = ["LinearAlgebra", "PackageExtensionCompat", "TupleTools", "VectorInterface"]
+git-tree-sha1 = "e06f75c460c16b5b08bf93324b7db2cf9c1fe831"
+uuid = "a9a3c162-d163-4c15-8926-b8794fbefed2"
+version = "0.3.1"
+
+    [deps.SparseArrayKit.extensions]
+    SparseArrayKitSparseArrays = "SparseArrays"
+    SparseArrayKitTensorOperations = "TensorOperations"
+
+    [deps.SparseArrayKit.weakdeps]
+    SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+    TensorOperations = "6aa20fa7-93e2-5fca-9bc0-fbd0db3c71a2"
+
 [[deps.SparseArrays]]
 deps = ["Libdl", "LinearAlgebra", "Random", "Serialization", "SuiteSparse_jll"]
 uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
@@ -297,6 +317,12 @@
 uuid = "98d24dd4-01ad-11ea-1b02-c9a08f80db04"
 version = "3.0.0"
 
+[[deps.Tokens]]
+deps = ["SparseArrayKit", "SparseArrays"]
+git-tree-sha1 = "c4f40125383ce3bfcfcd49a1b206080b7afd9a34"
+uuid = "040c2ec2-8d69-4aca-bf03-7d3a7092f2f6"
+version = "0.1.1"
+
 [[deps.Tullio]]
 deps = ["DiffRules", "LinearAlgebra", "Requires"]
 git-tree-sha1 = "6d476962ba4e435d7f4101a403b1d3d72afe72f3"
@@ -315,6 +341,11 @@
     FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b"
     Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c"
 
+[[deps.TupleTools]]
+git-tree-sha1 = "41d61b1c545b06279871ef1a4b5fcb2cac2191cd"
+uuid = "9d95972d-f1c8-5527-a6e0-b4b365fa01f6"
+version = "1.5.0"
+
 [[deps.UUIDs]]
 deps = ["Random", "SHA"]
 uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4"
@@ -322,6 +353,12 @@
 [[deps.Unicode]]
 uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5"
 
+[[deps.VectorInterface]]
+deps = ["LinearAlgebra"]
+git-tree-sha1 = "7aff7d62bffad9bba9928eb6ab55226b32a351eb"
+uuid = "409d34a3-91d5-4945-b6ec-7529ddf182d8"
+version = "0.4.6"
+
 [[deps.Zlib_jll]]
 deps = ["Libdl"]
 uuid = "83775a58-1f1d-513f-b197-d71354ab007a"
@@ -330,7 +367,7 @@
 [[deps.libblastrampoline_jll]]
 deps = ["Artifacts", "Libdl"]
 uuid = "8e850b90-86db-534c-a0d3-1478176c7d93"
-version = "5.8.0+1"
+version = "5.11.0+0"
 
 [[deps.nghttp2_jll]]
 deps = ["Artifacts", "Libdl"]
--- a/test/Project.toml	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/Project.toml	Tue Sep 10 21:59:10 2024 +0200
@@ -2,8 +2,11 @@
 BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
 Glob = "c27321d9-0574-5035-807b-f59d2c89b15c"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+SparseArrayKit = "a9a3c162-d163-4c15-8926-b8794fbefed2"
+SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 TOML = "fa267f1f-6049-4f14-aa54-33bafae1ed76"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 TestSetExtensions = "98d24dd4-01ad-11ea-1b02-c9a08f80db04"
+Tokens = "040c2ec2-8d69-4aca-bf03-7d3a7092f2f6"
 Tullio = "bc48ee85-29a4-5162-ae0b-a64e1601d4bc"
--- a/test/RegionIndices/RegionIndices_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/RegionIndices/RegionIndices_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,3 +1,3 @@
-using Sbplib.RegionIndices
+using Diffinitive.RegionIndices
 using Test
 
--- a/test/SbpOperators/boundary_conditions/boundary_condition_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/boundary_conditions/boundary_condition_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,8 +1,8 @@
 using Test
 
-using Sbplib.Grids
-using Sbplib.RegionIndices
-using Sbplib.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.RegionIndices
+using Diffinitive.SbpOperators
 
 @testset "BoundaryCondition" begin
     grid_1d = equidistant_grid(0.0, 1.0, 11)
@@ -15,8 +15,8 @@
     g = 3.14
     f(x,y,z) = x^2+y^2+z^2
     @testset "Constructors" begin
-        @test DirichletCondition(g,id_l) isa DirichletCondition{Float64,Lower}
-        @test NeumannCondition(f,id_b) isa NeumannCondition{<:Function,CartesianBoundary{3,Lower}}
+        @test DirichletCondition(g,id_l) isa DirichletCondition{Float64,LowerBoundary}
+        @test NeumannCondition(f,id_b) isa NeumannCondition{<:Function,CartesianBoundary{3,LowerBoundary}}
     end
 
     @testset "boundary" begin
--- a/test/SbpOperators/boundary_conditions/sat_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/boundary_conditions/sat_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,8 +1,8 @@
 using Test
 
-using Sbplib.Grids
-using Sbplib.LazyTensors
-using Sbplib.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
+using Diffinitive.SbpOperators
 
 stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order = 4)
 
--- a/test/SbpOperators/boundaryops/boundary_operator_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/boundaryops/boundary_operator_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,11 +1,11 @@
 using Test
 
-using Sbplib.LazyTensors
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.RegionIndices
-import Sbplib.SbpOperators.Stencil
-import Sbplib.SbpOperators.BoundaryOperator
+using Diffinitive.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.RegionIndices
+import Diffinitive.SbpOperators.Stencil
+import Diffinitive.SbpOperators.BoundaryOperator
 
 
 @testset "BoundaryOperator" begin
@@ -13,12 +13,12 @@
     g_1D = EquidistantGrid(range(0,1,length=11))
 
     @testset "Constructors" begin
-        @test BoundaryOperator(g_1D, closure_stencil, Lower()) isa LazyTensor{T,0,1} where T
-        @test BoundaryOperator(g_1D, closure_stencil, Upper()) isa LazyTensor{T,0,1} where T
+        @test BoundaryOperator(g_1D, closure_stencil, LowerBoundary()) isa LazyTensor{T,0,1} where T
+        @test BoundaryOperator(g_1D, closure_stencil, UpperBoundary()) isa LazyTensor{T,0,1} where T
     end
 
-    op_l = BoundaryOperator(g_1D, closure_stencil, Lower())
-    op_r = BoundaryOperator(g_1D, closure_stencil, Upper())
+    op_l = BoundaryOperator(g_1D, closure_stencil, LowerBoundary())
+    op_r = BoundaryOperator(g_1D, closure_stencil, UpperBoundary())
 
     @testset "Sizes" begin
         @test domain_size(op_l) == (11,)
--- a/test/SbpOperators/boundaryops/boundary_restriction_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/boundaryops/boundary_restriction_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,10 @@
 using Test
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
-using Sbplib.RegionIndices
-using Sbplib.SbpOperators: BoundaryOperator, Stencil
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
+using Diffinitive.RegionIndices
+using Diffinitive.SbpOperators: BoundaryOperator, Stencil
 
 @testset "boundary_restriction" begin
 	stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order = 4)
@@ -14,19 +14,19 @@
 
     @testset "boundary_restriction" begin
         @testset "1D" begin
-            e_l = boundary_restriction(g_1D,stencil_set,Lower())
-            @test e_l == BoundaryOperator(g_1D,Stencil{Float64}(e_closure),Lower())
-            @test e_l isa BoundaryOperator{T,Lower} where T
+            e_l = boundary_restriction(g_1D,stencil_set,LowerBoundary())
+            @test e_l == BoundaryOperator(g_1D,Stencil{Float64}(e_closure),LowerBoundary())
+            @test e_l isa BoundaryOperator{T,LowerBoundary} where T
             @test e_l isa LazyTensor{T,0,1} where T
 
-            e_r = boundary_restriction(g_1D,stencil_set,Upper())
-            @test e_r == BoundaryOperator(g_1D,Stencil{Float64}(e_closure),Upper())
-            @test e_r isa BoundaryOperator{T,Upper} where T
+            e_r = boundary_restriction(g_1D,stencil_set,UpperBoundary())
+            @test e_r == BoundaryOperator(g_1D,Stencil{Float64}(e_closure),UpperBoundary())
+            @test e_r isa BoundaryOperator{T,UpperBoundary} where T
             @test e_r isa LazyTensor{T,0,1} where T
         end
 
         @testset "2D" begin
-            e_w = boundary_restriction(g_2D,stencil_set,CartesianBoundary{1,Upper}())
+            e_w = boundary_restriction(g_2D,stencil_set,CartesianBoundary{1,UpperBoundary}())
             @test e_w isa InflatedTensor
             @test e_w isa LazyTensor{T,1,2} where T
         end
--- a/test/SbpOperators/boundaryops/normal_derivative_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/boundaryops/normal_derivative_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,10 @@
 using Test
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
-using Sbplib.RegionIndices
-import Sbplib.SbpOperators.BoundaryOperator
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
+using Diffinitive.RegionIndices
+import Diffinitive.SbpOperators.BoundaryOperator
 
 @testset "normal_derivative" begin
     g_1D = equidistant_grid(0.0, 1.0, 11)
@@ -12,19 +12,19 @@
     @testset "normal_derivative" begin
     	stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=4)
         @testset "1D" begin
-            d_l = normal_derivative(g_1D, stencil_set, Lower())
-            @test d_l == normal_derivative(g_1D, stencil_set, Lower())
-            @test d_l isa BoundaryOperator{T,Lower} where T
+            d_l = normal_derivative(g_1D, stencil_set, LowerBoundary())
+            @test d_l == normal_derivative(g_1D, stencil_set, LowerBoundary())
+            @test d_l isa BoundaryOperator{T,LowerBoundary} where T
             @test d_l isa LazyTensor{T,0,1} where T
         end
         @testset "2D" begin
-            d_w = normal_derivative(g_2D, stencil_set, CartesianBoundary{1,Lower}())
-            d_n = normal_derivative(g_2D, stencil_set, CartesianBoundary{2,Upper}())
+            d_w = normal_derivative(g_2D, stencil_set, CartesianBoundary{1,LowerBoundary}())
+            d_n = normal_derivative(g_2D, stencil_set, CartesianBoundary{2,UpperBoundary}())
             Ix = IdentityTensor{Float64}((size(g_2D)[1],))
             Iy = IdentityTensor{Float64}((size(g_2D)[2],))
-            d_l = normal_derivative(g_2D.grids[1], stencil_set, Lower())
-            d_r = normal_derivative(g_2D.grids[2], stencil_set, Upper())
-            @test d_w == normal_derivative(g_2D, stencil_set, CartesianBoundary{1,Lower}())
+            d_l = normal_derivative(g_2D.grids[1], stencil_set, LowerBoundary())
+            d_r = normal_derivative(g_2D.grids[2], stencil_set, UpperBoundary())
+            @test d_w == normal_derivative(g_2D, stencil_set, CartesianBoundary{1,LowerBoundary}())
             @test d_w ==  d_l⊗Iy
             @test d_n ==  Ix⊗d_r
             @test d_w isa LazyTensor{T,1,2} where T
--- a/test/SbpOperators/stencil_set_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/stencil_set_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,10 @@
 using Test
 
 using TOML
-using Sbplib.SbpOperators
+using Diffinitive.SbpOperators
 
-import Sbplib.SbpOperators.Stencil
-import Sbplib.SbpOperators.NestedStencil
+import Diffinitive.SbpOperators.Stencil
+import Diffinitive.SbpOperators.NestedStencil
 
 @testset "readoperator" begin
     toml_str = """
--- a/test/SbpOperators/stencil_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/stencil_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,10 @@
 using Test
-using Sbplib.SbpOperators
+using Diffinitive.SbpOperators
 using StaticArrays
-import Sbplib.SbpOperators.Stencil
-import Sbplib.SbpOperators.NestedStencil
-import Sbplib.SbpOperators.scale
-import Sbplib.SbpOperators: apply_stencil, apply_stencil_backwards
+import Diffinitive.SbpOperators.Stencil
+import Diffinitive.SbpOperators.NestedStencil
+import Diffinitive.SbpOperators.scale
+import Diffinitive.SbpOperators: apply_stencil, apply_stencil_backwards
 
 @testset "Stencil" begin
     s = Stencil(-2:2, (1.,2.,2.,3.,4.))
--- a/test/SbpOperators/volumeops/constant_interior_scaling_operator_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/constant_interior_scaling_operator_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,9 +1,9 @@
 using Test
 
-using Sbplib.LazyTensors
-using Sbplib.SbpOperators
-import Sbplib.SbpOperators: ConstantInteriorScalingOperator
-using Sbplib.Grids
+using Diffinitive.LazyTensors
+using Diffinitive.SbpOperators
+import Diffinitive.SbpOperators: ConstantInteriorScalingOperator
+using Diffinitive.Grids
 
 @testset "ConstantInteriorScalingOperator" begin
     @test ConstantInteriorScalingOperator(1, (2,3), 10) isa ConstantInteriorScalingOperator{Int,2}
--- a/test/SbpOperators/volumeops/derivatives/dissipation_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/derivatives/dissipation_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,17 +1,17 @@
 using Test
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
 
-using Sbplib.SbpOperators: Stencil
+using Diffinitive.SbpOperators: Stencil
 
-using Sbplib.SbpOperators: dissipation_interior_weights
-using Sbplib.SbpOperators: dissipation_interior_stencil, dissipation_transpose_interior_stencil
-using Sbplib.SbpOperators: midpoint, midpoint_transpose
-using Sbplib.SbpOperators: dissipation_lower_closure_size, dissipation_upper_closure_size
-using Sbplib.SbpOperators: dissipation_lower_closure_stencils,dissipation_upper_closure_stencils
-using Sbplib.SbpOperators: dissipation_transpose_lower_closure_stencils, dissipation_transpose_upper_closure_stencils
+using Diffinitive.SbpOperators: dissipation_interior_weights
+using Diffinitive.SbpOperators: dissipation_interior_stencil, dissipation_transpose_interior_stencil
+using Diffinitive.SbpOperators: midpoint, midpoint_transpose
+using Diffinitive.SbpOperators: dissipation_lower_closure_size, dissipation_upper_closure_size
+using Diffinitive.SbpOperators: dissipation_lower_closure_stencils,dissipation_upper_closure_stencils
+using Diffinitive.SbpOperators: dissipation_transpose_lower_closure_stencils, dissipation_transpose_upper_closure_stencils
 
 
 @testset "undivided_skewed04" begin
--- a/test/SbpOperators/volumeops/derivatives/first_derivative_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/derivatives/first_derivative_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,11 +1,11 @@
 using Test
 
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
 
-using Sbplib.SbpOperators: closure_size, Stencil, VolumeOperator
+using Diffinitive.SbpOperators: closure_size, Stencil, VolumeOperator
 
 @testset "first_derivative" begin
     @testset "Constructors" begin
--- a/test/SbpOperators/volumeops/derivatives/second_derivative_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/derivatives/second_derivative_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,10 @@
 using Test
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
 
-import Sbplib.SbpOperators.VolumeOperator
+import Diffinitive.SbpOperators.VolumeOperator
 
 # TODO: Refactor these test to look more like the tests in first_derivative_test.jl.
 
--- a/test/SbpOperators/volumeops/derivatives/second_derivative_variable_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/derivatives/second_derivative_variable_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,10 @@
 using Test
 
-using Sbplib.Grids
-using Sbplib.LazyTensors
-using Sbplib.SbpOperators
-using Sbplib.RegionIndices
-using Sbplib.SbpOperators: NestedStencil, CenteredNestedStencil, SecondDerivativeVariable
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.RegionIndices
+using Diffinitive.SbpOperators: NestedStencil, CenteredNestedStencil, SecondDerivativeVariable
 
 using LinearAlgebra
 
--- a/test/SbpOperators/volumeops/inner_products/inner_product_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/inner_products/inner_product_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,10 @@
 using Test
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
 
-import Sbplib.SbpOperators.ConstantInteriorScalingOperator
+import Diffinitive.SbpOperators.ConstantInteriorScalingOperator
 
 @testset "Diagonal-stencil inner_product" begin
     Lx = π/2.
--- a/test/SbpOperators/volumeops/inner_products/inverse_inner_product_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/inner_products/inverse_inner_product_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,10 +1,10 @@
 using Test
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
 
-import Sbplib.SbpOperators.ConstantInteriorScalingOperator
+import Diffinitive.SbpOperators.ConstantInteriorScalingOperator
 
 @testset "Diagonal-stencil inverse_inner_product" begin
     Lx = π/2.
--- a/test/SbpOperators/volumeops/laplace/laplace_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/laplace/laplace_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,8 +1,8 @@
 using Test
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
 
 @testset "Laplace" begin
     # Default stencils (4th order)
--- a/test/SbpOperators/volumeops/stencil_operator_distinct_closures_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/stencil_operator_distinct_closures_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,11 +1,11 @@
 using Test
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.LazyTensors
 
-import Sbplib.SbpOperators.Stencil
-import Sbplib.SbpOperators.StencilOperatorDistinctClosures
+import Diffinitive.SbpOperators.Stencil
+import Diffinitive.SbpOperators.StencilOperatorDistinctClosures
 
 @testset "StencilOperatorDistinctClosures" begin
     g = equidistant_grid(0., 1., 11)
--- a/test/SbpOperators/volumeops/volume_operator_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/SbpOperators/volumeops/volume_operator_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,14 +1,14 @@
 using Test
 
-using Sbplib.SbpOperators
-using Sbplib.Grids
-using Sbplib.RegionIndices
-using Sbplib.LazyTensors
+using Diffinitive.SbpOperators
+using Diffinitive.Grids
+using Diffinitive.RegionIndices
+using Diffinitive.LazyTensors
 
-import Sbplib.SbpOperators.Stencil
-import Sbplib.SbpOperators.VolumeOperator
-import Sbplib.SbpOperators.odd
-import Sbplib.SbpOperators.even
+import Diffinitive.SbpOperators.Stencil
+import Diffinitive.SbpOperators.VolumeOperator
+import Diffinitive.SbpOperators.odd
+import Diffinitive.SbpOperators.even
 
 
 @testset "VolumeOperator" begin
--- a/test/StaticDicts/StaticDicts_test.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/StaticDicts/StaticDicts_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -1,5 +1,5 @@
 using Test
-using Sbplib.StaticDicts
+using Diffinitive.StaticDicts
 
 @testset "StaticDicts" begin
 
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test/ext/sparse_array_kit_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -0,0 +1,35 @@
+using Test
+
+using Diffinitive
+using Diffinitive.Grids
+using Diffinitive.SbpOperators
+
+using SparseArrayKit
+using Tokens
+using Tullio
+
+
+@testset "SparseArray" begin
+    g = equidistant_grid((0,0),(1,2), 20,30)
+    stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=4)
+
+
+    @testset let Δ = laplace(g, stencil_set), M = SparseArray(Δ)
+        @test ndims(M) == 4
+        @test size(M) == (20,30,20,30)
+
+        v = rand(size(g)...)
+        @tullio Mv[i,j] := M[i,j,k,l]*v[k,l]
+
+        @test Mv ≈ Δ*v
+    end
+
+    @testset let dₙ = normal_derivative(g, stencil_set,CartesianBoundary{1,LowerBoundary}()), M = SparseArray(dₙ)
+        @test ndims(M) == 3
+        @test size(M) == (30,20,30)
+
+        v = rand(size(g)...)
+        @tullio Mv[i] := M[i,j,k]*v[j,k]
+        @test Mv ≈ dₙ*v
+    end
+end
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test/ext/sparse_arrays_test.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -0,0 +1,34 @@
+using Test
+
+using Diffinitive
+using Diffinitive.Grids
+using Diffinitive.SbpOperators
+
+using SparseArrays
+using Tokens
+
+
+@testset "SparseArray" begin
+    g = equidistant_grid((0,0),(1,2), 20,30)
+    stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=4)
+
+
+    @testset let Δ = laplace(g, stencil_set), M = sparse(Δ)
+        @test ndims(M) == 2
+        @test size(M) == (20*30,20*30)
+
+        v = rand(size(g)...)
+
+        Mv = M*reshape(v,:)
+        @test Mv ≈ reshape(Δ*v,:)
+    end
+
+    @testset let dₙ = normal_derivative(g, stencil_set,CartesianBoundary{1,LowerBoundary}()), M = sparse(dₙ)
+        @test ndims(M) == 2
+        @test size(M) == (30,20*30)
+
+        v = rand(size(g)...)
+        Mv = M*reshape(v,:)
+        @test Mv ≈ reshape(dₙ*v,:)
+    end
+end
--- a/test/runtests.jl	Tue Sep 10 21:35:55 2024 +0200
+++ b/test/runtests.jl	Tue Sep 10 21:59:10 2024 +0200
@@ -45,7 +45,7 @@
     end
 end
 
-testsetname = isempty(ARGS) ? "Sbplib.jl" : "["*join(ARGS, ", ")*"]"
+testsetname = isempty(ARGS) ? "Diffinitive.jl" : "["*join(ARGS, ", ")*"]"
 
 @testset "$testsetname" begin
     run_testfiles(ARGS)