SimpleMesh: Implement multidimensional indexing (#671)

* SimpleMesh: Implement multidimensional indexing

* Fix typo

* Apply suggestions

src/mesh/simplemesh.jl

@@ -48,3 +48,12 @@ vertex(m::SimpleMesh, ind::Int) = m.vertices[ind]
 vertices(m::SimpleMesh) = m.vertices
 
 nvertices(m::SimpleMesh) = length(m.vertices)
+
+function Base.getindex(m::SimpleMesh{Dim,T,V,GridTopology{Dim}}, I::CartesianIndices{Dim}) where {Dim,T,V}
+  dims = size(I)
+  odims = size(m)
+  cinds = first(I):CartesianIndex(Tuple(last(I)) .+ 1)
+  inds = vec(LinearIndices(odims .+ 1)[cinds])
+  periodic = isperiodic(topology(m)) .&& dims .== odims
+  SimpleMesh(m.vertices[inds], GridTopology(dims, periodic))
+end

src/mesh/structuredgrid.jl

@@ -37,7 +37,7 @@ XYZ(g::StructuredGrid) = g.XYZ
 
 function Base.getindex(g::StructuredGrid{Dim}, I::CartesianIndices{Dim}) where {Dim}
   dims = size(I)
-  range = first(I):CartesianIndex(Tuple(last(I)) .+ 1)
-  XYZ = ntuple(i -> g.XYZ[i][range], Dim)
+  cinds = first(I):CartesianIndex(Tuple(last(I)) .+ 1)
+  XYZ = ntuple(i -> g.XYZ[i][cinds], Dim)
   StructuredGrid(XYZ, GridTopology(dims))
 end

test/mesh.jl

@@ -432,6 +432,14 @@
     @test vertex(mesh, nvertices(mesh)) == vertex(mesh, size(mesh) .+ 1)
     @test mesh[1, 1] == mesh[1]
     @test mesh[10, 10] == mesh[100]
+    sub = mesh[2:4, 3:7]
+    @test size(sub) == (3, 5)
+    @test minimum(sub) == P2(1, 2)
+    @test maximum(sub) == P2(4, 7)
+    sub = mesh[2, 3:7]
+    @test size(sub) == (1, 5)
+    @test minimum(sub) == P2(1, 2)
+    @test maximum(sub) == P2(2, 7)
 
     # test for https://github.com/JuliaGeometry/Meshes.jl/issues/261
     points = rand(P2, 5)