Add adjoint testing to test utilities

Author: gaurav-arya

ext/AbstractFFTsTestUtilsExt.jl

@@ -51,6 +51,26 @@ const TEST_CASES = (
                              dims=3)),
         )
 
+# Perform generic adjoint plan tests 
+function _adjoint_test(P, x; real_plan=false)
+    y = rand(eltype(P * x), size(P * x))
+    # test basic properties
+    @test_broken eltype(P') === typeof(y) # (AbstactFFTs.jl#110)
+    @test fftdims(P') == fftdims(P)
+    @test (P')' === P # test adjoint of adjoint
+    @test size(P') == AbstractFFTs.output_size(P) # test size of adjoint 
+    # test correctness of adjoint and its inverse via the dot test
+    if !real_plan
+        @test dot(y, P * x) ≈ dot(P' * y, x)
+        @test dot(y, P \ x) ≈ dot(P' \ y, x) 
+    else
+        _component_dot(x, y) = dot(real.(x), real.(y)) + dot(imag.(x), imag.(y))
+        @test _component_dot(y, P * copy(x)) ≈ _component_dot(P' * copy(y), x)
+        @test _component_dot(x, P \ copy(y)) ≈ _component_dot(P' \ copy(x), y) 
+    end
+    @test_throws MethodError mul!(x, P', y)
+end
+
 """
     TestUtils.test_complex_fft(ArrayType=Array; test_real=true, test_inplace=true) 
 
@@ -63,8 +83,9 @@ The backend implementation is assumed to be loaded prior to calling this functio
   which the correctness tests are run. Arrays are constructed via
   `convert(ArrayType, ...)`.
 - `test_inplace=true`: whether to test in-place plans. 
+- `test_adjoint=true`: whether to test adjoints of plans. 
 """
-function TestUtils.test_complex_fft(ArrayType=Array; test_inplace=true) 
+function TestUtils.test_complex_fft(ArrayType=Array; test_inplace=true, test_adjoint=true) 
     @testset "correctness of fft, bfft, ifft" begin
         for test_case in TEST_CASES
             _x, dims, _x_fft = test_case.x, test_case.dims, test_case.x_fft
@@ -89,6 +110,9 @@ function TestUtils.test_complex_fft(ArrayType=Array; test_inplace=true)
                 _x_out = similar(x_fft)
                 @test mul!(_x_out, P, x_complexf) ≈ x_fft
                 @test _x_out ≈ x_fft
+                if test_adjoint
+                    _adjoint_test(P, x_complexf)
+                end
             end
             if test_inplace
                 # test IIP plans
@@ -120,6 +144,9 @@ function TestUtils.test_complex_fft(ArrayType=Array; test_inplace=true)
                 _x_complexf = similar(x_complexf)
                 @test mul!(_x_complexf, P, x_fft) ≈ x_scaled
                 @test _x_complexf ≈ x_scaled
+                if test_adjoint
+                    _adjoint_test(P, x_complexf)
+                end
             end
             # test IIP plans
             for P in (plan_bfft!(similar(x_fft), dims),)
@@ -148,6 +175,9 @@ function TestUtils.test_complex_fft(ArrayType=Array; test_inplace=true)
                 _x_complexf = similar(x_complexf)
                 @test mul!(_x_complexf, P, x_fft) ≈ x
                 @test _x_complexf ≈ x
+                if test_adjoint
+                    _adjoint_test(P, x_complexf)
+                end
             end
             # test IIP plans
             if test_inplace
@@ -177,10 +207,11 @@ The backend implementation is assumed to be loaded prior to calling this functio
   which the correctness tests are run. Arrays are constructed via
   `convert(ArrayType, ...)`.
 - `test_inplace=true`: whether to test in-place plans. 
+- `test_adjoint=true`: whether to test adjoints of plans. 
 """
-function TestUtils.test_real_fft(ArrayType=Array; test_inplace=true)
+function TestUtils.test_real_fft(ArrayType=Array; test_inplace=true, test_adjoint=true)
     @testset "correctness of rfft, brfft, irfft" begin
-        for test_case in TEST_CASES[5:5]
+        for test_case in TEST_CASES
             _x, dims, _x_fft = test_case.x, test_case.dims, test_case.x_fft
             x = convert(ArrayType, _x) # dummy array that will be passed to plans
             x_real = float.(x) # for testing mutating real FFTs
@@ -202,6 +233,9 @@ function TestUtils.test_real_fft(ArrayType=Array; test_inplace=true)
                 _x_rfft = similar(x_rfft)
                 @test mul!(_x_rfft, P, copy(x_real)) ≈ x_rfft
                 @test _x_rfft ≈ x_rfft
+                if test_adjoint
+                    _adjoint_test(P, x_real; real_plan=true)
+                end
             end
 
             # BRFFT

src/definitions.jl

@@ -619,6 +619,7 @@ Base.adjoint(p::ScaledPlan) = ScaledPlan(p.p', p.scale)
 
 size(p::AdjointPlan) = output_size(p.p)
 output_size(p::AdjointPlan) = size(p.p)
+fftdims(p::AdjointPlan) = fftdims(p.p)
 
 Base.:*(p::AdjointPlan, x::AbstractArray) = _mul(p, x, ProjectionStyle(p.p))
 

test/runtests.jl

@@ -153,55 +153,6 @@ end
     end
 end
 
-@testset "adjoint" begin
-    @testset "complex fft adjoint" begin
-        for x_shape in ((3,), (3, 4), (3, 4, 5))
-            N = length(x_shape)
-            real_x = randn(x_shape)
-            complex_x = randn(ComplexF64, x_shape)
-            y = randn(ComplexF64, x_shape)
-            for x in (real_x, complex_x)
-                for dims in unique((1, 1:N, N))
-                    P = plan_fft(x, dims)
-                    @test (P')' === P # test adjoint of adjoint
-                    @test size(P') == AbstractFFTs.output_size(P) # test size of adjoint 
-                    @test dot(y, P * x) ≈ dot(P' * y, x) # test validity of adjoint
-                    @test dot(y, P \ x) ≈ dot(P' \ y, x) # test inv of adjoint
-                    @test dot(y, P \ x) ≈ dot(AbstractFFTs.plan_inv(P') * y, x) # test plan_inv of adjoint 
-                    Pinv = plan_ifft(y)
-                    @test (Pinv')' * y == Pinv * y 
-                    @test size(Pinv') == AbstractFFTs.output_size(Pinv) 
-                    @test dot(x, Pinv * y) ≈ dot(Pinv' * x, y)
-                    @test dot(x, Pinv \ y) ≈ dot(Pinv' \ x, y)
-                    @test dot(x, Pinv \ y) ≈ dot(AbstractFFTs.plan_inv(Pinv') * x, y)
-                    @test_throws MethodError mul!(x, P', y)
-                end
-            end
-        end
-    end
-    @testset "real fft adjoint" begin
-        for x in (randn(3), randn(4), randn(3, 4), randn(3, 4, 5)) # test odd and even lengths
-            N = ndims(x)
-            for dims in unique((1, 1:N, N))
-                P = plan_rfft(x, dims)        
-                y = randn(ComplexF64, size(P * x))
-                @test (P')' * x == P * x
-                @test size(P') == AbstractFFTs.output_size(P)
-                @test dot(real.(y), real.(P * x)) + dot(imag.(y), imag.(P * x)) ≈ dot(P' * y, x)
-                @test dot(real.(y), real.(P' \ x)) + dot(imag.(y), imag.(P' \ x)) ≈ dot(P \ y, x)
-                @test dot(real.(y), real.(AbstractFFTs.plan_inv(P') * x)) + 
-                      dot(imag.(y), imag.(AbstractFFTs.plan_inv(P') * x)) ≈ dot(P \ y, x)
-                Pinv = plan_irfft(y, size(x)[first(dims)], dims)
-                @test (Pinv')' * y == Pinv * y
-                @test size(Pinv') == AbstractFFTs.output_size(Pinv) 
-                @test dot(x, Pinv * y) ≈ dot(real.(y), real.(Pinv' * x)) + dot(imag.(y), imag.(Pinv' * x))
-                @test dot(x, Pinv' \ y) ≈ dot(real.(y), real.(Pinv \ x)) + dot(imag.(y), imag.(Pinv \ x))
-                @test dot(x, AbstractFFTs.plan_inv(Pinv') * y) ≈ dot(real.(y), real.(Pinv \ x)) + dot(imag.(y), imag.(Pinv \ x))
-            end
-        end
-    end
-end
-
 # Test that dims defaults to 1:ndims for fft-like functions
 @testset "Default dims" begin
     for x in (randn(3), randn(3, 4), randn(3, 4, 5))