FFT was at single precision only
This commit is contained in:
parent
4d1568aa51
commit
1f160dfd72
|
@ -20,7 +20,7 @@ def curlFFT(geomdim,field):
|
||||||
dataType = 'tensor'
|
dataType = 'tensor'
|
||||||
|
|
||||||
field_fourier=np.fft.fftpack.rfftn(field,axes=(0,1,2))
|
field_fourier=np.fft.fftpack.rfftn(field,axes=(0,1,2))
|
||||||
curl_fourier=np.zeros(field_fourier.shape,'c8')
|
curl_fourier=np.zeros(field_fourier.shape,'c16')
|
||||||
|
|
||||||
# differentiation in Fourier space
|
# differentiation in Fourier space
|
||||||
k_s=np.zeros([3],'i')
|
k_s=np.zeros([3],'i')
|
||||||
|
@ -34,7 +34,7 @@ def curlFFT(geomdim,field):
|
||||||
for k in xrange(grid[2]/2+1):
|
for k in xrange(grid[2]/2+1):
|
||||||
k_s[2] = k
|
k_s[2] = k
|
||||||
if(k > grid[2]/2 ): k_s[2] = k_s[2] - grid[2]
|
if(k > grid[2]/2 ): k_s[2] = k_s[2] - grid[2]
|
||||||
xi=np.array([k_s[2]/geomdim[2]+0.0j,k_s[1]/geomdim[1]+0.j,k_s[0]/geomdim[0]+0.j],'c8')
|
xi=np.array([k_s[2]/geomdim[2]+0.0j,k_s[1]/geomdim[1]+0.j,k_s[0]/geomdim[0]+0.j],'c16')
|
||||||
if dataType == 'tensor':
|
if dataType == 'tensor':
|
||||||
for l in xrange(3):
|
for l in xrange(3):
|
||||||
curl_fourier[i,j,k,0,l] = ( field_fourier[i,j,k,l,2]*xi[1]\
|
curl_fourier[i,j,k,0,l] = ( field_fourier[i,j,k,l,2]*xi[1]\
|
||||||
|
|
|
@ -18,10 +18,10 @@ def divFFT(geomdim,field):
|
||||||
|
|
||||||
if len(np.shape(field)) == 4:
|
if len(np.shape(field)) == 4:
|
||||||
dataType = 'vector'
|
dataType = 'vector'
|
||||||
div_fourier=np.zeros(field_fourier.shape[0:3],'c8') # div is a scalar
|
div_fourier=np.zeros(field_fourier.shape[0:3],'c16') # div is a scalar
|
||||||
elif len(np.shape(field)) == 5:
|
elif len(np.shape(field)) == 5:
|
||||||
dataType = 'tensor'
|
dataType = 'tensor'
|
||||||
div_fourier=np.zeros(field_fourier.shape[0:4],'c8') # div is a vector
|
div_fourier=np.zeros(field_fourier.shape[0:4],'c16') # div is a vector
|
||||||
|
|
||||||
# differentiation in Fourier space
|
# differentiation in Fourier space
|
||||||
k_s=np.zeros([3],'i')
|
k_s=np.zeros([3],'i')
|
||||||
|
@ -35,7 +35,7 @@ def divFFT(geomdim,field):
|
||||||
for k in xrange(grid[2]/2+1):
|
for k in xrange(grid[2]/2+1):
|
||||||
k_s[2] = k
|
k_s[2] = k
|
||||||
if(k > grid[2]/2 ): k_s[2] = k_s[2] - grid[2]
|
if(k > grid[2]/2 ): k_s[2] = k_s[2] - grid[2]
|
||||||
xi=np.array([k_s[2]/geomdim[2]+0.0j,k_s[1]/geomdim[1]+0.j,k_s[0]/geomdim[0]+0.j],'c8')
|
xi=np.array([k_s[2]/geomdim[2]+0.0j,k_s[1]/geomdim[1]+0.j,k_s[0]/geomdim[0]+0.j],'c16')
|
||||||
if dataType == 'tensor':
|
if dataType == 'tensor':
|
||||||
for l in xrange(3):
|
for l in xrange(3):
|
||||||
div_fourier[i,j,k,l] = sum(field_fourier[i,j,k,l,0:3]*xi) *TWOPIIMG
|
div_fourier[i,j,k,l] = sum(field_fourier[i,j,k,l,0:3]*xi) *TWOPIIMG
|
||||||
|
|
Loading…
Reference in New Issue