DiagOp not working with @view on linear operator

[SebastianFlassbeck]

Hi,

We had the issue that we wanted to use DiagOp with a view on a list of linear operators, see the example below.
This unfortunately created a type inconsistency between typeof(ops) and typeof([ops...]. The first one is of type SubArray{LinearOperator...} while the second one is of type Vector{LinearOperator...}.

If I'm not mistaken, [ops...] creates a copy of the operator which can be quite costly if the operators are large. So it would be better to avoid this.

I'm unsure though if this change causes issues somewhere else especially since we do not create a copy.

Best,
Sebastian

Example:

N=32
K=4
arrayType = Array
x = arrayType(rand(ComplexF64, K*N))
block = arrayType(rand(ComplexF64, N, N))

F = arrayType(zeros(ComplexF64, K*N, K*N))
for k = 1:K
    start = (k-1)*N + 1
    stop = k*N
    F[start:stop,start:stop] = block
end

blocks = [block for k = 1:K]
op1 = DiagOp(blocks[1:2])
op2 = DiagOp(@view blocks[1:2]) # this line doesn't work with the current version

[nHackel]

Hello,

we definitely should have support for views, thanks for catching that missing feature.

You are correct that [ops...] creates a copy and is a bit inefficient. However, there is a small caveat to that and that is that the call creates copies of the references/pointers to the elements of ops and not the ops themselves. The operator size itself has no impact here, more the number of elements in the view.

We have two options that I can see atm. One is the one you presented here. The other one would be to switch out [ops...] with ops = collect(ops) . This returns again a vector of ops and has to happen before we create the LinearOperator to avoid the typeof(ops) mismatch.

From a mul! standpoint both versions are equal. The only difference I can see is if someone does something special based on the parameters of the DiagOp{T, vecT, vecO}. If there exists a function that expects vecO to be a dense vector, this would not work with the first approach. But I don't think that is a strong argument against your solution and we can just leave it as you proposed.

Could you add a test case to the DiagOps tests here? The tests should essentially do your MWE and see if the normal operator for it also works

We can then merge this and do a patch release. Alternatively, I can merge this now and add the tests myself once I find time and then do a release

[SebastianFlassbeck]

Thanks for your detailed response!

Sounds good, I added the tests for @view.

Best,
Sebastian