--- a/test/testLazyTensors.jl	Sat Nov 07 13:31:55 2020 +0100
+++ b/test/testLazyTensors.jl	Thu Nov 26 21:56:33 2020 +0100
@@ -340,61 +340,115 @@
     B = LazyLinearMap(B̃,(1,2),(3,))
     C = LazyLinearMap(C̃,(1,),(2,3))
 
-    @test InflatedTensorMapping(I(3,2), A, I(4)) isa TensorMapping{Float64, 4, 4}
-    @test InflatedTensorMapping(I(3,2), B, I(4)) isa TensorMapping{Float64, 5, 4}
-    @test InflatedTensorMapping(I(3), C, I(2,3)) isa TensorMapping{Float64, 4, 5}
-    @test InflatedTensorMapping(C, I(2,3)) isa TensorMapping{Float64, 3, 4}
-    @test InflatedTensorMapping(I(3), C) isa TensorMapping{Float64, 2, 3}
-    @test InflatedTensorMapping(I(3), I(2,3)) isa TensorMapping{Float64, 3, 3}
-
-    @test range_size(InflatedTensorMapping(I(3,2), A, I(4))) == (3,2,4,4)
-    @test domain_size(InflatedTensorMapping(I(3,2), A, I(4))) == (3,2,2,4)
-
-    @test range_size(InflatedTensorMapping(I(3,2), B, I(4))) == (3,2,4,2,4)
-    @test domain_size(InflatedTensorMapping(I(3,2), B, I(4))) == (3,2,3,4)
-
-    @test range_size(InflatedTensorMapping(I(3), C, I(2,3))) == (3,4,2,3)
-    @test domain_size(InflatedTensorMapping(I(3), C, I(2,3))) == (3,2,3,2,3)
-
-    @inferred range_size(InflatedTensorMapping(I(3,2), A, I(4))) == (3,2,4,4)
-    @inferred domain_size(InflatedTensorMapping(I(3,2), A, I(4))) == (3,2,2,4)
+    @testset "Constructors" begin
+        @test InflatedTensorMapping(I(3,2), A, I(4)) isa TensorMapping{Float64, 4, 4}
+        @test InflatedTensorMapping(I(3,2), B, I(4)) isa TensorMapping{Float64, 5, 4}
+        @test InflatedTensorMapping(I(3), C, I(2,3)) isa TensorMapping{Float64, 4, 5}
+        @test InflatedTensorMapping(C, I(2,3)) isa TensorMapping{Float64, 3, 4}
+        @test InflatedTensorMapping(I(3), C) isa TensorMapping{Float64, 2, 3}
+        @test InflatedTensorMapping(I(3), I(2,3)) isa TensorMapping{Float64, 3, 3}
+    end
 
-    # Test InflatedTensorMapping mapping w. before and after
-    tm = InflatedTensorMapping(I(3,2), A, I(4))
-    v = rand(domain_size(tm)...)
-    @tullio IAIv[a,b,c,d] := Ã[c,i]*v[a,b,i,d]
-    @test tm*v ≈ IAIv rtol=1e-14
-    @inferred LazyTensors.split_index(tm,1,1,1,1)
+    @testset "Range and domain size" begin
+        @test range_size(InflatedTensorMapping(I(3,2), A, I(4))) == (3,2,4,4)
+        @test domain_size(InflatedTensorMapping(I(3,2), A, I(4))) == (3,2,2,4)
 
-    # Test InflatedTensorMapping mapping w. before
-    tm = InflatedTensorMapping(I(3,2), A)
-    v = rand(domain_size(tm)...)
-    @tullio IAIv[a,b,c] := Ã[c,i]*v[a,b,i]
-    @test tm*v ≈ IAIv rtol=1e-14
-    @inferred LazyTensors.split_index(tm,1,1,1)
+        @test range_size(InflatedTensorMapping(I(3,2), B, I(4))) == (3,2,4,2,4)
+        @test domain_size(InflatedTensorMapping(I(3,2), B, I(4))) == (3,2,3,4)
 
-    # Test InflatedTensorMapping mapping w. after
-    tm = InflatedTensorMapping(A,I(4))
-    v = rand(domain_size(tm)...)
-    @tullio IAIv[c,d] := Ã[c,i]*v[i,d]
-    @test tm*v ≈ IAIv rtol=1e-14
-    @inferred LazyTensors.split_index(tm,1,1)
+        @test range_size(InflatedTensorMapping(I(3), C, I(2,3))) == (3,4,2,3)
+        @test domain_size(InflatedTensorMapping(I(3), C, I(2,3))) == (3,2,3,2,3)
 
