Added integration weights for Fourier discretization

Author: brownbaerchen

pySDC/helpers/spectral_helper.py

@@ -384,6 +384,7 @@ def get_integration_matrix(self, lbnd=0):
         return S
 
     def get_integration_weights(self):
+        """Weights for integration across entire domain"""
         n = self.xp.arange(self.N, dtype=float)
 
         weights = (-1) ** n + 1
@@ -821,6 +822,12 @@ def get_integration_matrix(self, p=1):
         k[0] = 1j * self.L
         return self.linalg.matrix_power(self.sparse_lib.diags(1 / (1j * k)), p)
 
+    def get_integration_weights(self):
+        """Weights for integration across entire domain"""
+        weights = self.xp.zeros(self.N)
+        weights[0] = self.L / self.N
+        return weights
+
     def get_plan(self, u, forward, *args, **kwargs):
         if self.fft_lib.__name__ == 'mpi4py_fft.fftw':
             if 'axes' in kwargs.keys():

pySDC/tests/test_helpers/test_spectral_helper_1d_fft.py

@@ -112,6 +112,36 @@ def test_integration_matrix(N, plot=False):
     assert np.allclose(expect, Du)
 
 
+@pytest.mark.base
+@pytest.mark.parametrize('x0', [-1, 0])
+@pytest.mark.parametrize('x1', [0.789, 1])
+@pytest.mark.parametrize('N', [32, 45])
+def test_integral_whole_interval(x0, x1, N):
+    import numpy as np
+    from pySDC.helpers.spectral_helper import FFTHelper
+    from qmat.lagrange import LagrangeApproximation
+
+    helper = FFTHelper(N, x0=x0, x1=x1)
+    x = helper.get_1dgrid()
+
+    u = np.zeros_like(x)
+
+    num_coef = N // 2 - 1
+    coeffs = np.random.random((2, N))
+    u += coeffs[0, 0]
+    for i in range(1, num_coef + 1):
+        u += coeffs[0, i] * np.sin(2 * np.pi * i * x / helper.L)
+        u += coeffs[1, i] * np.cos(2 * np.pi * i * x / helper.L)
+
+    u_hat = helper.transform(u)
+
+    weights = helper.get_integration_weights()
+    integral = weights @ u_hat
+    integral_ref = coeffs[0, 0] * helper.L
+
+    assert np.isclose(integral, integral_ref, atol=1e-7), abs(integral_ref - integral)
+
+
 @pytest.mark.base
 @pytest.mark.parametrize('N', [4, 32])
 @pytest.mark.parametrize('v', [0, 4.78])
@@ -141,4 +171,5 @@ def test_tau_method(N, v):
     # test_tau_method(6, 1)
     # test_transform(True)
     # test_transform(False)
-    test_transform_cupy(4)
+    # test_transform_cupy(4)
+    test_integral_whole_interval(0, 2, 90)