--- a/src/Grids/curvilinear_grid.jl	Fri Aug 25 08:57:27 2023 +0200
+++ b/src/Grids/curvilinear_grid.jl	Wed Jul 26 21:00:44 2023 +0200
@@ -1,1 +1,84 @@
-struct CurvilinearGrid end
+struct CurvilinearGrid{T,D, GT<:Grid{<:Any,D}, CT<:AbstractArray{T,D}, JT<:AbstractArray{<:AbstractArray{<:Any, 2}, D}} <: Grid{T,D}
+    logicalgrid::GT
+    physicalcoordinates::CT
+    Jacobian::JT
+end
+
+
+# Indexing interface
+Base.getindex(g::CurvilinearGrid, I...) = g.physicalcoordinates[I...]
+
+# function Base.getindex(g::TensorGrid, I...)
+#     szs = ndims.(g.grids)
+
+#     Is = LazyTensors.split_tuple(I, szs)
+#     ps = map((g,I)->SVector(g[I...]), g.grids, Is)
+
+#     return vcat(ps...)
+# end
+
+# function Base.eachindex(g::TensorGrid)
+#     szs = LazyTensors.concatenate_tuples(size.(g.grids)...)
+#     return CartesianIndices(szs)
+# end
+
+# # Iteration interface
+# Base.iterate(g::TensorGrid) = iterate(Iterators.product(g.grids...)) |> _iterate_combine_coords
+# Base.iterate(g::TensorGrid, state) = iterate(Iterators.product(g.grids...), state) |> _iterate_combine_coords
+# _iterate_combine_coords(::Nothing) = nothing
+# _iterate_combine_coords((next,state)) = combine_coordinates(next...), state
+
+# Base.IteratorSize(::Type{<:TensorGrid{<:Any, D}}) where D = Base.HasShape{D}()
+# Base.eltype(::Type{<:TensorGrid{T}}) where T = T
+# Base.length(g::TensorGrid) = sum(length, g.grids)
+# Base.size(g::TensorGrid) = LazyTensors.concatenate_tuples(size.(g.grids)...)
+
+
+# refine(g::TensorGrid, r::Int) = mapreduce(g->refine(g,r), TensorGrid, g.grids)
+# coarsen(g::TensorGrid, r::Int) = mapreduce(g->coarsen(g,r), TensorGrid, g.grids)
+
+# """
+#     TensorGridBoundary{N, BID} <: BoundaryIdentifier
+
+# A boundary identifier for a tensor grid. `N` Specifies which grid in the
+# tensor product and `BID` which boundary on that grid.
+# """
+# struct TensorGridBoundary{N, BID} <: BoundaryIdentifier end
+# grid_id(::TensorGridBoundary{N, BID}) where {N, BID} = N
+# boundary_id(::TensorGridBoundary{N, BID}) where {N, BID} = BID()
+
+# """
+#     boundary_identifiers(g::TensorGrid)
+
+# Returns a tuple containing the boundary identifiers of `g`.
+# """
+# function boundary_identifiers(g::TensorGrid)
+#     per_grid = map(eachindex(g.grids)) do i
+#         return map(bid -> TensorGridBoundary{i, typeof(bid)}(), boundary_identifiers(g.grids[i]))
+#     end
+#     return LazyTensors.concatenate_tuples(per_grid...)
+# end
+
+
+# """
+#     boundary_grid(g::TensorGrid, id::TensorGridBoundary)
+
+# The grid for the boundary of `g` specified by `id`.
+# """
+# function boundary_grid(g::TensorGrid, id::TensorGridBoundary)
+#     local_boundary_grid = boundary_grid(g.grids[grid_id(id)], boundary_id(id))
+#     new_grids = Base.setindex(g.grids, local_boundary_grid, grid_id(id))
+#     return TensorGrid(new_grids...)
+# end
+
+
+
+
+
+
+
+
+# Do we add a convenience function `curvilinear_grid`? It could help with
+# creating the logical grid, evaluating functions and possibly calculating the
+# entries in the jacobian.
+