--- a/src/Grids/grid.jl	Tue Jun 27 09:36:32 2023 +0200
+++ b/src/Grids/grid.jl	Wed Aug 23 15:51:32 2023 +0200
@@ -22,6 +22,8 @@
 Base.ndims(::Grid{T,D}) where {T,D} = D
 Base.eltype(::Type{<:Grid{T}}) where T = T
 
+Base.getindex(g::Grid, I::CartesianIndex) = g[Tuple(I)...]
+
 """
     coordinate_size(g)
 
@@ -74,7 +76,7 @@
 """
     eval_on(g::Grid, f)
 
-Lazy evaluation `f` on the grid. `f` can either be on the form `f(x,y,...)`
+Lazy evaluation of `f` on the grid. `f` can either be on the form `f(x,y,...)`
 with each coordinate as an argument, or on the form `f(x̄)` taking a
 coordinate vector.
 
@@ -89,5 +91,12 @@
     end
 end
 
+"""
+    eval_on(g::Grid, f::Number)
+
+Lazy evaluation of a scalar `f` on the grid.
+"""
+eval_on(g::Grid, f::Number) = return LazyTensors.LazyConstantArray(f, size(g))
+
 _ncomponents(::Type{<:Number}) = 1
 _ncomponents(T::Type{<:SVector}) = length(T)

--- a/src/Grids/tensor_grid.jl	Tue Jun 27 09:36:32 2023 +0200
+++ b/src/Grids/tensor_grid.jl	Wed Aug 23 15:51:32 2023 +0200
@@ -17,7 +17,7 @@
 end
 
 # Indexing interface
-function Base.getindex(g::TensorGrid, I...)
+function Base.getindex(g::TensorGrid, I::Vararg{Int})
     szs = ndims.(g.grids)
 
     Is = LazyTensors.split_tuple(I, szs)
@@ -26,8 +26,6 @@
     return vcat(ps...)
 end
 
-Base.getindex(g::TensorGrid, I::CartesianIndex) = g[Tuple(I)...]
-
 function Base.eachindex(g::TensorGrid)
     szs = LazyTensors.concatenate_tuples(size.(g.grids)...)
     return CartesianIndices(szs)

--- a/test/Grids/grid_test.jl	Tue Jun 27 09:36:32 2023 +0200
+++ b/test/Grids/grid_test.jl	Wed Aug 23 15:51:32 2023 +0200
@@ -38,6 +38,9 @@
     @test eval_on(ZeroDimGrid(1.), x̄->2x̄) isa LazyArray
     @test eval_on(ZeroDimGrid(1.), x̄->2x̄) == fill(2.)
 
+    @test eval_on(ZeroDimGrid(@SVector[1.,2.]), π) isa LazyArray
+    @test eval_on(ZeroDimGrid(@SVector[1.,2.]), π) == fill(π)
+
     @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)
 
@@ -49,6 +52,8 @@
 
     @test eval_on(g, x̄ -> sin(x̄[1])*cos(x̄[2])) == map(x̄->sin(x̄[1])*cos(x̄[2]), g)
 
+    @test eval_on(g, π) == fill(π, (5,3))
+
     # 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)