--- a/test/SbpOperators/volumeops/laplace/laplace_test.jl	Tue Feb 07 21:55:07 2023 +0100
+++ b/test/SbpOperators/volumeops/laplace/laplace_test.jl	Wed Jul 26 21:35:50 2023 +0200
@@ -4,40 +4,40 @@
 using Sbplib.Grids
 using Sbplib.LazyTensors
 
-# Default stencils (4th order)
-operator_path = sbp_operators_path()*"standard_diagonal.toml"
-stencil_set = read_stencil_set(operator_path; order=4)
-inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"])
-closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"])
-g_1D = EquidistantGrid(101, 0.0, 1.)
-g_3D = EquidistantGrid((51,101,52), (0.0, -1.0, 0.0), (1., 1., 1.))
+@testset "Laplace" begin
+    # Default stencils (4th order)
+    operator_path = sbp_operators_path()*"standard_diagonal.toml"
+    stencil_set = read_stencil_set(operator_path; order=4)
+    g_1D = equidistant_grid(101, 0.0, 1.)
+    g_3D = equidistant_grid((51,101,52), (0.0, -1.0, 0.0), (1., 1., 1.))
 
-@testset "Laplace" begin
     @testset "Constructors" begin
         @testset "1D" begin
-            Δ = laplace(g_1D, inner_stencil, closure_stencils)
-            @test Laplace(g_1D, stencil_set) == Laplace(Δ, stencil_set)
-            @test Laplace(g_1D, stencil_set) isa LazyTensor{T,1,1}  where T
+            @test Laplace(g_1D, stencil_set) == Laplace(laplace(g_1D, stencil_set), stencil_set)
+            @test Laplace(g_1D, stencil_set) isa LazyTensor{Float64,1,1}
         end
         @testset "3D" begin
-            Δ = laplace(g_3D, inner_stencil, closure_stencils)
-            @test Laplace(g_3D, stencil_set) == Laplace(Δ,stencil_set)
-            @test Laplace(g_3D, stencil_set) isa LazyTensor{T,3,3} where T
+            @test Laplace(g_3D, stencil_set) == Laplace(laplace(g_3D, stencil_set),stencil_set)
+            @test Laplace(g_3D, stencil_set) isa LazyTensor{Float64,3,3}
         end
     end
 
     # Exact differentiation is measured point-wise. In other cases
     # the error is measured in the l2-norm.
     @testset "Accuracy" begin
-        l2(v) = sqrt(prod(spacing(g_3D))*sum(v.^2));
+        l2(v) = sqrt(prod(spacing.(g_3D.grids))*sum(v.^2));
         polynomials = ()
         maxOrder = 4;
         for i = 0:maxOrder-1
             f_i(x,y,z) = 1/factorial(i)*(y^i + x^i + z^i)
-            polynomials = (polynomials...,evalOn(g_3D,f_i))
+            polynomials = (polynomials...,eval_on(g_3D,f_i))
         end
-        v = evalOn(g_3D, (x,y,z) -> sin(x) + cos(y) + exp(z))
-        Δv = evalOn(g_3D,(x,y,z) -> -sin(x) - cos(y) + exp(z))
+        # v = eval_on(g_3D, (x,y,z) -> sin(x) + cos(y) + exp(z))
+        # Δv = eval_on(g_3D,(x,y,z) -> -sin(x) - cos(y) + exp(z))
+
+        v =  eval_on(g_3D, x̄ -> sin(x̄[1]) + cos(x̄[2]) + exp(x̄[3]))
+        Δv = eval_on(g_3D, x̄ -> -sin(x̄[1]) - cos(x̄[2]) + exp(x̄[3]))
+        @inferred v[1,2,3]
 
         # 2nd order interior stencil, 1st order boundary stencil,
         # implies that L*v should be exact for binomials up to order 2.
@@ -67,19 +67,24 @@
 end
 
 @testset "laplace" begin
+    operator_path = sbp_operators_path()*"standard_diagonal.toml"
+    stencil_set = read_stencil_set(operator_path; order=4)
+    g_1D = equidistant_grid(101, 0.0, 1.)
+    g_3D = equidistant_grid((51,101,52), (0.0, -1.0, 0.0), (1., 1., 1.))
+
     @testset "1D" begin
-        Δ = laplace(g_1D, inner_stencil, closure_stencils)
-        @test Δ == second_derivative(g_1D, inner_stencil, closure_stencils, 1)
-        @test Δ isa LazyTensor{T,1,1}  where T
+        Δ = laplace(g_1D, stencil_set)
+        @test Δ == second_derivative(g_1D, stencil_set)
+        @test Δ isa LazyTensor{Float64,1,1}
     end
     @testset "3D" begin
-        Δ = laplace(g_3D, inner_stencil, closure_stencils)
-        @test Δ isa LazyTensor{T,3,3} where T
-        Dxx = second_derivative(g_3D, inner_stencil, closure_stencils, 1)
-        Dyy = second_derivative(g_3D, inner_stencil, closure_stencils, 2)
-        Dzz = second_derivative(g_3D, inner_stencil, closure_stencils, 3)
+        Δ = laplace(g_3D, stencil_set)
+        @test Δ isa LazyTensor{Float64,3,3}
+        Dxx = second_derivative(g_3D, stencil_set, 1)
+        Dyy = second_derivative(g_3D, stencil_set, 2)
+        Dzz = second_derivative(g_3D, stencil_set, 3)
         @test Δ == Dxx + Dyy + Dzz
-        @test Δ isa LazyTensor{T,3,3} where T
+        @test Δ isa LazyTensor{Float64,3,3}
     end
 end