Add fallback when resize! does not work (#31)

* Add fallback for when resize fails

* Only catch shared data error from ` resize!`

* Test composition of several maps

Author: cako

src/composition.jl

@@ -90,7 +90,15 @@ function A_mul_B!(y::AbstractVector, A::CompositeMap, x::AbstractVector)
             dest = Array{T}(undef, size(A.maps[2], 1))
         end
         for n=2:N-1
-            resize!(dest, size(A.maps[n], 1))
+            try
+                resize!(dest, size(A.maps[n], 1))
+            catch err
+                if err == ErrorException("cannot resize array with shared data")
+                    dest = Array{T}(undef, size(A.maps[n], 1))
+                else
+                    rethrow(err)
+                end
+            end
             A_mul_B!(dest, A.maps[n], source)
             dest, source = source, dest # alternate dest and source
         end
@@ -112,7 +120,15 @@ function At_mul_B!(y::AbstractVector, A::CompositeMap, x::AbstractVector)
             dest = Array{T}(undef, size(A.maps[N-1], 2))
         end
         for n = N-1:-1:2
-            resize!(dest, size(A.maps[n], 2))
+            try
+                resize!(dest, size(A.maps[n], 2))
+            catch err
+                if err == ErrorException("cannot resize array with shared data")
+                    dest = Array{T}(undef, size(A.maps[n], 2))
+                else
+                    rethrow(err)
+                end
+            end
             At_mul_B!(dest, A.maps[n], source)
             dest, source = source, dest # alternate dest and source
         end
@@ -134,7 +150,15 @@ function Ac_mul_B!(y::AbstractVector, A::CompositeMap, x::AbstractVector)
             dest = Array{T}(undef, size(A.maps[N-1], 2))
         end
         for n = N-1:-1:2
-            resize!(dest, size(A.maps[n], 2))
+            try
+                resize!(dest, size(A.maps[n], 2))
+            catch err
+                if err == ErrorException("cannot resize array with shared data")
+                    dest = Array{T}(undef, size(A.maps[n], 2))
+                else
+                    rethrow(err)
+                end
+            end
             Ac_mul_B!(dest, A.maps[n], source)
             dest, source = source, dest # alternate dest and source
         end

test/runtests.jl

@@ -179,3 +179,42 @@ mul!(w, F, v)
 @test_throws MethodError transpose(F) * v
 @test_throws MethodError mul!(w, adjoint(F), v)
 @test_throws MethodError mul!(w, transpose(F), v)
+
+# test composition of several maps with shared data #31
+sizes = ( (5, 2), (3, 3), (3, 2), (2, 2), (9, 2), (7, 1) )
+N = length(sizes) - 1
+Lf = []
+Lt = []
+Lc = []
+
+# build list of operators [LN, ..., L2, L1] for each mode
+for (fi, i) in [ (Symbol("f$i"), i) for i in 1:N]
+    @eval begin
+        function ($fi)(source)
+            dest = ones(prod(sizes[$i + 1]))
+            tmp = reshape(source, sizes[$i])
+            return conj.($i * dest)
+        end
+        insert!(Lf, 1, LinearMap($fi, prod(sizes[$i + 1]), prod(sizes[$i])))
+        insert!(Lt, 1, LinearMap(x -> x, $fi, prod(sizes[$i]), prod(sizes[$i + 1])))
+        insert!(Lc, 1, LinearMap{ComplexF64}(x -> x, $fi, prod(sizes[$i]), prod(sizes[$i + 1])))
+    end
+end
+
+# multiply as composition and as recursion
+v1 = ones(prod(sizes[1]))
+u1 = ones(prod(sizes[1]))
+w1 = im.*ones(ComplexF64, prod(sizes[1]))
+for i = N:-1:1
+    v2 = prod(Lf[i:N]) * ones(prod(sizes[1]))
+    u2 = transpose(LinearMap(prod(Lt[N:-1:i]))) * ones(prod(sizes[1]))
+    w2 = adjoint(LinearMap(prod(Lc[N:-1:i]))) * ones(prod(sizes[1]))
+
+    global v1 = Lf[i] * v1
+    global u1 = transpose(Lt[i]) * u1
+    global w1 = adjoint(Lc[i]) * w1
+
+    @test v1 == v2
+    @test u1 == u2
+    @test w1 == w2
+end