complete chebyshev_filter

Author: yuehhua

src/polynomials.jl

@@ -165,25 +165,23 @@ function chebyshev_filter(k)
     Φ1 = zeros(k, k)
     ϕ, ψ1, ψ2 = chebyshev_ϕ_ψ(k)
 
-    # ----------------------------------------------------------
-
-    # x = Symbol('x')
-    # kUse = 2*k
-    # roots = Poly(chebyshevt(kUse, 2*x-1)).all_roots()
-    # x_m = np.array([rt.evalf(20) for rt in roots]).astype(np.float64)
-    # # x_m[x_m==0.5] = 0.5 + 1e-8 # add small noise to avoid the case of 0.5 belonging to both phi(2x) and phi(2x-1)
-    # # not needed for our purpose here, we use even k always to avoid
-    # wm = np.pi / kUse / 2
-
-    # for ki in range(k):
-    #     for kpi in range(k):
-    #         H0[ki, kpi] = 1/np.sqrt(2) * (wm * phi[ki](x_m/2) * phi[kpi](x_m)).sum()
-    #         G0[ki, kpi] = 1/np.sqrt(2) * (wm * psi(psi1, psi2, ki, x_m/2) * phi[kpi](x_m)).sum()
-    #         H1[ki, kpi] = 1/np.sqrt(2) * (wm * phi[ki]((x_m+1)/2) * phi[kpi](x_m)).sum()
-    #         G1[ki, kpi] = 1/np.sqrt(2) * (wm * psi(psi1, psi2, ki, (x_m+1)/2) * phi[kpi](x_m)).sum()
-
-    #         PHI0[ki, kpi] = (wm * phi[ki](2*x_m) * phi[kpi](2*x_m)).sum() * 2
-    #         PHI1[ki, kpi] = (wm * phi[ki](2*x_m-1) * phi[kpi](2*x_m-1)).sum() * 2
+    k_use = 2k
+    c = convert(Polynomial, gen_poly(Chebyshev, k_use))
+    x_m = roots(c(Polynomial([-1, 2])))  # 2x-1
+    # x_m[x_m==0.5] = 0.5 + 1e-8 # add small noise to avoid the case of 0.5 belonging to both phi(2x) and phi(2x-1)
+    # not needed for our purpose here, we use even k always to avoid
+    wm = π / k_use / 2
+
+    for ki in 0:(k-1)
+        for kpi in 0:(k-1)
+            H0[ki+1, kpi+1] = 1/sqrt(2) * sum(wm .* ϕ[ki+1].(x_m/2) .* ϕ[kpi+1].(x_m))
+            H1[ki+1, kpi+1] = 1/sqrt(2) * sum(wm .* ϕ[ki+1].((x_m.+1)/2) .* ϕ[kpi+1].(x_m))
+            G0[ki+1, kpi+1] = 1/sqrt(2) * sum(wm .* ψ(ψ1, ψ2, ki, x_m/2) .* ϕ[kpi+1].(x_m))
+            G1[ki+1, kpi+1] = 1/sqrt(2) * sum(wm .* ψ(ψ1, ψ2, ki, (x_m.+1)/2) .* ϕ[kpi+1].(x_m))
+            Φ0[ki+1, kpi+1] = 2*sum(wm .* ϕ[ki+1].(2x_m) .* ϕ[kpi+1].(2x_m))
+            Φ1[ki+1, kpi+1] = 2*sum(wm .* ϕ[ki+1].(2 .* x_m .- 1) .* ϕ[kpi+1].(2 .* x_m .- 1))
+        end
+    end
 
     zero_out!(H0)
     zero_out!(H1)

test/polynomials.jl

@@ -85,25 +85,25 @@
     end
 
     @testset "chebyshev_filter" begin
-        # H0, H1, G0, G1, Φ1, Φ2 = NeuralOperators.chebyshev_filter(3)
+        H0, H1, G0, G1, Φ0, Φ1 = NeuralOperators.chebyshev_filter(3)
 
-        # @test H0 ≈ [0.70710678  0.          0.       ;
-        #            -0.5         0.35355339  0.       ;
-        #            -0.25       -0.70710678  0.1767767]
-        # @test H1 ≈ [0.70710678  0.          0.       ;
-        #             0.5         0.35355339  0.       ;
-        #            -0.25        0.70710678  0.1767767]
-        # @test G0 ≈ [0.60944614 0.77940383 0.        ;
-        #             0.66325172 1.02726613 1.14270252;
-        #             0.61723435 0.90708619 1.1562954 ]
-        # @test G1 ≈ [-0.60944614  0.77940383  0.       ;
-        #             0.66325172 -1.02726613  1.14270252;
-        #            -0.61723435  0.90708619 -1.1562954 ]
-        # @test Φ1 ≈ [1.         -0.40715364 -0.21440101;
-        #            -0.40715364  0.84839559 -0.44820615;
-        #            -0.21440101 -0.44820615  0.84002127]
-        # @test Φ2 ≈ [1.          0.40715364 -0.21440101;
-        #             0.40715364  0.84839559  0.44820615;
-        #            -0.21440101  0.44820615  0.84002127]
+        @test H0 ≈ [0.70710678  0.          0.       ;
+                   -0.5         0.35355339  0.       ;
+                   -0.25       -0.70710678  0.1767767]
+        @test H1 ≈ [0.70710678  0.          0.       ;
+                    0.5         0.35355339  0.       ;
+                   -0.25        0.70710678  0.1767767]
+        @test G0 ≈ [0.60944614 0.77940383 0.        ;
+                    0.66325172 1.02726613 1.14270252;
+                    0.61723435 0.90708619 1.1562954 ]
+        @test G1 ≈ [-0.60944614  0.77940383  0.       ;
+                    0.66325172 -1.02726613  1.14270252;
+                   -0.61723435  0.90708619 -1.1562954 ]
+        @test Φ0 ≈ [1.         -0.40715364 -0.21440101;
+                   -0.40715364  0.84839559 -0.44820615;
+                   -0.21440101 -0.44820615  0.84002127]
+        @test Φ1 ≈ [1.          0.40715364 -0.21440101;
+                    0.40715364  0.84839559  0.44820615;
+                   -0.21440101  0.44820615  0.84002127]
     end
 end