diff --git a/pySDC/helpers/spectral_helper.py b/pySDC/helpers/spectral_helper.py index 6cb247606..97170381c 100644 --- a/pySDC/helpers/spectral_helper.py +++ b/pySDC/helpers/spectral_helper.py @@ -117,12 +117,12 @@ def get_1dgrid(self): raise NotImplementedError -class ChebychovHelper(SpectralHelper1D): +class ChebychevHelper(SpectralHelper1D): """ - The Chebychov base consists of special kinds of polynomials, with the main advantage that you can easily transform + The Chebychev base consists of special kinds of polynomials, with the main advantage that you can easily transform between physical and spectral space by discrete cosine transform. - The differentiation in the Chebychov T base is dense, but can be preconditioned to yield a differentiation operator - that moves to Chebychov U basis during differentiation, which is sparse. When using this technique, problems need to + The differentiation in the Chebychev T base is dense, but can be preconditioned to yield a differentiation operator + that moves to Chebychev U basis during differentiation, which is sparse. When using this technique, problems need to be formulated in first order formulation. This implementation is largely based on the Dedalus paper (arXiv:1905.10388). @@ -142,9 +142,9 @@ def __init__(self, *args, transform_type='fft', x0=-1, x1=1, **kwargs): def get_1dgrid(self): ''' - Generates a 1D grid with Chebychov points. These are clustered at the boundary. You need this kind of grid to - use discrete cosine transformation (DCT) to get the Chebychov representation. If you want a different grid, you - need to do an affine transformation before any Chebychov business. + Generates a 1D grid with Chebychev points. These are clustered at the boundary. You need this kind of grid to + use discrete cosine transformation (DCT) to get the Chebychev representation. If you want a different grid, you + need to do an affine transformation before any Chebychev business. Returns: numpy.ndarray: 1D grid @@ -157,8 +157,8 @@ def get_wavenumbers(self): def get_conv(self, name, N=None): ''' Get conversion matrix between different kinds of polynomials. The supported kinds are - - T: Chebychov polynomials of first kind - - U: Chebychov polynomials of second kind + - T: Chebychev polynomials of first kind + - U: Chebychev polynomials of second kind - D: Dirichlet recombination. You get the desired matrix by choosing a name as ``A2B``. I.e. ``T2U`` for the conversion matrix from T to U. @@ -420,11 +420,11 @@ def get_Dirichlet_recombination_matrix(self): return sp.eye(N) - sp.diags(xp.ones(N - 2), offsets=+2) -class Ultraspherical(ChebychovHelper): +class Ultraspherical(ChebychevHelper): """ This implementation follows https://doi.org/10.1137/120865458. - The ultraspherical method works in Chebychov polynomials as well, but also uses various Gegenbauer polynomials. - The idea is that for every derivative of Chebychov T polynomials, there is a basis of Gegenbauer polynomials where the differentiation matrix is a single off-diagonal. + The ultraspherical method works in Chebychev polynomials as well, but also uses various Gegenbauer polynomials. + The idea is that for every derivative of Chebychev T polynomials, there is a basis of Gegenbauer polynomials where the differentiation matrix is a single off-diagonal. There are also conversion operators from one derivative basis to the next that are sparse. This basis is used like this: For every equation that you have, look for the highest derivative and bump all matrices to the correct basis. If your highest derivative is 2 and you have an identity, it needs to get bumped from 0 to 1 and from 1 to 2. If you have a first derivative as well, it needs to be bumped from 1 to 2. @@ -672,7 +672,7 @@ def add_axis(self, base, *args, **kwargs): if base.lower() in ['chebychov', 'chebychev', 'cheby', 'chebychovhelper']: kwargs['transform_type'] = kwargs.get('transform_type', 'fft') - self.axes.append(ChebychovHelper(*args, **kwargs)) + self.axes.append(ChebychevHelper(*args, **kwargs)) elif base.lower() in ['fft', 'fourier', 'ffthelper']: self.axes.append(FFTHelper(*args, **kwargs)) elif base.lower() in ['ultraspherical', 'gegenbauer']: @@ -1146,7 +1146,7 @@ def transform_single_component(self, u, axes=None, padding=None): Transform a single component of the solution """ trfs = { - ChebychovHelper: self._transform_dct, + ChebychevHelper: self._transform_dct, Ultraspherical: self._transform_dct, FFTHelper: self._transform_fft, } @@ -1283,7 +1283,7 @@ def _transform_idct(self, u, axes, padding=None, **kwargs): def itransform_single_component(self, u, axes=None, padding=None): trfs = { FFTHelper: self._transform_ifft, - ChebychovHelper: self._transform_idct, + ChebychevHelper: self._transform_idct, Ultraspherical: self._transform_idct, } @@ -1504,10 +1504,10 @@ def get_Id(self): def get_Dirichlet_recombination_matrix(self, axis=-1, **kwargs): """ - Get Dirichlet recombination matrix along axis. Not that it only makes sense in directions discretized with variations of Chebychov bases. + Get Dirichlet recombination matrix along axis. Not that it only makes sense in directions discretized with variations of Chebychev bases. Args: - axis (int): Axis you discretized with Chebychov + axis (int): Axis you discretized with Chebychev Returns: sparse matrix diff --git a/pySDC/implementations/problem_classes/Burgers.py b/pySDC/implementations/problem_classes/Burgers.py index 1788d7937..2eb569cee 100644 --- a/pySDC/implementations/problem_classes/Burgers.py +++ b/pySDC/implementations/problem_classes/Burgers.py @@ -6,9 +6,17 @@ class Burgers1D(GenericSpectralLinear): """ - See https://en.wikipedia.org/wiki/Burgers'_equation - Discretization is done with a Chebychov method, which requires a first order derivative formulation. + See https://en.wikipedia.org/wiki/Burgers'_equation for the equation that is solved. + Discretization is done with a Chebychev method, which requires a first order derivative formulation. Feel free to do a more efficient implementation using an ultraspherical method to avoid the first order business. + + Parameters: + N (int): Spatial resolution + epsilon (float): viscosity + BCl (float): Value at left boundary + BCr (float): Value at right boundary + f (int): Frequency of the initial conditions + mode (str): 'T2U' or 'T2T'. Use 'T2U' to get sparse differentiation matrices """ dtype_u = mesh @@ -131,7 +139,18 @@ def plot(self, u, t=None, fig=None, comp='u'): # pragma: no cover class Burgers2D(GenericSpectralLinear): """ - See documentation of `Burgers1D`. + See https://en.wikipedia.org/wiki/Burgers'_equation for the equation that is solved. + This implementation is discretized with FFTs in x and Chebychev in z. + + Parameters: + nx (int): Spatial resolution in x direction + nz (int): Spatial resolution in z direction + epsilon (float): viscosity + BCl (float): Value at left boundary + BCr (float): Value at right boundary + fux (int): Frequency of the initial conditions in x-direction + fuz (int): Frequency of the initial conditions in z-direction + mode (str): 'T2U' or 'T2T'. Use 'T2U' to get sparse differentiation matrices """ dtype_u = mesh diff --git a/pySDC/implementations/problem_classes/HeatEquation_Chebychov.py b/pySDC/implementations/problem_classes/HeatEquation_Chebychev.py similarity index 97% rename from pySDC/implementations/problem_classes/HeatEquation_Chebychov.py rename to pySDC/implementations/problem_classes/HeatEquation_Chebychev.py index 8442eeb28..139ef3ad2 100644 --- a/pySDC/implementations/problem_classes/HeatEquation_Chebychov.py +++ b/pySDC/implementations/problem_classes/HeatEquation_Chebychev.py @@ -6,14 +6,14 @@ from pySDC.implementations.problem_classes.generic_spectral import GenericSpectralLinear -class Heat1DChebychov(GenericSpectralLinear): +class Heat1DChebychev(GenericSpectralLinear): dtype_u = mesh dtype_f = mesh def __init__(self, nvars=128, a=0, b=0, f=1, nu=1.0, mode='T2U', **kwargs): self._makeAttributeAndRegister(*locals().keys(), localVars=locals(), readOnly=True) - bases = [{'base': 'chebychov', 'N': nvars}] + bases = [{'base': 'chebychev', 'N': nvars}] components = ['u', 'ux'] super().__init__(bases, components, **kwargs) @@ -153,11 +153,11 @@ def u_exact(self, t, noise=0): return u -class Heat2DChebychov(GenericSpectralLinear): +class Heat2DChebychev(GenericSpectralLinear): dtype_u = mesh dtype_f = mesh - def __init__(self, nx=128, ny=128, base_x='fft', base_y='chebychov', a=0, b=0, c=0, fx=1, fy=1, nu=1.0, **kwargs): + def __init__(self, nx=128, ny=128, base_x='fft', base_y='chebychev', a=0, b=0, c=0, fx=1, fy=1, nu=1.0, **kwargs): assert nx % 2 == 1 or base_x == 'fft' assert ny % 2 == 1 or base_y == 'fft' self._makeAttributeAndRegister(*locals().keys(), localVars=locals(), readOnly=True) @@ -184,20 +184,20 @@ def __init__(self, nx=128, ny=128, base_x='fft', base_y='chebychov', a=0, b=0, c self.setup_M(M_lhs) for base in [base_x, base_y]: - assert base in ['chebychov', 'fft'] + assert base in ['chebychev', 'fft'] alpha = (self.b - self.a) / 2.0 beta = (self.c - self.b) / 2.0 gamma = (self.c + self.a) / 2.0 - if base_x == 'chebychov': + if base_x == 'chebychev': y = self.Y[0, :] if self.base_y == 'fft': self.add_BC(component='u', equation='u', axis=0, x=-1, v=beta * y - alpha + gamma, kind='Dirichlet') self.add_BC(component='ux', equation='ux', axis=0, v=2 * alpha, kind='integral') else: assert a == b, f'Need periodic boundary conditions in x for {base_x} method!' - if base_y == 'chebychov': + if base_y == 'chebychev': x = self.X[:, 0] self.add_BC(component='u', equation='u', axis=1, x=-1, v=alpha * x - beta + gamma, kind='Dirichlet') self.add_BC(component='uy', equation='uy', axis=1, v=2 * beta, kind='integral') diff --git a/pySDC/implementations/problem_classes/RayleighBenard.py b/pySDC/implementations/problem_classes/RayleighBenard.py index a27a6cefe..4a1399652 100644 --- a/pySDC/implementations/problem_classes/RayleighBenard.py +++ b/pySDC/implementations/problem_classes/RayleighBenard.py @@ -9,6 +9,19 @@ class RayleighBenard(GenericSpectralLinear): + """ + Rayleigh-Benard Convection is a variation of incompressible Navier-Stokes. See, for instance https://doi.org/10.1007/s00791-020-00332-3. + + Parameters: + Prandl (float): Prandl number + Rayleigh (float): Rayleigh number + nx (int): Horizontal resolution + nz (int): Vertical resolution + BCs (dict): Can specify boundary conditions here + dealiasing (float): Dealiasing factor for evaluating the non-linear part + comm (mpi4py.Intracomm): Space communicator + """ + dtype_u = mesh dtype_f = imex_mesh @@ -21,7 +34,6 @@ def __init__( BCs=None, dealiasing=3 / 2, comm=None, - debug=False, **kwargs, ): BCs = {} if BCs is None else BCs @@ -50,7 +62,6 @@ def __init__( 'BCs', 'dealiasing', 'comm', - 'debug', localVars=locals(), readOnly=True, ) diff --git a/pySDC/implementations/problem_classes/generic_spectral.py b/pySDC/implementations/problem_classes/generic_spectral.py index 4c4bcdc99..2f922d533 100644 --- a/pySDC/implementations/problem_classes/generic_spectral.py +++ b/pySDC/implementations/problem_classes/generic_spectral.py @@ -150,7 +150,7 @@ def setup_preconditioner(self, Dirichlet_recombination=True, left_preconditioner self.Pl = self.spectral.sparse_lib.csc_matrix(R) - if Dirichlet_recombination and type(self.axes[-1]).__name__ in ['ChebychovHelper, Ultraspherical']: + if Dirichlet_recombination and type(self.axes[-1]).__name__ in ['ChebychevHelper, Ultraspherical']: _Pr = self.spectral.get_Dirichlet_recombination_matrix(axis=-1) else: _Pr = Id diff --git a/pySDC/tests/test_helpers/test_spectral_helper_1d_chebychov.py b/pySDC/tests/test_helpers/test_spectral_helper_1d_chebychev.py similarity index 89% rename from pySDC/tests/test_helpers/test_spectral_helper_1d_chebychov.py rename to pySDC/tests/test_helpers/test_spectral_helper_1d_chebychev.py index e94d7da48..37a1ee717 100644 --- a/pySDC/tests/test_helpers/test_spectral_helper_1d_chebychov.py +++ b/pySDC/tests/test_helpers/test_spectral_helper_1d_chebychev.py @@ -5,9 +5,9 @@ @pytest.mark.parametrize('N', [4, 32]) def test_D2T_conversion_matrices(N): import numpy as np - from pySDC.helpers.spectral_helper import ChebychovHelper + from pySDC.helpers.spectral_helper import ChebychevHelper - cheby = ChebychovHelper(N) + cheby = ChebychevHelper(N) x = np.linspace(-1, 1, N) D2T = cheby.get_conv('D2T') @@ -31,9 +31,9 @@ def test_D2T_conversion_matrices(N): def test_T_U_conversion(N): import numpy as np from scipy.special import chebyt, chebyu - from pySDC.helpers.spectral_helper import ChebychovHelper + from pySDC.helpers.spectral_helper import ChebychevHelper - cheby = ChebychovHelper(N) + cheby = ChebychevHelper(N) T2U = cheby.get_conv('T2U') U2T = cheby.get_conv('U2T') @@ -61,11 +61,11 @@ def eval_poly(poly, coeffs, x): @pytest.mark.base @pytest.mark.parametrize('name', ['T2U', 'T2D', 'T2T']) def test_conversion_inverses(name): - from pySDC.helpers.spectral_helper import ChebychovHelper + from pySDC.helpers.spectral_helper import ChebychevHelper import numpy as np N = 8 - cheby = ChebychovHelper(N) + cheby = ChebychevHelper(N) P = cheby.get_conv(name) Pinv = cheby.get_conv(name[::-1]) assert np.allclose((P @ Pinv).toarray(), np.diag(np.ones(N))) @@ -76,9 +76,9 @@ def test_conversion_inverses(name): def test_differentiation_matrix(N): import numpy as np import scipy - from pySDC.helpers.spectral_helper import ChebychovHelper + from pySDC.helpers.spectral_helper import ChebychevHelper - cheby = ChebychovHelper(N) + cheby = ChebychevHelper(N) x = np.cos(np.pi / N * (np.arange(N) + 0.5)) coeffs = np.random.random(N) norm = cheby.get_norm() @@ -95,9 +95,9 @@ def test_differentiation_matrix(N): @pytest.mark.parametrize('N', [4, 32]) def test_integration_matrix(N): import numpy as np - from pySDC.helpers.spectral_helper import ChebychovHelper + from pySDC.helpers.spectral_helper import ChebychevHelper - cheby = ChebychovHelper(N) + cheby = ChebychevHelper(N) coeffs = np.random.random(N) coeffs[-1] = 0 @@ -116,9 +116,9 @@ def test_integration_matrix(N): def test_transform(N, d, transform_type): import scipy import numpy as np - from pySDC.helpers.spectral_helper import ChebychovHelper + from pySDC.helpers.spectral_helper import ChebychevHelper - cheby = ChebychovHelper(N, transform_type=transform_type) + cheby = ChebychevHelper(N, transform_type=transform_type) u = np.random.random((d, N)) norm = cheby.get_norm() x = cheby.get_1dgrid() @@ -137,10 +137,10 @@ def test_transform(N, d, transform_type): @pytest.mark.base @pytest.mark.parametrize('N', [8, 32]) def test_integration_BC(N): - from pySDC.helpers.spectral_helper import ChebychovHelper + from pySDC.helpers.spectral_helper import ChebychevHelper import numpy as np - cheby = ChebychovHelper(N) + cheby = ChebychevHelper(N) coef = np.random.random(N) BC = cheby.get_integ_BC_row() @@ -155,11 +155,11 @@ def test_integration_BC(N): @pytest.mark.base @pytest.mark.parametrize('N', [4, 32]) def test_norm(N): - from pySDC.helpers.spectral_helper import ChebychovHelper + from pySDC.helpers.spectral_helper import ChebychevHelper import numpy as np import scipy - cheby = ChebychovHelper(N) + cheby = ChebychevHelper(N) coeffs = np.random.random(N) x = cheby.get_1dgrid() norm = cheby.get_norm() @@ -185,10 +185,10 @@ def test_tau_method(bc, N, bc_val): The test verifies that the solution satisfies the perturbed equation and the boundary condition. ''' - from pySDC.helpers.spectral_helper import ChebychovHelper + from pySDC.helpers.spectral_helper import ChebychevHelper import numpy as np - cheby = ChebychovHelper(N) + cheby = ChebychevHelper(N) x = cheby.get_1dgrid() coef = np.append(np.zeros(N - 1), [1]) @@ -221,13 +221,13 @@ def test_tau_method(bc, N, bc_val): def test_tau_method2D(bc, nz, nx, bc_val, plotting=False): ''' solve u_z - 0.1u_xx -u_x + tau P = 0, u(bc) = sin(bc_val*x) -> space-time discretization of advection-diffusion - problem. We do FFT in x-direction and Chebychov in z-direction. + problem. We do FFT in x-direction and Chebychev in z-direction. ''' - from pySDC.helpers.spectral_helper import ChebychovHelper, FFTHelper + from pySDC.helpers.spectral_helper import ChebychevHelper, FFTHelper import numpy as np import scipy.sparse as sp - cheby = ChebychovHelper(nz) + cheby = ChebychevHelper(nz) fft = FFTHelper(nx) # generate grid @@ -239,7 +239,7 @@ def test_tau_method2D(bc, nz, nx, bc_val, plotting=False): bcs = np.sin(bc_val * x) rhs = np.zeros_like(X) rhs[:, -1] = bcs - rhs_hat = fft.transform(rhs, axis=-2) # the rhs is already in Chebychov spectral space + rhs_hat = fft.transform(rhs, axis=-2) # the rhs is already in Chebychev spectral space # generate matrices Dx = fft.get_differentiation_matrix(p=2) * 1e-1 + fft.get_differentiation_matrix() @@ -299,11 +299,11 @@ def test_tau_method2D_diffusion(nz, nx, bc_val, plotting=False): ''' Solve a Poisson problem with funny Dirichlet BCs in z-direction and periodic in x-direction. ''' - from pySDC.helpers.spectral_helper import ChebychovHelper, FFTHelper + from pySDC.helpers.spectral_helper import ChebychevHelper, FFTHelper import numpy as np import scipy.sparse as sp - cheby = ChebychovHelper(nz) + cheby = ChebychevHelper(nz) fft = FFTHelper(nx) # generate grid @@ -315,7 +315,7 @@ def test_tau_method2D_diffusion(nz, nx, bc_val, plotting=False): rhs = np.zeros((2, nx, nz)) # components u and u_x rhs[0, :, -1] = np.sin(bc_val * x) + 1 rhs[1, :, -1] = 3 * np.exp(-((x - 3.6) ** 2)) + np.cos(x) - rhs_hat = fft.transform(rhs, axis=-2) # the rhs is already in Chebychov spectral space + rhs_hat = fft.transform(rhs, axis=-2) # the rhs is already in Chebychev spectral space # generate 1D matrices Dx = fft.get_differentiation_matrix() diff --git a/pySDC/tests/test_problems/test_heat_chebychov.py b/pySDC/tests/test_problems/test_heat_chebychev.py similarity index 83% rename from pySDC/tests/test_problems/test_heat_chebychov.py rename to pySDC/tests/test_problems/test_heat_chebychev.py index 4bd90230e..66564c2b1 100644 --- a/pySDC/tests/test_problems/test_heat_chebychov.py +++ b/pySDC/tests/test_problems/test_heat_chebychev.py @@ -7,13 +7,13 @@ @pytest.mark.parametrize('f', [0, 1]) @pytest.mark.parametrize('noise', [0, 1e-3]) @pytest.mark.parametrize('use_ultraspherical', [True, False]) -def test_heat1d_chebychov(a, b, f, noise, use_ultraspherical, nvars=2**4): +def test_heat1d_chebychev(a, b, f, noise, use_ultraspherical, nvars=2**4): import numpy as np if use_ultraspherical: - from pySDC.implementations.problem_classes.HeatEquation_Chebychov import Heat1DUltraspherical as problem_class + from pySDC.implementations.problem_classes.HeatEquation_Chebychev import Heat1DUltraspherical as problem_class else: - from pySDC.implementations.problem_classes.HeatEquation_Chebychov import Heat1DChebychov as problem_class + from pySDC.implementations.problem_classes.HeatEquation_Chebychev import Heat1DChebychev as problem_class P = problem_class( nvars=nvars, a=a, b=b, f=f, nu=1e-2, left_preconditioner=False, right_preconditioning='T2T', debug=True @@ -42,17 +42,17 @@ def test_heat1d_chebychov(a, b, f, noise, use_ultraspherical, nvars=2**4): @pytest.mark.parametrize('c', [0, 3.1415]) @pytest.mark.parametrize('fx', [2, 1]) @pytest.mark.parametrize('fy', [2, 1]) -@pytest.mark.parametrize('base_x', ['fft', 'chebychov', 'ultraspherical']) -@pytest.mark.parametrize('base_y', ['fft', 'chebychov', 'ultraspherical']) -def test_heat2d_chebychov(a, b, c, fx, fy, base_x, base_y, nx=2**5 + 1, ny=2**5 + 1): +@pytest.mark.parametrize('base_x', ['fft', 'chebychev', 'ultraspherical']) +@pytest.mark.parametrize('base_y', ['fft', 'chebychev', 'ultraspherical']) +def test_heat2d_chebychev(a, b, c, fx, fy, base_x, base_y, nx=2**5 + 1, ny=2**5 + 1): import numpy as np if base_x == 'ultraspherical' or base_y == 'ultraspherical': - from pySDC.implementations.problem_classes.HeatEquation_Chebychov import Heat2DUltraspherical as problem_class + from pySDC.implementations.problem_classes.HeatEquation_Chebychev import Heat2DUltraspherical as problem_class else: - from pySDC.implementations.problem_classes.HeatEquation_Chebychov import Heat2DChebychov as problem_class + from pySDC.implementations.problem_classes.HeatEquation_Chebychev import Heat2DChebychev as problem_class - if base_x == 'chebychov' and base_y == 'ultraspherical' or base_y == 'chebychov' and base_x == 'ultraspherical': + if base_x == 'chebychev' and base_y == 'ultraspherical' or base_y == 'chebychev' and base_x == 'ultraspherical': return None if base_y == 'fft' and (b != c): @@ -87,7 +87,7 @@ def test_heat2d_chebychov(a, b, c, fx, fy, base_x, base_y, nx=2**5 + 1, ny=2**5 def test_SDC(): import numpy as np from pySDC.implementations.sweeper_classes.generic_implicit import generic_implicit - from pySDC.implementations.problem_classes.HeatEquation_Chebychov import Heat1DChebychov + from pySDC.implementations.problem_classes.HeatEquation_Chebychev import Heat1DChebychev from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI from pySDC.implementations.problem_classes.generic_spectral import compute_residual_DAE from pySDC.helpers.stats_helper import get_sorted @@ -117,7 +117,7 @@ def test_SDC(): controller_params['hook_class'] = [LogSDCIterations] description = {} - description['problem_class'] = Heat1DChebychov + description['problem_class'] = Heat1DChebychev description['problem_params'] = problem_params description['sweeper_class'] = generic_implicit description['sweeper_params'] = sweeper_params @@ -139,5 +139,5 @@ def test_SDC(): if __name__ == '__main__': # test_SDC() - # test_heat1d_chebychov(1, 0, 1, 1e-3, True, 2**6) - test_heat2d_chebychov(0, 0, 0, 2, 2, 'ultraspherical', 'fft', 2**6, 2**6) + # test_heat1d_chebychev(1, 0, 1, 1e-3, True, 2**6) + test_heat2d_chebychev(0, 0, 0, 2, 2, 'ultraspherical', 'fft', 2**6, 2**6)