From 873c52ccebdd306ee862a4a8af4f3553c3d3d405 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 7 Nov 2016 09:49:53 +0100 Subject: [PATCH] using 3 way merge to have syntax as similar as possible --- processing/post/addCurl.py | 16 ++++++++-------- processing/post/addDivergence.py | 16 ++++++++-------- processing/post/addGradient.py | 14 +++++++------- 3 files changed, 23 insertions(+), 23 deletions(-) diff --git a/processing/post/addCurl.py b/processing/post/addCurl.py index bec3e3a81..98a00197c 100755 --- a/processing/post/addCurl.py +++ b/processing/post/addCurl.py @@ -12,8 +12,8 @@ scriptID = ' '.join([scriptName,damask.version]) def curlFFT(geomdim,field): shapeFFT = np.array(np.shape(field))[0:3] grid = np.array(np.shape(field)[2::-1]) - N = grid.prod() # field size - n = np.array(np.shape(field)[3:]).prod() # data size + N = grid.prod() # field size + n = np.array(np.shape(field)[3:]).prod() # data size if n == 3: dataType = 'vector' elif n == 9: dataType = 'tensor' @@ -24,23 +24,23 @@ def curlFFT(geomdim,field): # differentiation in Fourier space TWOPIIMG = 2.0j*math.pi k_sk = np.where(np.arange(grid[2])>grid[2]//2,np.arange(grid[2])-grid[2],np.arange(grid[2]))/geomdim[0] - if grid[2]%2 == 0: k_sk[grid[2]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) + if grid[2]%2 == 0: k_sk[grid[2]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) k_sj = np.where(np.arange(grid[1])>grid[1]//2,np.arange(grid[1])-grid[1],np.arange(grid[1]))/geomdim[1] - if grid[1]%2 == 0: k_sj[grid[1]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) - + if grid[1]%2 == 0: k_sj[grid[1]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) + k_si = np.arange(grid[0]//2+1)/geomdim[2] kk, kj, ki = np.meshgrid(k_sk,k_sj,k_si,indexing = 'ij') k_s = np.concatenate((ki[:,:,:,None],kj[:,:,:,None],kk[:,:,:,None]),axis = 3).astype('c16') e = np.zeros((3, 3, 3)) - e[0, 1, 2] = e[1, 2, 0] = e[2, 0, 1] = 1.0 # Levi-Civita symbols + e[0, 1, 2] = e[1, 2, 0] = e[2, 0, 1] = 1.0 # Levi-Civita symbols e[0, 2, 1] = e[2, 1, 0] = e[1, 0, 2] = -1.0 - if dataType == 'tensor': # tensor, 3x3 -> 3x3 + if dataType == 'tensor': # tensor, 3x3 -> 3x3 curl_fourier = np.einsum('slm,ijkl,ijknm->ijksn',e,k_s,field_fourier)*TWOPIIMG - elif dataType == 'vector': # vector, 3 -> 3 + elif dataType == 'vector': # vector, 3 -> 3 curl_fourier = np.einsum('slm,ijkl,ijkm->ijks',e,k_s,field_fourier)*TWOPIIMG return np.fft.irfftn(curl_fourier,axes=(0,1,2),s=shapeFFT).reshape([N,n]) diff --git a/processing/post/addDivergence.py b/processing/post/addDivergence.py index 84cdf4593..232b5bc21 100755 --- a/processing/post/addDivergence.py +++ b/processing/post/addDivergence.py @@ -12,30 +12,30 @@ scriptID = ' '.join([scriptName,damask.version]) def divFFT(geomdim,field): shapeFFT = np.array(np.shape(field))[0:3] grid = np.array(np.shape(field)[2::-1]) - N = grid.prod() # field size - n = np.array(np.shape(field)[3:]).prod() # data size + N = grid.prod() # field size + n = np.array(np.shape(field)[3:]).prod() # data size if n == 3: dataType = 'vector' elif n == 9: dataType = 'tensor' - + field_fourier = np.fft.rfftn(field,axes=(0,1,2),s=shapeFFT) - div_fourier = np.empty(field_fourier.shape[0:len(np.shape(field))-1],'c16') # size depents on whether tensor or vector + div_fourier = np.empty(field_fourier.shape[0:len(np.shape(field))-1],'c16') # differentiation in Fourier space TWOPIIMG = 2.0j*math.pi k_sk = np.where(np.arange(grid[2])>grid[2]//2,np.arange(grid[2])-grid[2],np.arange(grid[2]))/geomdim[0] - if grid[2]%2 == 0: k_sk[grid[2]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) + if grid[2]%2 == 0: k_sk[grid[2]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) k_sj = np.where(np.arange(grid[1])>grid[1]//2,np.arange(grid[1])-grid[1],np.arange(grid[1]))/geomdim[1] - if grid[1]%2 == 0: k_sj[grid[1]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) + if grid[1]%2 == 0: k_sj[grid[1]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) k_si = np.arange(grid[0]//2+1)/geomdim[2] kk, kj, ki = np.meshgrid(k_sk,k_sj,k_si,indexing = 'ij') k_s = np.concatenate((ki[:,:,:,None],kj[:,:,:,None],kk[:,:,:,None]),axis = 3).astype('c16') - if dataType == 'tensor': # tensor, 3x3 -> 3 + if dataType == 'tensor': # tensor, 3x3 -> 3 div_fourier = np.einsum('ijklm,ijkm->ijkl',field_fourier,k_s)*TWOPIIMG - elif dataType == 'vector': # vector, 3 -> 1 + elif dataType == 'vector': # vector, 3 -> 1 div_fourier = np.einsum('ijkl,ijkl->ijk',field_fourier,k_s)*TWOPIIMG return np.fft.irfftn(div_fourier,axes=(0,1,2),s=shapeFFT).reshape([N,n/3]) diff --git a/processing/post/addGradient.py b/processing/post/addGradient.py index 5ac784c01..fefe8f84e 100755 --- a/processing/post/addGradient.py +++ b/processing/post/addGradient.py @@ -12,8 +12,8 @@ scriptID = ' '.join([scriptName,damask.version]) def gradFFT(geomdim,field): shapeFFT = np.array(np.shape(field))[0:3] grid = np.array(np.shape(field)[2::-1]) - N = grid.prod() # field size - n = np.array(np.shape(field)[3:]).prod() # data size + N = grid.prod() # field size + n = np.array(np.shape(field)[3:]).prod() # data size if n == 3: dataType = 'vector' elif n == 1: dataType = 'scalar' @@ -24,18 +24,18 @@ def gradFFT(geomdim,field): # differentiation in Fourier space TWOPIIMG = 2.0j*math.pi k_sk = np.where(np.arange(grid[2])>grid[2]//2,np.arange(grid[2])-grid[2],np.arange(grid[2]))/geomdim[0] - if grid[2]%2 == 0: k_sk[grid[2]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) + if grid[2]%2 == 0: k_sk[grid[2]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) k_sj = np.where(np.arange(grid[1])>grid[1]//2,np.arange(grid[1])-grid[1],np.arange(grid[1]))/geomdim[1] - if grid[1]%2 == 0: k_sj[grid[1]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) - + if grid[1]%2 == 0: k_sj[grid[1]//2] = 0 # for even grid, set Nyquist freq to 0 (Johnson, MIT, 2011) + k_si = np.arange(grid[0]//2+1)/geomdim[2] kk, kj, ki = np.meshgrid(k_sk,k_sj,k_si,indexing = 'ij') k_s = np.concatenate((ki[:,:,:,None],kj[:,:,:,None],kk[:,:,:,None]),axis = 3).astype('c16') - if dataType == 'vector': # vector, 3 -> 3x3 + if dataType == 'vector': # vector, 3 -> 3x3 grad_fourier = np.einsum('ijkl,ijkm->ijklm',field_fourier,k_s)*TWOPIIMG - elif dataType == 'scalar': # scalar, 1 -> 3 + elif dataType == 'scalar': # scalar, 1 -> 3 grad_fourier = np.einsum('ijkl,ijkl->ijkl',field_fourier,k_s)*TWOPIIMG return np.fft.irfftn(grad_fourier,axes=(0,1,2),s=shapeFFT).reshape([N,3*n])