-    struct ScalingOperator{T,D} <: TensorMapping{T,D,D}
-        λ::T
-        size::NTuple{D,Int}
+        @inferred range_size(InflatedTensorMapping(I(3,2), A, I(4))) == (3,2,4,4)
+        @inferred domain_size(InflatedTensorMapping(I(3,2), A, I(4))) == (3,2,2,4)
     end
 
-    LazyTensors.apply(m::ScalingOperator{T,D}, v, I::Vararg{Index,D}) where {T,D} = m.λ*v[I]
-    LazyTensors.range_size(m::ScalingOperator) = m.size
-    LazyTensors.domain_size(m::ScalingOperator) = m.size
+    @testset "Application" begin
+        # Testing regular application and transposed application with inflation "before", "after" and "before and after".
+        # The inflated tensor mappings are chosen to preserve, reduce and increase the dimension of the result compared to the input.
+        tests = [
+            (
+                InflatedTensorMapping(I(3,2), A, I(4)),
+                (v-> @tullio res[a,b,c,d] := Ã[c,i]*v[a,b,i,d]), # Expected result of apply
+                (v-> @tullio res[a,b,c,d] := Ã[i,c]*v[a,b,i,d]), # Expected result of apply_transpose
+            ),
+            (
+                InflatedTensorMapping(I(3,2), B, I(4)),
+                (v-> @tullio res[a,b,c,d,e] := B̃[c,d,i]*v[a,b,i,e]),
+                (v-> @tullio res[a,b,c,d] := B̃[i,j,c]*v[a,b,i,j,d]),
+            ),
+            (
+                InflatedTensorMapping(I(3,2), C, I(4)),
+                (v-> @tullio res[a,b,c,d] := C̃[c,i,j]*v[a,b,i,j,d]),
+                (v-> @tullio res[a,b,c,d,e] := C̃[i,c,d]*v[a,b,i,e]),
+            ),
+            (
+                InflatedTensorMapping(I(3,2), A),
+                (v-> @tullio res[a,b,c] := Ã[c,i]*v[a,b,i]),
+                (v-> @tullio res[a,b,c] := Ã[i,c]*v[a,b,i]),
+            ),
+            (
+                InflatedTensorMapping(I(3,2), B),
+                (v-> @tullio res[a,b,c,d] := B̃[c,d,i]*v[a,b,i]),
+                (v-> @tullio res[a,b,c] := B̃[i,j,c]*v[a,b,i,j]),
+            ),
+            (
+                InflatedTensorMapping(I(3,2), C),
+                (v-> @tullio res[a,b,c] := C̃[c,i,j]*v[a,b,i,j]),
+                (v-> @tullio res[a,b,c,d] := C̃[i,c,d]*v[a,b,i]),
+            ),
+            (
+                InflatedTensorMapping(A,I(4)),
+                (v-> @tullio res[a,b] := Ã[a,i]*v[i,b]),
+                (v-> @tullio res[a,b] := Ã[i,a]*v[i,b]),
+            ),
+            (
+                InflatedTensorMapping(B,I(4)),
+                (v-> @tullio res[a,b,c] := B̃[a,b,i]*v[i,c]),
+                (v-> @tullio res[a,b] := B̃[i,j,a]*v[i,j,b]),
+            ),
+            (
+                InflatedTensorMapping(C,I(4)),
+                (v-> @tullio res[a,b] := C̃[a,i,j]*v[i,j,b]),
+                (v-> @tullio res[a,b,c] := C̃[i,a,b]*v[i,c]),
+            ),
+        ]
 
-    tm = InflatedTensorMapping(I(2,3),ScalingOperator(2.0, (3,2)),I(3,4))
-    v = rand(domain_size(tm)...)
+        @testset "apply" begin
+            for i ∈ 1:length(tests)
+                tm = tests[i][1]
+                v = rand(domain_size(tm)...)
+                true_value = tests[i][2](v)
+                @test tm*v ≈ true_value rtol=1e-14
+            end
+        end
+
+        @testset "apply_transpose" begin
+            for i ∈ 1:length(tests)
+                tm = tests[i][1]
+                v = rand(range_size(tm)...)
+                true_value = tests[i][3](v)
+                @test tm'*v ≈ true_value rtol=1e-14
+            end
+        end
 
