research/incomp-flame-post
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

pycompact - more informative test messages

Browse source
This commit is contained in:
ignis 2021-03-05 10:27:18 +09:00
parent ca3920c8f9
commit 52c2d1f574
1 changed files with 33 additions and 52 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

85
code/pycompact/pycompact.py
View file

@ -114,9 +114,9 @@ class CompactScheme:
dx_fortran = compact.rhs1np(self.hx, x) dx_fortran = compact.rhs1np(self.hx, x)
print (np.linalg.norm((dx - dx_fortran) / dx_fortran)) print ("RelError Norm: ", np.linalg.norm((dx - dx_fortran) / dx_fortran))
print (((dx - dx_fortran) / dx_fortran).min()) print ("RelError Min : ", ((dx - dx_fortran) / dx_fortran).min())
print (((dx - dx_fortran) / dx_fortran).max()) print ("RelError Min : ", ((dx - dx_fortran) / dx_fortran).max())
@ -130,16 +130,14 @@ class CompactScheme:
compact.tdslv(dx_fortran,l1) compact.tdslv(dx_fortran,l1)
print ("dx - exact") print ("dx - exact")
print (np.linalg.norm((dx - exact) / exact)) print ("RelError Norm: ", np.linalg.norm((dx - exact) / exact))
print ("RelError Min : ", ((dx - exact) / exact).min())
print (((dx - exact) / exact).min()) print ("RelError Min : ", ((dx - exact) / exact).max())
print (((dx - exact) / exact).max())
print ("dx_fortran - exact") print ("dx_fortran - exact")
print (np.linalg.norm((dx_fortran - exact) / exact)) print ("RelError Norm: ", np.linalg.norm((dx_fortran - exact) / exact))
print ("RelError Min : ", ((dx_fortran - exact) / exact).min())
print (((dx_fortran - exact) / exact).min()) print ("RelError Min : ", ((dx_fortran - exact) / exact).max())
print (((dx_fortran - exact) / exact).max())
''' '''
import pprint import pprint
@ -165,7 +163,7 @@ class CompactScheme:
coef_verify = self.py_stdlu(aa) coef_verify = self.py_stdlu(aa)
print (np.linalg.norm((coef - coef_verify)/coef_verify)) print ("RelError Norm: ", np.linalg.norm((coef - coef_verify)/coef_verify))
def verify_nonp_lud2(self): def verify_nonp_lud2(self):
@ -183,7 +181,7 @@ class CompactScheme:
coef_verify = self.py_stdlu(aa) coef_verify = self.py_stdlu(aa)
print (np.linalg.norm((coef - coef_verify)/coef_verify)) print ("RelError Norm: ", np.linalg.norm((coef - coef_verify)/coef_verify))
def py_stdlu(self, aa): def py_stdlu(self, aa):
@ -211,13 +209,11 @@ class CompactScheme:
if self.px: # Periodic BC if self.px: # Periodic BC
for i in range(nz): for i in range(nz):
xsrc[:] = src[i] dst[i] = compact.dfp(self.hx, src[i], 1)
dst[i] = compact.dfp(self.hx, xsrc, 1)
else: else:
for i in range(nz): for i in range(nz):
xsrc[:] = src[i] dst[i] = compact.dfnonp(self.hx, src[i], 1)
dst[i] = compact.dfnonp(self.hx, xsrc, 1)
# return np.swapaxes(dst, 1, 2) # return np.swapaxes(dst, 1, 2)
return dst return dst
@ -236,13 +232,11 @@ class CompactScheme:
if self.py: # Periodic BC if self.py: # Periodic BC
for i in range(nz): for i in range(nz):
xysrc = src[i].reshape(-1).reshape((nx,ny), order='F') dst[i] = compact.dfp(self.hx, src[i].T, 2).T
dst[i] = compact.dfp(self.hx, xysrc, 2).reshape(-1).reshape((ny,nx), order='C')
else: else:
for i in range(nz): for i in range(nz):
xysrc = src[i].reshape(-1).reshape((nx,ny), order='F') dst[i] = compact.dfnonp(self.hx, src[i].T, 2).T
dst[i] = compact.dfnonp(self.hx, xysrc, 2).reshape(-1).reshape((ny,nx), order='C')
# return np.swapaxes(dst, 1, 2) # return np.swapaxes(dst, 1, 2)
return dst return dst
@ -254,20 +248,17 @@ class CompactScheme:
nz, ny, nx = self.shape nz, ny, nx = self.shape
xsrc = np.zeros((nz, nx,), dtype=np.float64) # , order="F")
# dst = np.zeros((nx, ny, nz,), order="F") # dst = np.zeros((nx, ny, nz,), order="F")
dst = np.zeros((nz, ny, nx,), dtype=np.float64,) dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if self.pz: # Periodic BC if self.pz: # Periodic BC
for i in range(ny): for i in range(ny):
xsrc[:] = src[:,i,:] dst[:,i,:] = compact.dfp(self.hx, src[:,i,:], 3)
dst[i] = compact.dfp(self.hx, xsrc, 3)
else: else:
for i in range(ny): for i in range(ny):
xsrc[:] = src[:,i,:] dst[:,i,:] = compact.dfnonp(self.hx, src[:,i,:], 3)
dst[i] = compact.dfnonp(self.hx, xsrc, 3)
# return np.swapaxes(dst, 1, 2) # return np.swapaxes(dst, 1, 2)
return dst return dst
@ -377,14 +368,8 @@ def read_data (fname):
V = np.fromfile(f1, dtype=np.float64, count=(3*count)).reshape((3,nz,ny,nx)) V = np.fromfile(f1, dtype=np.float64, count=(3*count)).reshape((3,nz,ny,nx))
s = np.fromfile(f1, dtype=np.float64, count=(2*count)).reshape((2,nz,ny,nx)) s = np.fromfile(f1, dtype=np.float64, count=(2*count)).reshape((2,nz,ny,nx))
print (V.order) print (V.flags)
print (s.order) print (s.flags)
print (V.shape)
print (s.shape)
V.order="F"
s.order="F"
print (V.shape) print (V.shape)
print (s.shape) print (s.shape)
@ -443,11 +428,11 @@ def validate_trigonometric():
cos_fortran = compact.dfnonp(hxp, np.sin(1.1*XX).reshape((1,-1)), 1) cos_fortran = compact.dfnonp(hxp, np.sin(1.1*XX).reshape((1,-1)), 1)
cos_exact = 1.1 * np.cos(1.1*XX).reshape((1,-1)) cos_exact = 1.1 * np.cos(1.1*XX).reshape((1,-1))
print (cos_fortran.min(), cos_fortran.max()) print ("Compact Scheme: ", cos_fortran.min(), cos_fortran.max())
print (cos_exact.min(), cos_exact.max()) print ("Exact : ",cos_exact.min(), cos_exact.max())
print (np.linalg.norm((cos_fortran - cos_exact)/cos_exact)) print ("Norm of relative errors: ", np.linalg.norm((cos_fortran - cos_exact)/cos_exact))
# print (((cos_fortran - cos_exact)/cos_exact)) # print (((cos_fortran - cos_exact)/cos_exact))
@ -510,30 +495,26 @@ def validate_trigonometric():
def test_dns_data(): def test_dns_data():
import sys import sys
file_name = sys.argv[1] file_name = "./fort.1000"
t, nx, ny, nz, u, v, w, Y0, Y1 = read_data(file_name) answer = "export-ddxc.dat"
t, nx, ny, nz, u, v, w, Y0, Y1 = read_data(file_name)
ddx_answer = np.fromfile(answer, dtype=np.float64, count=Y1.size).reshape(Y1.shape)
cs = CompactScheme(nx, ny, nz, False, True, True, 4, 2, 2) cs = CompactScheme(nx, ny, nz, False, True, True, 4, 2, 2)
y = np.memmap("yr", dtype=np.float64, mode="w+", shape=cs.shape) ddx = cs.ddx(1.-Y1)
y[:] = Y1[:]
print (y.min(), y.max()) print (ddx.min(), ddx.max())
print (ddx_answer.min(), ddx_answer.max())
dydx = np.memmap("dyr", dtype=np.float64, mode="w+", shape=cs.shape)
dydx[:] = cs.ddx(Y1)[:]
dydx.flush()
print (dydx.min(), dydx.max())
if __name__ == "__main__": if __name__ == "__main__":
validate_trigonometric() validate_trigonometric()
# test_dns_data()
test_dns_data()