JuliaSmoothOptimizers / LinearOperators.jl

Linear Operators for Julia
Other
150 stars 32 forks source link

Bug in derived adjoint operation? #302

Closed JakobAsslaender closed 10 months ago

JakobAsslaender commented 10 months ago

It would well be that my math is failing me here, but I wonder if code is correct: https://github.com/JakobAsslaender/LinearOperators.jl/blob/1f0cfb4c6d44d00fdb2d326a255a6b8697b1123e/src/adjtrans.jl#L135C1-L135C1

My understanding is that the code aims to perform the following:

y = A^H b = (A^T b^*)^*

where H denotes the adjoint, T the transpose, and * the conjugate. However, I think the conj(b) part (and also the conj!(y) part before the 5-arg mul!) are missing. It is correctly implemented when inferring the transpose from the adjoint.

Is my math off at 9pm? Thanks for thinking this through with me!

dpo commented 10 months ago

The conj!(y) part is here but the conj(b) part indeed seems to be missing.

@geoffroyleconte Do you agree? Somehow that part of the code is missing a unit test.

JakobAsslaender commented 10 months ago

I think the unit test only uses a matrix as an input and if I understand the corresponding constructor correctly, it defines the transpose and adjoint explicitly. I added a test to the PR. I am sure one could be more exhaustive by, e.g. testing the derived multiplication for hermitian and symmetric operators, but it is at least a start.