--- a/grid_refactor.md	Fri Feb 24 21:54:39 2023 +0100
+++ b/grid_refactor.md	Sat Feb 25 22:22:50 2023 +0100
@@ -100,10 +100,24 @@
 ## Should Grid have function for the target manifold dimension?
 Where would it be used?
     In the constructor for TensorGrid
+    In eval on if we want to allow multiargument functions
     Elsewhere?
 
 An alternative is to analyze T in Grid{T,D} to find the answer. (See combined_coordinate_vector_type in tensor_grid.jl)
 
+## Lazy version of map for our needs?
+Could be used to
+ * evaulate functions on grids
+ * pick out components of grid functions
+ * More?
+
+Maybe this:
+```julia
+struct LazyMappedArray <: LazyArray
+    f::F
+    v::AT
+end
+```
 
 ## Notes from pluto notebook
 - Är det dåligt att använda ndims om antalet index inte matchar?

--- a/src/Grids/grid.jl	Fri Feb 24 21:54:39 2023 +0100
+++ b/src/Grids/grid.jl	Sat Feb 25 22:22:50 2023 +0100
@@ -44,11 +44,17 @@
 Enumerate the dimensions of the grid.
 """
 dims(grid::Grid) = 1:ndims(grid)
-
+# TBD: Is this function needed? Where is it used?
 
 # TBD: New file grid_functions.jl?
+"""
+TODO:
 
-function eval_on(::Grid) end # TODO: Should return a LazyArray and index the grid
+* Mention map(f,g) if you want a concrete array
+"""
+eval_on(g::Grid, f) = eval_on(g, f, Base.IteratorSize(g)) # TBD: Borde f vara först som i alla map, sum, och dylikt
+eval_on(g::Grid, f, ::Base.HasShape) = LazyTensors.LazyFunctionArray((I...)->f(g[I...]), size(g))
+
 
 """
     getcomponent(gfun, I::Vararg{Int})

--- a/test/Grids/grid_test.jl	Fri Feb 24 21:54:39 2023 +0100
+++ b/test/Grids/grid_test.jl	Sat Feb 25 22:22:50 2023 +0100
@@ -1,15 +1,42 @@
+using Test
 using Sbplib.Grids
-using Test
+using Sbplib.LazyTensors
+using StaticArrays
 
 @testset "Grid" begin
-    @test_broken false
+    struct DummyGrid{T,D} <: Grid{T,D} end
+
+    @test eltype(DummyGrid{Int, 2}) == Int
+    @test eltype(DummyGrid{Int, 2}()) == Int
+
+    @test ndims(DummyGrid{Int, 2}()) == 2
+    @test dims(DummyGrid{Int, 2}()) == 1:2
 end
 
 @testset "eval_on" begin
+    @test eval_on(ZeroDimGrid(@SVector[1.,2.]), x̄->x̄[1]+x̄[2]) isa LazyArray
+    @test eval_on(ZeroDimGrid(@SVector[1.,2.]), x̄->x̄[1]+x̄[2]) == fill(3.)
+    @test eval_on(ZeroDimGrid(@SVector[3.,2.]), x̄->x̄[1]+x̄[2]) == fill(5.)
+
+    @test eval_on(EquidistantGrid(range(0,1,length=4)), x->2x) isa LazyArray
+    @test eval_on(EquidistantGrid(range(0,1,length=4)), x->2x) == 2 .* range(0,1,length=4)
+
+
     g = equidistant_grid((5,3), (0.0,0.0), (2.0,1.0))
 
-    @test_broken eval_on(g, (x,y) -> 0.) isa LazyArray
-    @test_broken eval_on(g, (x,y) -> 0.) == fill(0., (5,3))
+
+    # Splat for only one dim, controllef by type specification in function.
+
+    @test eval_on(g, x̄ -> 0.) isa LazyArray
+    @test eval_on(g, x̄ -> 0.) == fill(0., (5,3))
+
+    @test eval_on(g, x̄ -> sin(x̄[1])*cos(x̄[2])) == map(x̄->sin(x̄[1])*cos(x̄[2]), g)
+
+    # Vector valued function
+    @test eval_on(g, x̄ -> @SVector[x̄[2], x̄[1]]) isa LazyArray{SVector{2,Float64}}
+    @test eval_on(g, x̄ -> @SVector[x̄[2], x̄[1]]) == map(x̄ -> @SVector[x̄[2], x̄[1]], g)
+
+
 
     f(x,y) = sin(x)*cos(y)
     @test_broken eval_on(g, f) == map(p->f(p...), points(g))