Skip to content
Snippets Groups Projects
Commit a2ca4a0019fd authored by Pierre Augier's avatar Pierre Augier
Browse files

Prepare fluidfft/microbench

parent d75b3a04d8c0
No related branches found
No related tags found
No related merge requests found
......@@ -18,6 +18,8 @@
*.bbl
*.blg
*.so
*/*.pstats
*/.ipynb_checkpoints
......
perf: perfpythran
clean: cleanfortran
rm -f *.so
cleanfortran:
rm -f bench_proj_fortran.out
code = 'func(c0, c1, c2, a0, a1, a2, a3)'
end=as func, c0, c1, c2, a0, a1, a2, a3
perfnumpy: proj_pythran.so
# numpy proj
@python -m perf timeit -s 'from bench import proj0 $(end)' $(code)
@# python -m perf timeit -s 'from bench import proj_loop0 $(end)' $(code)
# numpy proj_inplace
@python -m perf timeit -s 'from bench import proj_inplace0 $(end)' $(code)
@# python -m perf timeit -s 'from bench import proj_inplace_loop0 $(end)' $(code)
perfpythran: proj_pythran.so
# pythran proj_loop
@python -m perf timeit -s 'from bench import proj_loop $(end)' $(code)
# pythran proj_inplace_loop
@python -m perf timeit -s 'from bench import proj_inplace_loop $(end)' $(code)
# pythran proj
@python -m perf timeit -s 'from bench import proj $(end)' $(code)
# pythran proj_inplace
@python -m perf timeit -s 'from bench import proj_inplace $(end)' $(code)
perfnumba:
# numba proj
@python -m perf timeit -s 'from bench import proj_numba $(end)' $(code)
# numba proj_inplace
@python -m perf timeit -s 'from bench import proj_inplace_numba $(end)' $(code)
# numba proj_loop
@python -m perf timeit -s 'from bench import proj_loop_numba $(end)' $(code)
# numba proj_inplace_loop
@python -m perf timeit -s 'from bench import proj_inplace_loop_numba $(end)' $(code)
proj_pythran.so: proj.py
pythran -v proj.py -march=native -o proj_pythran.so
perffortran: bench_proj_fortran.out bench_proj_fortran_inplace.out
./bench_proj_fortran.out
./bench_proj_fortran_inplace.out
bench_proj_fortran.out: bench_proj_fortran.f90
gfortran bench_proj_fortran.f90 -O3 -o bench_proj_fortran.out
bench_proj_fortran_inplace.out: bench_proj_fortran_inplace.f90
gfortran bench_proj_fortran_inplace.f90 -O3 -o bench_proj_fortran_inplace.out
Done on my (pa) laptop (see sys_info_laptop-pierre-kth.xml).
I used:
```
sudo python -m perf system tune
```
## outplace
### fortran with allocate/deallocate
Mean Time = 22.1 ms
### fortran without allocate/deallocate
Mean Time = 11.3 ms
### numpy proj
Mean +- std dev: 66.6 ms +- 0.7 ms
### pythran proj
Mean +- std dev: 38.8 ms +- 3.1 ms
### pythran proj_loop
Mean +- std dev: 19.8 ms +- 1.4 ms
### numba proj
Mean +- std dev: 78.3 ms +- 2.0 ms
### numba proj_loop
Mean +- std dev: 26.5 ms +- 2.5 ms
## inplace
### fortran inplace
Mean Time = 9.0 ms
### numpy proj_inplace
Mean +- std dev: 54.2 ms +- 2.4 ms
### pythran proj_inplace
Mean +- std dev: 18.7 ms +- 0.8 ms
### pythran proj_inplace_loop
Mean +- std dev: 8.60 ms +- 0.08 ms
### numba proj_inplace
Mean +- std dev: 65.9 ms +- 1.2 ms
### numba proj_inplace_loop
Mean +- std dev: 16.3 ms +- 1.5 ms
import numpy as np
from proj import (
proj as proj0,
proj_loop as proj_loop0,
proj_inplace as proj_inplace0,
proj_inplace_loop as proj_inplace_loop0
)
from proj_pythran import (
proj,
proj_loop,
proj_inplace,
proj_inplace_loop
)
import proj_pythran
from proj_numba import (
proj as proj_numba,
proj_loop as proj_loop_numba,
proj_inplace as proj_inplace_numba,
proj_inplace_loop as proj_inplace_loop_numba
)
assert hasattr(proj_pythran, '__pythran__')
__all__ = [
'proj0',
'proj_loop0',
'proj_inplace0',
'proj_inplace_loop0',
'proj',
'proj_loop',
'proj_inplace',
'proj_inplace_loop',
'proj_numba',
'proj_loop_numba',
'proj_inplace_numba',
'proj_inplace_loop_numba',
]
n0 = n1 = n2 = 128
shape = (n0, n1, n2//2+1)
c0 = 1.3j + np.ones(shape, dtype=np.complex128)
c1 = 2.3j + np.ones(shape, dtype=np.complex128)
c2 = 3.3j + np.ones(shape, dtype=np.complex128)
a0 = np.ones(shape)
a1 = np.ones(shape)
a2 = np.ones(shape)
a3 = np.ones(shape)
# JIT compiler
for func in (proj_numba,
proj_loop_numba,
proj_inplace_numba,
proj_inplace_loop_numba):
func(c0, c1, c2, a0, a1, a2, a3)
! benchmark of the function proj
program benchmark_proj
implicit none
integer, parameter :: N0=128, N1=128, N2=64, N=100
double precision, allocatable :: vx(:,:,:,:), vy(:,:,:,:), vz(:,:,:,:)
double precision, allocatable :: kx(:,:,:), ky(:,:,:), kz(:,:,:)
double precision, allocatable :: inv_k_square_nozero(:,:,:)
double precision, allocatable :: res(:,:,:,:,:)
real :: start, finish, cumtime
integer :: i
allocate(vx(2, N2, N1, N0), vy(2, N2, N1, N0), vz(2, N2, N1, N0))
allocate(kx(N2, N1, N0), ky(N2, N1, N0), kz(N2, N1, N0))
allocate(inv_k_square_nozero(N2, N1, N0))
call random_number(vx)
call random_number(vy)
call random_number(vz)
call random_number(kx)
call random_number(ky)
call random_number(kz)
call random_number(inv_k_square_nozero)
cumtime = 0
print*, "This program make some calculations."
do i = 1, N
call cpu_time(start)
allocate(res(2, 3, N2, N1, N0))
call proj(res, vx, vy, vz, kx, ky, kz, inv_k_square_nozero, N0, N1, N2)
deallocate(res)
call cpu_time(finish)
cumtime = cumtime + finish - start
enddo
print '("Mean Time = ",f6.3," ms")', 1000*cumtime/N
cumtime = 0
print*, "without allocate/deallocate."
allocate(res(2, 3, N2, N1, N0))
do i = 1, N
call cpu_time(start)
call proj(res, vx, vy, vz, kx, ky, kz, inv_k_square_nozero, N0, N1, N2)
call cpu_time(finish)
cumtime = cumtime + finish - start
enddo
print '("Mean Time = ",f6.3," ms")', 1000*cumtime/N
deallocate(res)
deallocate(vx, vy, vz)
deallocate(kx, ky, kz)
deallocate(inv_k_square_nozero)
end program benchmark_proj
subroutine proj(res, vx, vy, vz, kx, ky, kz, inv_k_square_nozero, N0, N1, N2)
implicit none
! Input/Output
integer, intent(in) :: N0, N1, N2
double precision, intent(in) :: vx(2, N2, N1, N0), vy(2, N2, N1, N0), vz(2, N2, N1, N0)
double precision, intent(in) :: kx(N2, N1, N0), ky(N2, N1, N0), kz(N2, N1, N0)
double precision, intent(in) :: inv_k_square_nozero(N2, N1, N0)
double precision, intent(out) :: res(2, 3, N2, N1, N0)
! Locals
double precision :: tmp(2)
integer:: i, j, k
do k = 1, N0
do j = 1, N1
do i = 1, N2
tmp(1:2) = (kx(i,j,k) * vx(1:2,i,j,k) &
+ ky(i,j,k) * vy(1:2,i,j,k) &
+ kz(i,j,k) * vz(1:2,i,j,k)) * inv_k_square_nozero(i,j,k)
res(1:2,1,i,j,k) = vx(1:2,i,j,k) - kx(i,j,k) * tmp(1:2)
res(1:2,2,i,j,k) = vy(1:2,i,j,k) - ky(i,j,k) * tmp(1:2)
res(1:2,3,i,j,k) = vz(1:2,i,j,k) - kz(i,j,k) * tmp(1:2)
enddo
enddo
enddo
end subroutine proj
! benchmark of the function proj
program benchmark_proj
implicit none
integer, parameter :: N0=128, N1=128, N2=64, N=100
double precision, allocatable :: vx(:,:,:,:), vy(:,:,:,:), vz(:,:,:,:)
double precision, allocatable :: kx(:,:,:), ky(:,:,:), kz(:,:,:)
double precision, allocatable :: inv_k_square_nozero(:,:,:)
real :: start, finish, cumtime
integer :: i
allocate(vx(2, N2, N1, N0), vy(2, N2, N1, N0), vz(2, N2, N1, N0))
allocate(kx(N2, N1, N0), ky(N2, N1, N0), kz(N2, N1, N0))
allocate(inv_k_square_nozero(N2, N1, N0))
call random_number(vx)
call random_number(vy)
call random_number(vz)
call random_number(kx)
call random_number(ky)
call random_number(kz)
call random_number(inv_k_square_nozero)
cumtime = 0
print*, "This program make some calculations."
do i = 1, N
call cpu_time(start)
call proj(vx, vy, vz, kx, ky, kz, inv_k_square_nozero, N0, N1, N2)
call cpu_time(finish)
cumtime = cumtime + finish - start
enddo
print '("Mean Time = ",f6.3," ms")', 1000*cumtime/N
deallocate(vx, vy, vz)
deallocate(kx, ky, kz)
deallocate(inv_k_square_nozero)
end program benchmark_proj
subroutine proj(vx, vy, vz, kx, ky, kz, inv_k_square_nozero, N0, N1, N2)
implicit none
! Input/Output
integer, intent(in) :: N0, N1, N2
double precision, intent(inout) :: vx(2, N2, N1, N0), vy(2, N2, N1, N0), vz(2, N2, N1, N0)
double precision, intent(in) :: kx(N2, N1, N0), ky(N2, N1, N0), kz(N2, N1, N0)
double precision, intent(in) :: inv_k_square_nozero(N2, N1, N0)
! Locals
double precision :: tmp(2)
integer:: i, j, k
do k = 1, N0
do j = 1, N1
do i = 1, N2
tmp(1:2) = (kx(i,j,k) * vx(1:2,i,j,k) &
+ ky(i,j,k) * vy(1:2,i,j,k) &
+ kz(i,j,k) * vz(1:2,i,j,k) &
) * inv_k_square_nozero(i,j,k)
vx(1:2,i,j,k) = vx(1:2,i,j,k) - kx(i,j,k) * tmp
vy(1:2,i,j,k) = vy(1:2,i,j,k) - ky(i,j,k) * tmp
vz(1:2,i,j,k) = vz(1:2,i,j,k) - kz(i,j,k) * tmp
enddo
enddo
enddo
end subroutine proj
import numpy as np
# pythran export proj(
# complex128[][][], complex128[][][], complex128[][][],
# float64[][][], float64[][][], float64[][][], float64[][][])
def proj(vx, vy, vz, kx, ky, kz, inv_k_square_nozero):
tmp = (kx * vx + ky * vy + kz * vz) * inv_k_square_nozero
return (vx - kx * tmp,
vy - ky * tmp,
vz - kz * tmp)
# pythran export proj_loop(
# complex128[][][], complex128[][][], complex128[][][],
# float64[][][], float64[][][], float64[][][], float64[][][])
def proj_loop(vx, vy, vz, kx, ky, kz, inv_k_square_nozero):
rx = np.empty_like(vx)
ry = np.empty_like(vx)
rz = np.empty_like(vx)
n0, n1, n2 = kx.shape
for i0 in range(n0):
for i1 in range(n1):
for i2 in range(n2):
tmp = (kx[i0, i1, i2] * vx[i0, i1, i2]
+ ky[i0, i1, i2] * vy[i0, i1, i2]
+ kz[i0, i1, i2] * vz[i0, i1, i2]
) * inv_k_square_nozero[i0, i1, i2]
rx[i0, i1, i2] = vx[i0, i1, i2] - kx[i0, i1, i2] * tmp
ry[i0, i1, i2] = vz[i0, i1, i2] - kx[i0, i1, i2] * tmp
rz[i0, i1, i2] = vy[i0, i1, i2] - kx[i0, i1, i2] * tmp
return rx, ry, rz
# pythran export proj_inplace(
# complex128[][][], complex128[][][], complex128[][][],
# float64[][][], float64[][][], float64[][][], float64[][][])
def proj_inplace(vx, vy, vz, kx, ky, kz, inv_k_square_nozero):
tmp = (kx * vx + ky * vy + kz * vz) * inv_k_square_nozero
vx -= kx * tmp
vy -= ky * tmp
vz -= kz * tmp
# pythran export proj_inplace_loop(
# complex128[][][], complex128[][][], complex128[][][],
# float64[][][], float64[][][], float64[][][], float64[][][])
def proj_inplace_loop(vx, vy, vz, kx, ky, kz, inv_k_square_nozero):
n0, n1, n2 = kx.shape
for i0 in range(n0):
for i1 in range(n1):
for i2 in range(n2):
tmp = (kx[i0, i1, i2] * vx[i0, i1, i2]
+ ky[i0, i1, i2] * vy[i0, i1, i2]
+ kz[i0, i1, i2] * vz[i0, i1, i2]
) * inv_k_square_nozero[i0, i1, i2]
vx[i0, i1, i2] -= kx[i0, i1, i2] * tmp
vy[i0, i1, i2] -= ky[i0, i1, i2] * tmp
vz[i0, i1, i2] -= kz[i0, i1, i2] * tmp
import numpy as np
from numba import jit
@jit
def proj(vx, vy, vz, kx, ky, kz, inv_k_square_nozero):
tmp = (kx * vx + ky * vy + kz * vz) * inv_k_square_nozero
return (vx - kx * tmp,
vy - ky * tmp,
vz - kz * tmp)
@jit
def proj_loop(vx, vy, vz, kx, ky, kz, inv_k_square_nozero):
rx = np.empty_like(vx)
ry = np.empty_like(vx)
rz = np.empty_like(vx)
n0, n1, n2 = kx.shape
for i0 in range(n0):
for i1 in range(n1):
for i2 in range(n2):
tmp = (kx[i0, i1, i2] * vx[i0, i1, i2]
+ ky[i0, i1, i2] * vy[i0, i1, i2]
+ kz[i0, i1, i2] * vz[i0, i1, i2]
) * inv_k_square_nozero[i0, i1, i2]
rx[i0, i1, i2] = vx[i0, i1, i2] - kx[i0, i1, i2] * tmp
ry[i0, i1, i2] = vz[i0, i1, i2] - kx[i0, i1, i2] * tmp
rz[i0, i1, i2] = vy[i0, i1, i2] - kx[i0, i1, i2] * tmp
return rx, ry, rz
@jit
def proj_inplace(vx, vy, vz, kx, ky, kz, inv_k_square_nozero):
tmp = (kx * vx + ky * vy + kz * vz) * inv_k_square_nozero
vx -= kx * tmp
vy -= ky * tmp
vz -= kz * tmp
@jit
def proj_inplace_loop(vx, vy, vz, kx, ky, kz, inv_k_square_nozero):
n0, n1, n2 = kx.shape
for i0 in range(n0):
for i1 in range(n1):
for i2 in range(n2):
tmp = (kx[i0, i1, i2] * vx[i0, i1, i2]
+ ky[i0, i1, i2] * vy[i0, i1, i2]
+ kz[i0, i1, i2] * vz[i0, i1, i2]
) * inv_k_square_nozero[i0, i1, i2]
vx[i0, i1, i2] -= kx[i0, i1, i2] * tmp
vy[i0, i1, i2] -= ky[i0, i1, i2] * tmp
vz[i0, i1, i2] -= kz[i0, i1, i2] * tmp
<sys_info>
<software CC="gcc (Ubuntu 5.4.0-6ubuntu1~16.04.9) 5.4.0 20160609" MPI="mpirun
(Open MPI) 1.10.2" distro="Ubuntu 16.04 Xenial Xerus"
hostname="pierre-KTH" kernel="4.4.0-119-generic" system="Linux"/>
<hardware address_sizes="36 bits physical, 48 bits virtual" arch="x86_64"
bogomips="5388.02" cache_size="4096 KB" cpu_MHz_actual="[3363.187,
3341.144, 3340.617, 3334.289]" cpu_MHz_current="3344.809"
cpu_MHz_max="3400.000" cpu_MHz_min="3400.000" cpu_name="Intel(R)
Core(TM) i7-2620M CPU @ 2.70GHz" nb_cores="2" nb_procs="4"
nb_siblings="4"/>
<python compiler="GCC 4.8.2 20140120 (Red Hat 4.8.2-15)"
implementation="CPython" version="3.6.4">
<fluiddyn installed="True" local_path="/home/pierre/Dev/fluiddyn"
remote_path="https://paugier@bitbucket.org/fluiddyn/fluiddyn"
version="0.2.2"/>
<fluidsim installed="True" local_path="/home/pierre/Dev/fluidsim"
remote_path="calpe = https://paugier@bitbucket.org/calpe/fluidsim"
version="0.1.1"/>
<fluidlab installed="True" local_path="/home/pierre/Dev/fluidlab"
remote_path="https://paugier@bitbucket.org/fluiddyn/fluidlab"
version="0.0.3"/>
<fluidimage installed="True" local_path="/home/pierre/Dev/fluidimage"
remote_path="https://paugier@bitbucket.org/fluiddyn/fluidimage"
version="0.0.2"/>
<fluidfft installed="True" local_path="/home/pierre/Dev/fluidfft"
remote_path="https://paugier@bitbucket.org/fluiddyn/fluidfft"
version="0.2.2"/>
<fluidcoriolis installed="True" local_path="/home/pierre/Dev/fluidcoriolis"
remote_path="https://paugier@bitbucket.org/antcampagne/fluidcoriolis"
version="0.0.0a0"/>
<fluiddevops installed="False" local_path="" remote_path="" version=""/>
<numpy installed="True"
local_path="/home/pierre/opt/miniconda3/lib/python3.6/site-packages/numpy"
remote_path="" version="1.13.3">
<system>
<lapack_opt define_macros="[('SCIPY_MKL_H', None), ('HAVE_CBLAS',
None)]"
include_dirs="['/home/pierre/opt/miniconda3/include']"
libraries="['mkl_rt', 'pthread']"
library_dirs="['/home/pierre/opt/miniconda3/lib']"/>
<blas_opt define_macros="[('SCIPY_MKL_H', None), ('HAVE_CBLAS', None)]"
include_dirs="['/home/pierre/opt/miniconda3/include']"
libraries="['mkl_rt', 'pthread']"
library_dirs="['/home/pierre/opt/miniconda3/lib']"/>
<fftw/>
</system>
<build>
<lapack_opt define_macros="[('HAVE_CBLAS', None)]" language="c"
libraries="['openblas', 'openblas']"
library_dirs="['/home/pierre/opt/miniconda3/lib']"/>
<blas_opt define_macros="[('HAVE_CBLAS', None)]" language="c"
libraries="['openblas', 'openblas']"
library_dirs="['/home/pierre/opt/miniconda3/lib']"/>
</build>
</numpy>
<cython installed="True"
local_path="/home/pierre/opt/miniconda3/lib/python3.6/site-packages"
remote_path="" version="0.26"/>
<mpi4py installed="True"
local_path="/home/pierre/opt/miniconda3/lib/python3.6/site-packages/mpi4py"
remote_path="" version="3.0.0"/>
<pythran installed="True" local_path="/home/pierre/Dev/pythran"
remote_path="git@github.com:fluiddyn/pythran.git"
version="0.8.4post0"/>
<pyfftw installed="True"
local_path="/home/pierre/opt/miniconda3/lib/python3.6/site-packages/pyfftw"
remote_path="" version="0.10.3.dev0+e827cb5"/>
<matplotlib installed="True"
local_path="/home/pierre/opt/miniconda3/lib/python3.6/site-packages/matplotlib"
remote_path="" version="2.1.1"/>
<scipy installed="True"
local_path="/home/pierre/opt/miniconda3/lib/python3.6/site-packages/scipy"
remote_path="" version="0.19.1"/>
<skimage installed="True"
local_path="/home/pierre/opt/miniconda3/lib/python3.6/site-packages/skimage"
remote_path="" version="0.13.0"/>
<h5py installed="True"
local_path="/home/pierre/opt/miniconda3/lib/python3.6/site-packages/h5py"
remote_path="" version="2.7.1">
<config HDF5_version="1.10.1" MPI_enabled="False"
single_writer_multiple_reader_available="True"
virtual_dataset_available="True"/>
</h5py>
</python>
</sys_info>
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment