Coverage for pySDC/implementations/problem_classes/boussinesq_helpers/buildBoussinesq2DMatrix.py: 100%

1import numpy as np 

2import scipy.sparse as sp 

3 

4from pySDC.implementations.problem_classes.boussinesq_helpers.build2DFDMatrix import ( 

5 get2DMatrix, 

6 get2DUpwindMatrix, 

7) 

8 

9 

10def getBoussinesq2DUpwindMatrix(N, dx, u_adv, order): 

11 Dx = get2DUpwindMatrix(N, dx, order) 

12 

13 # Note: In the equations it is u_t + u_adv* D_x u = ... so in order to comply with the form u_t = M u, 

14 # add a minus sign in front of u_adv 

15 

16 Zero = np.zeros((N[0] * N[1], N[0] * N[1])) 

17 M1 = sp.hstack((-u_adv * Dx, Zero, Zero, Zero), format="csr") 

18 M2 = sp.hstack((Zero, -u_adv * Dx, Zero, Zero), format="csr") 

19 M3 = sp.hstack((Zero, Zero, -u_adv * Dx, Zero), format="csr") 

20 M4 = sp.hstack((Zero, Zero, Zero, -u_adv * Dx), format="csr") 

21 M = sp.vstack((M1, M2, M3, M4), format="csr") 

22 

23 return sp.csc_matrix(M) 

24 

25 

26def getBoussinesq2DMatrix(N, h, bc_hor, bc_ver, c_s, Nfreq, order): 

27 Dx_u, Dz_u = get2DMatrix(N, h, bc_hor[0], bc_ver[0], order) 

28 Dx_w, Dz_w = get2DMatrix(N, h, bc_hor[1], bc_ver[1], order) 

29 # Dx_b, Dz_b = get2DMatrix(N, h, bc_hor[2], bc_ver[2], order) 

30 Dx_p, Dz_p = get2DMatrix(N, h, bc_hor[3], bc_ver[3], order) 

31 

32 # Id_N = sp.eye(N[0] * N[1]) 

33 

34 Zero = np.zeros((N[0] * N[1], N[0] * N[1])) 

35 Id_w = sp.eye(N[0] * N[1]) 

36 

37 # Note: Bring all terms to right hand side, therefore a couple of minus signs 

38 # are needed 

39 

40 M1 = sp.hstack((Zero, Zero, Zero, -Dx_p), format="csr") 

41 M2 = sp.hstack((Zero, Zero, Id_w, -Dz_p), format="csr") 

42 M3 = sp.hstack((Zero, -(Nfreq**2) * Id_w, Zero, Zero), format="csr") 

43 M4 = sp.hstack((-(c_s**2) * Dx_u, -(c_s**2) * Dz_w, Zero, Zero), format="csr") 

44 M = sp.vstack((M1, M2, M3, M4), format="csr") 

45 

46 Id = sp.eye(4 * N[0] * N[1]) 

47 

48 return sp.csc_matrix(Id), sp.csc_matrix(M)