diff --git a/code/Compact.f90 b/code/Compact.f90 index 46b292b..0fd1a9e 100644 --- a/code/Compact.f90 +++ b/code/Compact.f90 @@ -155,6 +155,51 @@ END SUBROUTINE ludcmp_calculate + + SUBROUTINE test_nonp_lud1(xx, coef) + INTEGER :: xx + REAL(KIND=8), DIMENSION(xx) :: aa + REAL(KIND=8), DIMENSION(xx), INTENT(OUT) :: coef + aa=3. + aa(1)=0.5 ; aa(2)=4. + aa(xx-1)=4. ; aa(xx)=0.5 + + CALL stdlu(aa,xx,coef) + END SUBROUTINE test_nonp_lud1 + + + SUBROUTINE test_nonp_lud2(xx, coef) + INTEGER :: xx + REAL(KIND=8), DIMENSION(xx) :: aa + REAL(KIND=8), DIMENSION(xx), INTENT(OUT) :: coef + aa=5.5 + aa(1)=2./11. ; aa(2)=10. + aa(xx-1)=10. ; aa(xx)=2./11. + + CALL stdlu(aa,xx,coef) + END SUBROUTINE test_nonp_lud2 + + + SUBROUTINE test_p_lud1(xx, coef1, coef2) + INTEGER :: xx + REAL(KIND=8) :: a + REAL(KIND=8), DIMENSION(xx), INTENT(OUT) :: coef1, coef2 + a=3. ! first derivative + + CALL ptdlu(a,xx,coef1,coef2) ! x-direction + END SUBROUTINE test_p_lud1 + + + SUBROUTINE test_p_lud2(xx, coef1, coef2) + INTEGER :: xx + REAL(KIND=8) :: a + REAL(KIND=8), DIMENSION(xx), INTENT(OUT) :: coef1, coef2 + a=11./2. ! second derivative + + CALL ptdlu(a,xx,coef1,coef2) ! x-direction + END SUBROUTINE test_p_lud2 + + SUBROUTINE nonp_lud(xyz,xx) INTEGER :: i,xyz,xx REAL(KIND=8), DIMENSION(xx) :: aa @@ -313,7 +358,7 @@ SUBROUTINE ptdslv(r,n,l,w,nd) INTEGER,INTENT(IN) :: n,nd REAL(KIND=8),INTENT(INOUT),DIMENSION(nd,n) :: r - REAL(KIND=8),INTENT(IN),DIMENSION(:) :: l,w + REAL(KIND=8),INTENT(IN),DIMENSION(n) :: l,w INTEGER i,j REAL(KIND=8), DIMENSION(nd) :: sum DO j=1,nd @@ -406,7 +451,7 @@ SUBROUTINE tdslv(r,n,l,nd) INTEGER,INTENT(IN) :: n,nd REAL(KIND=8),INTENT(INOUT),DIMENSION(nd,n) :: r - REAL(KIND=8),INTENT(IN),DIMENSION(:) :: l + REAL(KIND=8),INTENT(IN),DIMENSION(n) :: l INTEGER i,j REAL(KIND=8) t1 DO j=1,nd diff --git a/code/pycompact/pycompact.py b/code/pycompact/pycompact.py index dd553a2..064f0b9 100644 --- a/code/pycompact/pycompact.py +++ b/code/pycompact/pycompact.py @@ -67,7 +67,7 @@ class CompactScheme: # dst = np.zeros((nx, ny, nz,), order="F") dst = np.zeros((nz, ny, nx,), dtype=np.float64,) - if self.px: + if self.px: # Periodic BC for i in range(nz): xsrc[:] = src[i] dst[i] = compact.dfp(self.hx, xsrc, 1) @@ -221,16 +221,18 @@ def validate_trigonometric(): pi8 = np.arccos(-1.) - print (pi8) + print ("calculated pi value = ", pi8) l_0 = 2.0 hyp=l_0*pi8/ny hxp=hyp hzp=hyp - Y1 = np.zeros(shape) - dY = np.zeros(shape) + cs = CompactScheme(nx, ny, nz, False, True, True, 4., 2., 2.) + + + Y1 = np.zeros(shape) true = np.zeros(shape) @@ -244,8 +246,6 @@ def validate_trigonometric(): Y1[:] = np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz)[:] true[:] = (1.1 * np.cos(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz))[:] - cs = CompactScheme(nx, ny, nz, False, True, True, 4., 2., 2.) - y = np.memmap("phi", dtype=np.float64, mode="w+", shape=cs.shape) y[:] = Y1[:]