Skip to content
Snippets Groups Projects
Commit 068de907 authored by Cristian Lalescu's avatar Cristian Lalescu
Browse files

remove ipython notebooks

parent 5051398d
Branches
Tags
No related merge requests found
{
"metadata": {
"name": "",
"signature": "sha256:70c4b065d524fe79b51e885d5e77aec8982c6ae0996146ba76aa979790c214b7"
},
"nbformat": 3,
"nbformat_minor": 0,
"worksheets": [
{
"cells": [
{
"cell_type": "code",
"collapsed": false,
"input": [
"import numpy as np\n",
"import subprocess\n",
"import pyfftw"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 1
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"def generate_data_3D(\n",
" n,\n",
" dtype = np.complex128,\n",
" p = 1.5):\n",
" \"\"\"\n",
" generate something that has the proper shape\n",
" \"\"\"\n",
" assert(n % 2 == 0)\n",
" a = np.zeros((n, n, n/2+1), dtype = dtype)\n",
" a[:] = np.random.randn(*a.shape) + 1j*np.random.randn(*a.shape)\n",
" k, j, i = np.mgrid[-n/2+1:n/2+1, -n/2+1:n/2+1, 0:n/2+1]\n",
" k = (k**2 + j**2 + i**2)**.5\n",
" k = np.roll(k, n//2+1, axis = 0)\n",
" k = np.roll(k, n//2+1, axis = 1)\n",
" a /= k**p\n",
" a[0, :, :] = 0\n",
" a[:, 0, :] = 0\n",
" a[:, :, 0] = 0\n",
" ii = np.where(k == 0)\n",
" a[ii] = 0\n",
" ii = np.where(k > n/3)\n",
" a[ii] = 0\n",
" return a\n",
"\n",
"n = 32\n",
"\n",
"Kdata0 = generate_data_3D(n, p = 2).astype(np.complex64)\n",
"Kdata0.tofile(\"Kdata0\")"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stderr",
"text": [
"-c:15: RuntimeWarning: divide by zero encountered in divide\n"
]
}
],
"prompt_number": 2
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"def compute_cpp_data(\n",
" branch = None):\n",
" if not (type(branch) == type(None)):\n",
" subprocess.call(['git', 'checkout', branch])\n",
" if subprocess.call(['make', 'main_fluid_solver.elf']) == 0:\n",
" subprocess.call(['time',\n",
" 'mpirun',\n",
" '-np',\n",
" '4',\n",
" './main_fluid_solver.elf'])\n",
" return np.fromfile('Kdata1')\n",
" else:\n",
" print ('compilation error')\n",
" return None\n",
"\n",
"Kdata1 = compute_cpp_data()"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 4
},
{
"cell_type": "code",
"collapsed": true,
"input": [
"distance = np.max(np.abs(Rdata_py - Rdata), axis = (1, 2, 3, 4))\n",
"print(np.max(distance))\n",
"if np.max(distance) > 1e-5:\n",
" ax = plt.figure(figsize=(6,2)).add_subplot(111)\n",
" ax.plot(distance)\n",
" i0 = np.random.randint(8)\n",
" i1 = np.random.randint(8)\n",
" i2 = np.random.randint(8)\n",
" z = cm.grid3D_to_zindex(np.array([i0, i1, i2]))\n",
" #z = 0\n",
" print(cm.zindex_to_grid3D(z))\n",
" s = np.max(np.abs(Rdata_py[None, z, :, :, :, 1] - Rdata[..., 1]),\n",
" axis = (1, 2, 3))\n",
" z1 = np.argmin(s)\n",
" print(z, z1, s[z1])\n",
" #print(Rdata[z1] - Rdata_py[z1])\n",
" ta0 = Rdata_py.ravel()\n",
" ta1 = Rdata.ravel()\n",
" print (Rdata_py[254:259, 7, 4, 3, 1])\n",
" print (Rdata[254:259, 7, 4, 3, 1])\n",
" print (ta0[ta0.shape[0]/2-1:ta0.shape[0]/2+7])\n",
" print (ta1[ta1.shape[0]/2-1:ta1.shape[0]/2+7])\n",
"else:\n",
" print('distance is small')\n",
"print(np.max(np.abs(Rdata)))"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"1.90735e-06\n",
"distance is small\n",
"15.316\n"
]
}
],
"prompt_number": 7
},
{
"cell_type": "code",
"collapsed": false,
"input": [],
"language": "python",
"metadata": {},
"outputs": []
}
],
"metadata": {}
}
]
}
\ No newline at end of file
{
"metadata": {
"name": "",
"signature": "sha256:46290caaf0d19fceaf2c8987f305c20f9083f9d3d3d0ee8ef3987689c3fccbf0"
},
"nbformat": 3,
"nbformat_minor": 0,
"worksheets": [
{
"cells": [
{
"cell_type": "code",
"collapsed": false,
"input": [
"import numpy as np\n",
"import subprocess\n",
"%matplotlib inline\n",
"import matplotlib.pyplot as plt\n",
"import pyfftw"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 1
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"def generate_data_3D(\n",
" n,\n",
" dtype = np.complex128,\n",
" p = 1.5):\n",
" \"\"\"\n",
" generate something that has the proper shape\n",
" \"\"\"\n",
" assert(n % 2 == 0)\n",
" a = np.zeros((n, n, n/2+1), dtype = dtype)\n",
" a[:] = np.random.randn(*a.shape) + 1j*np.random.randn(*a.shape)\n",
" k, j, i = np.mgrid[-n/2+1:n/2+1, -n/2+1:n/2+1, 0:n/2+1]\n",
" k = (k**2 + j**2 + i**2)**.5\n",
" k = np.roll(k, n//2+1, axis = 0)\n",
" k = np.roll(k, n//2+1, axis = 1)\n",
" a /= k**p\n",
" a[0, :, :] = 0\n",
" a[:, 0, :] = 0\n",
" a[:, :, 0] = 0\n",
" ii = np.where(k == 0)\n",
" a[ii] = 0\n",
" ii = np.where(k > n/3)\n",
" a[ii] = 0\n",
" return a\n",
"\n",
"n = 31*4\n",
"N = 256\n",
"\n",
"Kdata0 = generate_data_3D(n, p = 2).astype(np.complex64)\n",
"Kdata1 = generate_data_3D(n, p = 2).astype(np.complex64)\n",
"Kdata2 = generate_data_3D(n, p = 2).astype(np.complex64)\n",
"Kdata0.T.copy().astype('>c8').tofile(\"Kdata0\")\n",
"Kdata1.T.copy().astype('>c8').tofile(\"Kdata1\")\n",
"Kdata2.T.copy().astype('>c8').tofile(\"Kdata2\")"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stderr",
"text": [
"-c:15: RuntimeWarning: divide by zero encountered in divide\n"
]
}
],
"prompt_number": 2
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"def padd_with_zeros(\n",
" a,\n",
" n,\n",
" odtype = None):\n",
" if (type(odtype) == type(None)):\n",
" odtype = a.dtype\n",
" assert(a.shape[0] <= n)\n",
" b = np.zeros((n, n, n/2 + 1), dtype = odtype)\n",
" m = a.shape[0]\n",
" b[ :m/2, :m/2, :m/2+1] = a[ :m/2, :m/2, :m/2+1]\n",
" b[ :m/2, n-m/2: , :m/2+1] = a[ :m/2, m-m/2: , :m/2+1]\n",
" b[n-m/2: , :m/2, :m/2+1] = a[m-m/2: , :m/2, :m/2+1]\n",
" b[n-m/2: , n-m/2: , :m/2+1] = a[m-m/2: , m-m/2: , :m/2+1]\n",
" return b\n",
"\n",
"def transform_py(bla):\n",
" b = padd_with_zeros(bla, N)\n",
" c = np.zeros((N, N, N), np.float32)\n",
" t = pyfftw.FFTW(\n",
" b, c,\n",
" axes = (0, 1, 2),\n",
" direction = 'FFTW_BACKWARD',\n",
" flags = ('FFTW_ESTIMATE',),\n",
" threads = 2)\n",
" t.execute()\n",
" return c\n",
"\n",
"import chichi_misc as cm\n",
"\n",
"def array_to_8cubes(\n",
" a,\n",
" odtype = None):\n",
" assert(len(a.shape) >= 3)\n",
" assert((a.shape[0] % 8 == 0) and\n",
" (a.shape[1] % 8 == 0) and\n",
" (a.shape[2] % 8 == 0))\n",
" if type(odtype) == type(None):\n",
" odtype = a.dtype\n",
" c = np.zeros(\n",
" ((((a.shape[0] // 8)*(a.shape[1] // 8)*(a.shape[2] // 8)),) +\n",
" (8, 8, 8) +\n",
" a.shape[3:]),\n",
" dtype = odtype)\n",
" zi = np.zeros( c.shape[0], dtype = np.int)\n",
" ri = np.zeros((c.shape[0], 3, 2), dtype = np.int)\n",
" ii = 0\n",
" for k in range(a.shape[0]//8):\n",
" for j in range(a.shape[1]//8):\n",
" for i in range(a.shape[2]//8):\n",
" tindex = np.array([k, j, i])\n",
" zi[ii] = cm.grid3D_to_zindex(tindex)\n",
" ri[ii, 0] = np.array([8*tindex[0], 8*(tindex[0]+1)])\n",
" ri[ii, 1] = np.array([8*tindex[1], 8*(tindex[1]+1)])\n",
" ri[ii, 2] = np.array([8*tindex[2], 8*(tindex[2]+1)])\n",
" ii += 1\n",
" for ii in range(zi.shape[0]):\n",
" c[zi[ii]] = a[ri[ii, 0, 0]:ri[ii, 0, 1],\n",
" ri[ii, 1, 0]:ri[ii, 1, 1],\n",
" ri[ii, 2, 0]:ri[ii, 2, 1]]\n",
" return c\n",
"\n",
"d0 = transform_py(Kdata0)\n",
"d1 = transform_py(Kdata1)\n",
"d2 = transform_py(Kdata2)\n",
"\n",
"Rdata_py_tmp = np.array([d0, d1, d2]).transpose((1, 2, 3, 0))\n",
"\n",
"Rdata_py = array_to_8cubes(Rdata_py_tmp)\n",
"\n",
"# i0 = np.random.randint(16)\n",
"# i1 = np.random.randint(16)\n",
"# i2 = np.random.randint(16)\n",
"# z = cm.grid3D_to_zindex(np.array([i0, i1, i2]))"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 3
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"def compute_cpp_data(\n",
" branch = None,\n",
" nfiles = 16):\n",
" if not (type(branch) == type(None)):\n",
" subprocess.call(['git', 'checkout', branch])\n",
" if subprocess.call(['make', 'full.elf']) == 0:\n",
" subprocess.call([#'valgrind',\n",
" #'--tool=callgrind',\n",
" #'--callgrind-out-file=tmp.txt',\n",
" 'time',\n",
" 'mpirun',\n",
" '-np',\n",
" '4',\n",
" './full.elf',\n",
" '{0}'.format(n),\n",
" '{0}'.format(N),\n",
" '{0}'.format(nfiles),\n",
" '3'])\n",
" else:\n",
" print ('compilation error')\n",
" return None\n",
" \n",
"def get_cpp_data(\n",
" branch = None,\n",
" run = True,\n",
" nfiles = 16):\n",
" if run:\n",
" for nf in range(nfiles):\n",
" subprocess.call(\n",
" ['rm',\n",
" 'Rdata_z{0:0>7x}'.format(nf*Rdata_py.shape[0]//nfiles)])\n",
" compute_cpp_data(branch, nfiles = nfiles)\n",
" Rdata = []\n",
" for nf in range(nfiles):\n",
" Rdata.append(np.fromfile(\n",
" 'Rdata_z{0:0>7x}'.format(nf*Rdata_py.shape[0]//nfiles),\n",
" dtype = np.float32).reshape(-1, 8, 8, 8, 3))\n",
" return np.concatenate(Rdata)\n",
"\n",
"#Rdata = get_cpp_data(branch = 'develop')\n",
"# develop says 30 secs, inplace fft is 28 secs\n",
"#Rdata = get_cpp_data(branch = 'feature-inplace_fft')\n",
"Rdata = get_cpp_data(run = True, nfiles = 16)"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 6
},
{
"cell_type": "code",
"collapsed": true,
"input": [
"distance = np.max(np.abs(Rdata_py - Rdata), axis = (1, 2, 3, 4))\n",
"print(np.max(distance))\n",
"if np.max(distance) > 1e-5:\n",
" ax = plt.figure(figsize=(6,2)).add_subplot(111)\n",
" ax.plot(distance)\n",
" i0 = np.random.randint(8)\n",
" i1 = np.random.randint(8)\n",
" i2 = np.random.randint(8)\n",
" z = cm.grid3D_to_zindex(np.array([i0, i1, i2]))\n",
" #z = 0\n",
" print(cm.zindex_to_grid3D(z))\n",
" s = np.max(np.abs(Rdata_py[None, z, :, :, :, 1] - Rdata[..., 1]),\n",
" axis = (1, 2, 3))\n",
" z1 = np.argmin(s)\n",
" print(z, z1, s[z1])\n",
" #print(Rdata[z1] - Rdata_py[z1])\n",
" ta0 = Rdata_py.ravel()\n",
" ta1 = Rdata.ravel()\n",
" print (Rdata_py[254:259, 7, 4, 3, 1])\n",
" print (Rdata[254:259, 7, 4, 3, 1])\n",
" print (ta0[ta0.shape[0]/2-1:ta0.shape[0]/2+7])\n",
" print (ta1[ta1.shape[0]/2-1:ta1.shape[0]/2+7])\n",
"else:\n",
" print('distance is small')\n",
"print(np.max(np.abs(Rdata)))"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"1.90735e-06\n",
"distance is small\n",
"15.316\n"
]
}
],
"prompt_number": 7
},
{
"cell_type": "code",
"collapsed": false,
"input": [],
"language": "python",
"metadata": {},
"outputs": []
}
],
"metadata": {}
}
]
}
\ No newline at end of file
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment