pycompact - performance improvement, buffer reuse and x-dir tiling

This commit is contained in:
ignis 2021-08-12 15:48:23 +09:00
parent 7f6df479c1
commit 1d161dd900

View file

@ -7,6 +7,8 @@ class CompactScheme:
pi8 = np.arccos(-1., dtype=np.float64)
self.nb = 32
self.shape = (nz, ny, nx)
self.px = px
@ -195,34 +197,36 @@ class CompactScheme:
return coef
def ddx (self, src):
def ddx (self, src, dst=None):
if src.shape != self.shape:
print ("error")
if src.shape != self.shape: print ("error: src.shape wrong!")
nz, ny, nx = self.shape
nb = self.nb
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if dst is None or dst.shape != self.shape :
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if self.px: # Periodic BC
for i in range(nz):
dst[i] = compact.dfp(self.hx, src[i], 1)
for j in range(0, ny, nb):
dst[i,j:j+nb] = compact.dfp(self.hx, src[i,j:j+nb], 1)
else:
for i in range(nz):
dst[i] = compact.dfnonp(self.hx, src[i], 1)
for j in range(0, ny, nb):
dst[i,j:j+nb] = compact.dfnonp(self.hx, src[i,j:j+nb], 1)
# return np.swapaxes(dst, 1, 2)
return dst
def ddy (self, src):
def ddy (self, src, dst=None):
if src.shape != self.shape:
print ("error")
if src.shape != self.shape: print ("error: src.shape wrong!")
nz, ny, nx = self.shape
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if dst is None or dst.shape != self.shape :
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if self.py: # Periodic BC
for i in range(nz):
@ -232,53 +236,60 @@ class CompactScheme:
for i in range(nz):
dst[i] = compact.dfnonp(self.hx, src[i].T, 2).T
# return np.swapaxes(dst, 1, 2)
return dst
def ddz (self, src):
def ddz (self, src, dst=None):
if src.shape != self.shape:
print ("error")
if src.shape != self.shape: print ("error: src.shape wrong!")
nz, ny, nx = self.shape
# dst = np.zeros((nx, ny, nz,), order="F")
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if dst is None or dst.shape != self.shape :
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
zxsrc = np.zeros((nz, nx))
if self.pz: # Periodic BC
for i in range(ny):
dst[:,i,:] = compact.dfp(self.hx, src[:,i,:].T, 3).T
zxsrc[:,:] = src[:,i,:]
dst[:,i,:] = compact.dfp(self.hx, zxsrc.T, 3).T
else:
for i in range(ny):
dst[:,i,:] = compact.dfnonp(self.hx, src[:,i,:].T, 3).T
zxsrc[:,:] = src[:,i,:]
dst[:,i,:] = compact.dfnonp(self.hx, zxsrc.T, 3).T
# return np.swapaxes(dst, 1, 2)
return dst
def d2dx (self, src):
def d2dx (self, src, dst=None):
if src.shape != self.shape: print ("error")
if src.shape != self.shape: print ("error: src.shape wrong!")
nz, ny, nx = self.shape
nb = self.nb
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if dst is None or dst.shape != self.shape :
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if self.px: # Periodic BC
for i in range(nz):
dst[i] = compact.d2fp(self.hx, src[i], 1)
for j in range(0, ny, nb):
dst[i,j:j+nb] = compact.d2fp(self.hx, src[i,j:j+nb], 1)
else:
for i in range(nz):
dst[i] = compact.d2fnonp(self.hx, src[i], 1)
for j in range(0, ny, nb):
dst[i,j:j+nb] = compact.d2fnonp(self.hx, src[i,j:j+nb], 1)
return dst
def d2dy (self, src):
def d2dy (self, src, dst=None):
if src.shape != self.shape: print ("error")
if src.shape != self.shape: print ("error: src.shape wrong!")
nz, ny, nx = self.shape
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if dst is None or dst.shape != self.shape :
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if self.py: # Periodic BC
for i in range(nz):
@ -289,20 +300,25 @@ class CompactScheme:
return dst
def d2dz (self, src):
def d2dz (self, src, dst=None):
if src.shape != self.shape: print ("error")
if src.shape != self.shape: print ("error: src.shape wrong!")
nz, ny, nx = self.shape
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
if dst is None or dst.shape != self.shape :
dst = np.zeros((nz, ny, nx,), dtype=np.float64,)
zxsrc = np.zeros((nz, nx))
if self.pz: # Periodic BC
for i in range(ny):
dst[:,i,:] = compact.d2fp(self.hx, src[:,i,:].T, 3).T
zxsrc[:,:] = src[:,i,:]
dst[:,i,:] = compact.d2fp(self.hx, zxsrc.T, 3).T
else:
for i in range(ny):
dst[:,i,:] = compact.d2fnonp(self.hx, src[:,i,:].T, 3).T
zxsrc[:,:] = src[:,i,:]
dst[:,i,:] = compact.d2fnonp(self.hx, zxsrc.T, 3).T
return dst
@ -513,7 +529,7 @@ def validate_trigonometric():
print(" DDX Test ")
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))[:]
dY1 = cs.ddx(Y1)[:]
dY1 = cs.ddx(Y1, dY1)[:]
compare_3d_result(true, dY1)
@ -521,7 +537,7 @@ def validate_trigonometric():
print("D2DX Test ")
Y1[:] = np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz)[:]
true[:] = (- 1.1 * 1.1 * np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz))[:]
dY1 = cs.d2dx(Y1)[:]
dY1 = cs.d2dx(Y1, dY1)[:]
compare_3d_result(true, dY1)
@ -529,7 +545,7 @@ def validate_trigonometric():
print(" DDY Test ")
Y1[:] = np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz)[:]
true[:] = (3.0 * np.sin(1.1 * xx) * np.cos(3.0 * yy) * np.sin(2.0 * zz))[:]
dY1 = cs.ddy(Y1)[:]
dY1 = cs.ddy(Y1, dY1)[:]
compare_3d_result(true, dY1)
@ -537,7 +553,7 @@ def validate_trigonometric():
print(" DDZ Test ")
Y1[:] = np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz)[:]
true[:] = (2.0 * np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.cos(2.0 * zz))[:]
dY1 = cs.ddz(Y1)[:]
dY1 = cs.ddz(Y1, dY1)[:]
compare_3d_result(true, dY1)
@ -576,6 +592,7 @@ def performance_benchmark():
print ("Test ddx")
Y1 = np.zeros(shape)
dY1 = np.zeros(shape)
true = np.zeros(shape)
@ -605,7 +622,7 @@ def performance_benchmark():
true[:] = (1.1 * np.cos(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz))[:]
t0 = time.perf_counter()
for i in range(100):
dY1 = cs.ddx(Y1)[:]
cs.ddx(Y1, dY1)
t1 = time.perf_counter()
print ("Elapsed time in sec for 100 repetition: ", t1-t0)
compare_3d_result(true, dY1)
@ -614,7 +631,7 @@ def performance_benchmark():
true[:] = (- 1.1 * 1.1 * np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz))[:]
t0 = time.perf_counter()
for i in range(100):
dY1 = cs.d2dx(Y1)[:]
cs.d2dx(Y1, dY1)
t1 = time.perf_counter()
print ("Elapsed time in sec for 100 repetition: ", t1-t0)
compare_3d_result(true, dY1)
@ -623,7 +640,7 @@ def performance_benchmark():
true[:] = (3.0 * np.sin(1.1 * xx) * np.cos(3.0 * yy) * np.sin(2.0 * zz))[:]
t0 = time.perf_counter()
for i in range(100):
dY1 = cs.ddy(Y1)[:]
cs.ddy(Y1, dY1)
t1 = time.perf_counter()
print ("Elapsed time in sec for 100 repetition: ", t1-t0)
compare_3d_result(true, dY1)
@ -632,7 +649,7 @@ def performance_benchmark():
true[:] = (-3. * 3.0 * np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz))[:]
t0 = time.perf_counter()
for i in range(100):
dY1 = cs.d2dy(Y1)[:]
cs.d2dy(Y1, dY1)
t1 = time.perf_counter()
print ("Elapsed time in sec for 100 repetition: ", t1-t0)
compare_3d_result(true, dY1)
@ -641,7 +658,7 @@ def performance_benchmark():
true[:] = (2.0 * np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.cos(2.0 * zz))[:]
t0 = time.perf_counter()
for i in range(100):
dY1 = cs.ddz(Y1)[:]
cs.ddz(Y1, dY1)
t1 = time.perf_counter()
print ("Elapsed time in sec for 100 repetition: ", t1-t0)
compare_3d_result(true, dY1)
@ -651,7 +668,7 @@ def performance_benchmark():
true[:] = (-2. * 2.0 * np.sin(1.1 * xx) * np.sin(3.0 * yy) * np.sin(2.0 * zz))[:]
t0 = time.perf_counter()
for i in range(100):
dY1 = cs.d2dz(Y1)[:]
cs.d2dz(Y1, dY1)
t1 = time.perf_counter()
print ("Elapsed time in sec for 100 repetition: ", t1-t0)
compare_3d_result(true, dY1)