diff --git a/src/operators/multilinempo.jl b/src/operators/multilinempo.jl
index 92e3e791b..d9c09aa06 100644
--- a/src/operators/multilinempo.jl
+++ b/src/operators/multilinempo.jl
@@ -43,3 +43,9 @@ function Base.:*(mpo1::MultilineMPO, mpo2::MultilineMPO)
     size(mpo1) == size(mpo2) || throw(ArgumentError("dimension mismatch"))
     return Multiline(map(*, zip(mpo1, mpo2)))
 end
+
+for f_space in (:physicalspace, :left_virtualspace, :right_virtualspace)
+    @eval $f_space(t::MultilineMPO, i::Int, j::Int) = $f_space(t[i], j)
+    @eval $f_space(t::MultilineMPO, I::CartesianIndex{2}) = $f_space(t, Tuple(I)...)
+    @eval $f_space(t::MultilineMPO) = map(Base.Fix1($f_space, t), eachindex(t))
+end
diff --git a/test/operators/mpo.jl b/test/operators/mpo.jl
index f93d4a9c2..afbd99e09 100644
--- a/test/operators/mpo.jl
+++ b/test/operators/mpo.jl
@@ -86,10 +86,10 @@ end
 
 @testset "InfiniteMPO" begin
     P = ℂ^2
+    V = ℂ^2
     T = Float64
 
-    H1 = randn(T, P ← P)
-    H1 += H1'
+    H1 = randn(T, V ⊗ P ← P ⊗ V)
     H = InfiniteMPO([H1])
 
     @test !isfinite(H)
@@ -98,6 +98,16 @@ end
     @test GeometryStyle(H) == InfiniteChainStyle()
     @test OperatorStyle(typeof(H)) == MPOStyle()
     @test OperatorStyle(H) == MPOStyle()
+
+    @test physicalspace(H, 1) == P
+    @test left_virtualspace(H, 1) == V
+    @test left_virtualspace(H, 4) == V
+    @test right_virtualspace(H, 1) == V
+
+    multiH = MultilineMPO([H, H])
+    @test physicalspace(multiH, 1, 1) == P
+    @test left_virtualspace(multiH, 1, 1) == left_virtualspace(multiH, 2, 1) == V
+    @test right_virtualspace(multiH, CartesianIndex(1, 2)) == V
 end
 
 @testset "Adapt" for V in (ℂ^2, U1Space(-1 => 1, 0 => 1, 1 => 1))