-    @inferred LazyTensors.split_index(tm,1,2,3,2,2,4)
-    @inferred apply(tm,v,Index{Unknown}.((1,2,3,2,2,4))...)
-    @inferred (tm*v)[1,2,3,2,2,4]
+        @testset "Inference of application" begin
+            struct ScalingOperator{T,D} <: TensorMapping{T,D,D}
+                λ::T
+                size::NTuple{D,Int}
+            end
+
+            LazyTensors.apply(m::ScalingOperator{T,D}, v, I::Vararg{Index,D}) where {T,D} = m.λ*v[I]
+            LazyTensors.range_size(m::ScalingOperator) = m.size
+            LazyTensors.domain_size(m::ScalingOperator) = m.size
+
+            tm = InflatedTensorMapping(I(2,3),ScalingOperator(2.0, (3,2)),I(3,4))
+            v = rand(domain_size(tm)...)
+
+            @inferred apply(tm,v,Index{Unknown}.((1,2,3,2,2,4))...)
+            @inferred (tm*v)[1,2,3,2,2,4]
+        end
+    end
 
     @testset "InflatedTensorMapping of InflatedTensorMapping" begin
         A = ScalingOperator(2.0,(2,3))
@@ -405,7 +459,20 @@
 
         @test InflatedTensorMapping(I(2), I(2), I(2)) isa InflatedTensorMapping # The constructor should always return its type.
     end
+end
 
+@testset "split_index" begin
+    @test LazyTensors.split_index(Val(2),Val(1),Val(2),Val(2),1,2,3,4,5,6) == ((1,2,:,5,6),(3,4))
+    @test LazyTensors.split_index(Val(2),Val(3),Val(2),Val(2),1,2,3,4,5,6) == ((1,2,:,:,:,5,6),(3,4))
+    @test LazyTensors.split_index(Val(3),Val(1),Val(1),Val(2),1,2,3,4,5,6) == ((1,2,3,:,5,6),(4,))
+    @test LazyTensors.split_index(Val(3),Val(2),Val(1),Val(2),1,2,3,4,5,6) == ((1,2,3,:,:,5,6),(4,))
+    @test LazyTensors.split_index(Val(1),Val(1),Val(2),Val(3),1,2,3,4,5,6) == ((1,:,4,5,6),(2,3))
+    @test LazyTensors.split_index(Val(1),Val(2),Val(2),Val(3),1,2,3,4,5,6) == ((1,:,:,4,5,6),(2,3))
+
+    @test LazyTensors.split_index(Val(0),Val(1),Val(3),Val(3),1,2,3,4,5,6) == ((:,4,5,6),(1,2,3))
+    @test LazyTensors.split_index(Val(3),Val(1),Val(3),Val(0),1,2,3,4,5,6) == ((1,2,3,:),(4,5,6))
+
+    @inferred LazyTensors.split_index(Val(2),Val(3),Val(2),Val(2),1,2,3,2,2,4)
 end
 
 @testset "slice_tuple" begin
@@ -415,6 +482,32 @@
     @test LazyTensors.slice_tuple((1,2,3,4,5,6),Val(4), Val(6)) == (4,5,6)
 end
 
+@testset "split_tuple" begin
+    @testset "2 parts" begin
+        @test LazyTensors.split_tuple((),Val(0)) == ((),())
+        @test LazyTensors.split_tuple((1,),Val(0)) == ((),(1,))
+        @test LazyTensors.split_tuple((1,),Val(1)) == ((1,),())
+
+        @test LazyTensors.split_tuple((1,2,3,4),Val(0)) == ((),(1,2,3,4))
+        @test LazyTensors.split_tuple((1,2,3,4),Val(1)) == ((1,),(2,3,4))
+        @test LazyTensors.split_tuple((1,2,3,4),Val(2)) == ((1,2),(3,4))
+        @test LazyTensors.split_tuple((1,2,3,4),Val(3)) == ((1,2,3),(4,))
+        @test LazyTensors.split_tuple((1,2,3,4),Val(4)) == ((1,2,3,4),())
+
+        @inferred LazyTensors.split_tuple((1,2,3,4),Val(3))
+    end
+
+    @testset "3 parts" begin
+        @test LazyTensors.split_tuple((),Val(0),Val(0)) == ((),(),())
+        @test LazyTensors.split_tuple((1,2,3),Val(1), Val(1)) == ((1,),(2,),(3,))
+
+        @test LazyTensors.split_tuple((1,2,3,4,5,6),Val(1),Val(2)) == ((1,),(2,3),(4,5,6))
+        @test LazyTensors.split_tuple((1,2,3,4,5,6),Val(3),Val(2)) == ((1,2,3),(4,5),(6,))
+
+        @inferred LazyTensors.split_tuple((1,2,3,4,5,6),Val(3),Val(2))
+    end
+end
+
 @testset "flatten_tuple" begin
     @test LazyTensors.flatten_tuple((1,)) == (1,)
     @test LazyTensors.flatten_tuple((1,2,3,4,5,6)) == (1,2,3,4,5,6)