smyalygames/FiniteVolumeGPU
1
0
Fork 0
mirror of https://github.com/smyalygames/FiniteVolumeGPU.git synced 2026-01-13 23:28:41 +01:00
Code Issues Actions Packages Projects Releases Wiki Activity

refactor: replace string format with f strings and follow PEP8 formatting

Browse Source
This commit is contained in:
Anthony Berg
2025-06-24 18:27:35 +02:00
parent 1786ca979f
commit 2726a4dfce
13 changed files with 916 additions and 1453 deletions
Download Patch File Download Diff File Expand all files Collapse all files

75
Autotuning.ipynb
View File

@@ -37,7 +37,7 @@
" from StringIO import StringIO\n",
"except ImportError:\n",
" from io import StringIO\n",
" \n",
"\n",
"#Set large figure sizes\n",
"#Note, this prevents nice figures for articles...\n",
"rc('figure', figsize=(16.0, 12.0))\n",
@@ -51,7 +51,7 @@
"outputs": [],
"source": [
"from GPUSimulators import LxF, FORCE, HLL, HLL2, KP07, KP07_dimsplit, WAF, Autotuner\n",
"from GPUSimulators.common import common"
"from GPUSimulators.common import Timer"
]
},
{
@@ -93,15 +93,15 @@
"source": [
"with np.load(autotuner.filename) as data:\n",
" simulators = data['simulators']\n",
" \n",
"\n",
" for simulator in simulators:\n",
" megacells = data[simulator + \"_megacells\"]\n",
" block_widths = data[simulator + '_block_widths']\n",
" block_heights = data[simulator + '_block_heights']\n",
" arguments = data[simulator + '_arguments']\n",
" \n",
"\n",
" plt.figure()\n",
" plt.imshow(megacells, origin='lower')#, vmax=maximum, vmin=minimum)\n",
" plt.imshow(megacells, origin='lower') #, vmax=maximum, vmin=minimum)\n",
" plt.xlabel('Block Width')\n",
" plt.xticks(range(len(block_widths)), block_widths)\n",
" plt.ylabel('Block Height')\n",
@@ -135,7 +135,8 @@
"simulators = [LxF.LxF, FORCE.FORCE, HLL.HLL, HLL2.HLL2, KP07.KP07, KP07_dimsplit.KP07_dimsplit, WAF.WAF]\n",
"peak_performance = [autotuner.get_peak_performance(simulator) for simulator in simulators]\n",
"megacells = [performance['megacells'] for performance in peak_performance]\n",
"xlabels = [\"{:s}\\n[{:d}x{:d}]\".format(simulators[i].__name__, performance['block_width'], performance['block_height']) for i, performance in enumerate(peak_performance)]\n",
"xlabels = [f\"{simulators[i].__name__}\\n[{performance['block_width']}x{performance['block_height']}]\" for i, performance\n",
" in enumerate(peak_performance)]\n",
"\n",
"plt.figure()\n",
"plt.bar(range(len(simulators)), megacells)\n",
@@ -208,51 +209,51 @@
" dx = width / float(nx)\n",
" dy = height / float(ny)\n",
"\n",
" x_center = dx*nx/2.0\n",
" y_center = dy*ny/2.0\n",
" x_center = dx * nx / 2.0\n",
" y_center = dy * ny / 2.0\n",
"\n",
" #Create a gaussian \"dam break\" that will not form shocks\n",
" size = width / 5.0\n",
" dt = 10**10\n",
" dt = 10 ** 10\n",
"\n",
" h = np.zeros((ny, nx), dtype=np.float32); \n",
" hu = np.zeros((ny, nx), dtype=np.float32);\n",
" hv = np.zeros((ny, nx), dtype=np.float32);\n",
" h = np.zeros((ny, nx), dtype=np.float32)\n",
" hu = np.zeros((ny, nx), dtype=np.float32)\n",
" hv = np.zeros((ny, nx), dtype=np.float32)\n",
"\n",
" extent = 1.0/np.sqrt(2.0)\n",
" x = (dx*(np.arange(0, nx, dtype=np.float32)+0.5) - x_center) / size\n",
" y = (dy*(np.arange(0, ny, dtype=np.float32)+0.5) - y_center) / size\n",
" extent = 1.0 / np.sqrt(2.0)\n",
" x = (dx * (np.arange(0, nx, dtype=np.float32) + 0.5) - x_center) / size\n",
" y = (dy * (np.arange(0, ny, dtype=np.float32) + 0.5) - y_center) / size\n",
" xv, yv = np.meshgrid(x, y, sparse=False, indexing='xy')\n",
" r = np.minimum(1.0, np.sqrt(xv**2 + yv**2))\n",
" r = np.minimum(1.0, np.sqrt(xv ** 2 + yv ** 2))\n",
" xv = None\n",
" yv = None\n",
" gc.collect()\n",
"\n",
" #Generate highres\n",
" cos = np.cos(np.pi*r)\n",
" h = 0.5 + 0.1*0.5*(1.0 + cos)\n",
" hu = 0.1*0.5*(1.0 + cos)\n",
" cos = np.cos(np.pi * r)\n",
" h = 0.5 + 0.1 * 0.5 * (1.0 + cos)\n",
" hu = 0.1 * 0.5 * (1.0 + cos)\n",
" hv = hu.copy()\n",
"\n",
" scale = 0.7\n",
" max_h_estimate = 0.6\n",
" max_u_estimate = 0.1*np.sqrt(2.0)\n",
" dx = width/nx\n",
" dy = height/ny\n",
" dt = scale * min(dx, dy) / (max_u_estimate + np.sqrt(g*max_h_estimate))\n",
" max_u_estimate = 0.1 * np.sqrt(2.0)\n",
" dx = width / nx\n",
" dy = height / ny\n",
" dt = scale * min(dx, dy) / (max_u_estimate + np.sqrt(g * max_h_estimate))\n",
"\n",
" return h, hu, hv, dx, dy, dt\n",
"\n",
"\n",
"with Common.Timer(\"gen_data\", log_level=logging.INFO) as t:\n",
"with Timer(\"gen_data\", log_level=logging.INFO) as t:\n",
" h0, hu0, hv0, dx, dy, dt = gen_test_data(4096, 4096, 9.81)\n",
" \n",
"\n",
"plt.figure(figsize=(12, 8))\n",
"plt.subplot(1,3,1)\n",
"plt.subplot(1, 3, 1)\n",
"plt.imshow(h0, origin='lower', interpolation=\"none\")\n",
"plt.subplot(1,3,2)\n",
"plt.subplot(1, 3, 2)\n",
"plt.imshow(hu0, origin='lower', interpolation=\"none\")\n",
"plt.subplot(1,3,3)\n",
"plt.subplot(1, 3, 3)\n",
"plt.imshow(hv0, origin='lower', interpolation=\"none\")"
]
},
@@ -490,7 +491,7 @@
"run_simulation = True\n",
"sizes = list(range(64, 512, 64)) + list(range(512, 2048, 128)) + list(range(2048, 4096, 256)) + [4096]\n",
"simulators = [LxF.LxF, FORCE.FORCE, HLL.HLL, HLL2.HLL2, KP07.KP07, KP07_dimsplit.KP07_dimsplit, WAF.WAF]\n",
"if (run_simulation):\n",
"if run_simulation:\n",
" megacells = {}\n",
" for simulator in simulators:\n",
" print(simulator.__name__)\n",
@@ -508,9 +509,9 @@
" 'context': my_context,\n",
" 'h0': h0, 'hu0': hu0, 'hv0': hv0,\n",
" 'nx': nx, 'ny': ny,\n",
" 'dx': dx, 'dy': dy, 'dt': 0.9*dt,\n",
" 'dx': dx, 'dy': dy, 'dt': 0.9 * dt,\n",
" 'g': g\n",
" } \n",
" }\n",
"\n",
" sim = simulator(**arguments)\n",
"\n",
@@ -535,8 +536,8 @@
" gc.collect()\n",
"\n",
" #Compute megacells\n",
" gpu_elapsed = end.time_since(start)*1.0e-3\n",
" megacells[simulator.__name__][k] = (nx*ny*timesteps / (1000*1000)) / gpu_elapsed\n",
" gpu_elapsed = end.time_since(start) * 1.0e-3\n",
" megacells[simulator.__name__][k] = (nx * ny * timesteps / (1000 * 1000)) / gpu_elapsed\n",
" print(\"[{:d}x{:d}] => {:.1f} ({:2f})\".format(nx, ny, megacells[simulator.__name__][k], gpu_elapsed))"
]
},
@@ -555,7 +556,7 @@
],
"source": [
"datafilename = \"megacells.npz\"\n",
"if (not os.path.isfile(datafilename) and \"megacells\" in globals()):\n",
"if not os.path.isfile(datafilename) and \"megacells\" in globals():\n",
" print(\"Saving data to file\")\n",
" np.savez_compressed(datafilename, megacells=megacells)\n",
"else:\n",
@@ -633,7 +634,7 @@
"for simulator in simulators:\n",
" old = megacells[simulator.__name__][0:-1]\n",
" new = megacells[simulator.__name__][1:]\n",
" change = 100*(new-old)/old\n",
" change = 100 * (new - old) / old\n",
" plt.plot(sizes[0:-1], change, '.-', label=simulator.__name__)\n",
"plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)\n",
"plt.title(\"Relative performance change\")\n",
@@ -654,9 +655,7 @@
]
}
],
"source": [
"print(type(None) == None)"
]
"source": "print(type(None) is None)"
},
{
"cell_type": "code",

271
ConvergenceShock1D.ipynb
View File

@@ -57,7 +57,7 @@
"\n",
"#Finally, import our simulator\n",
"from GPUSimulators import LxF, FORCE, HLL, HLL2, KP07, KP07_dimsplit, WAF\n",
"from GPUSimulators.common import common"
"from GPUSimulators.common import Timer"
]
},
{
@@ -95,13 +95,14 @@
"outputs": [],
"source": [
"#Misc plotting setup\n",
"def setBwStyles(ax):\n",
"def set_bw_styles(ax):\n",
" from cycler import cycler\n",
"\n",
" ax.set_prop_cycle( cycler('marker', ['.', 'x', 4, '+', '*', '1', 5]) +\n",
" cycler('linestyle', ['-.', '--', ':', '-.', '--', ':', '-.']) +\n",
" #cycler('markersize', [5, 5, 5, 5, 5, 5]) +\n",
" cycler('color', ['k', 'k', 'k', 'k', 'k', 'k', 'k']) ) \n",
" ax.set_prop_cycle(cycler('marker', ['.', 'x', 4, '+', '*', '1', 5]) +\n",
" cycler('linestyle', ['-.', '--', ':', '-.', '--', ':', '-.']) +\n",
" #cycler('markersize', [5, 5, 5, 5, 5, 5]) +\n",
" cycler('color', ['k', 'k', 'k', 'k', 'k', 'k', 'k']))\n",
"\n",
"\n",
"#Set large figure sizes\n",
"plt.rcParams['figure.figsize'] = [12, 8]\n",
@@ -120,26 +121,26 @@
"outputs": [],
"source": [
"def save_figure(fig, stem):\n",
" if (not os.path.isdir(\"figures\")):\n",
" if not os.path.isdir(\"figures\"):\n",
" os.mkdir(\"figures\")\n",
" \n",
"\n",
" fig_filename = os.path.join(\"figures\", \"convergence_shock1d_\" + stem + \".pdf\")\n",
" \n",
"\n",
" metadata = {\n",
" 'created': time.strftime(\"%Y_%m_%d-%H_%M_%S\"),\n",
" 'hostname': socket.gethostname()\n",
" }\n",
" \n",
"\n",
" legend = fig.gca().legend_\n",
" if (legend != None):\n",
" fig.savefig(fig_filename, dpi=300,format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None, \n",
" metadata=metadata, \n",
" bbox_extra_artists=(legend, ), bbox_inches='tight')\n",
" if legend != None:\n",
" fig.savefig(fig_filename, dpi=300, format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None,\n",
" metadata=metadata,\n",
" bbox_extra_artists=(legend,), bbox_inches='tight')\n",
" else:\n",
" fig.savefig(fig_filename, dpi=300,format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None, \n",
" metadata=metadata)"
" fig.savefig(fig_filename, dpi=300, format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None,\n",
" metadata=metadata)"
]
},
{
@@ -152,24 +153,24 @@
" csv_filename = os.path.abspath(os.path.join(\"reference\", \"swashes_1_nx=\" + str(nx) + \".csv\"))\n",
"\n",
" #If we do not have the data, generate it \n",
" if (not os.path.isfile(csv_filename)):\n",
" if not os.path.isfile(csv_filename):\n",
" print(\"Generating new reference!\")\n",
" swashes_path = r'C:\\Users\\anbro\\Documents\\programs\\SWASHES-1.03.00_win\\bin\\swashes_win.exe'\n",
"\n",
" swashes_args = [\\\n",
" '1', # 1D problems \\\n",
" '3', # Dam breaks \\\n",
" '1', # Domain 1 \\\n",
" '1', # Wet domain no friction\n",
" str(nx) #Number of cells X\n",
" ]\n",
" swashes_args = [\n",
" '1', # 1D problems \\\n",
" '3', # Dam breaks \\\n",
" '1', # Domain 1 \\\n",
" '1', # Wet domain no friction\n",
" str(nx) #Number of cells X\n",
" ]\n",
"\n",
" with open(csv_filename, 'w') as csv_file:\n",
" p = subprocess.check_call([swashes_path] + swashes_args, stdout=csv_file)\n",
"\n",
" reference = np.genfromtxt(csv_filename, comments='#', delimiter='\\t', skip_header=0, usecols=(0, 1, 2))\n",
" x, h, u = reference[:, 0], reference[:, 1], reference[:, 2]\n",
" return x, h, h*u"
" return x, h, h * u"
]
},
{
@@ -179,17 +180,17 @@
"outputs": [],
"source": [
"def downsample(highres_solution, factor):\n",
" if (len(highres_solution.shape) == 1):\n",
" if len(highres_solution.shape) == 1:\n",
" highres_solution = highres_solution.reshape((1, highres_solution.size))\n",
"\n",
" assert(highres_solution.shape[1] % factor == 0)\n",
" ny = highres_solution.shape[0] \n",
" assert (highres_solution.shape[1] % factor == 0)\n",
" ny = highres_solution.shape[0]\n",
" nx = highres_solution.shape[1] / factor\n",
"\n",
" if (factor == 1):\n",
" if factor == 1:\n",
" return highres_solution\n",
" else:\n",
" \n",
"\n",
" return highres_solution.reshape([int(ny), int(1), int(nx), int(factor)]).mean(3).mean(1)"
]
},
@@ -200,10 +201,10 @@
"outputs": [],
"source": [
"def wall_boundary_conditions(data, num_ghost_cells):\n",
" data[0:num_ghost_cells,:] = data[2*num_ghost_cells-1:num_ghost_cells-1:-1,:]\n",
" data[-num_ghost_cells:,:] = data[-num_ghost_cells-1:-2*num_ghost_cells-1:-1,:]\n",
" data[:,0:num_ghost_cells] = data[:,2*num_ghost_cells-1:num_ghost_cells-1:-1]\n",
" data[:,-num_ghost_cells:] = data[:,-num_ghost_cells-1:-2*num_ghost_cells-1:-1]\n",
" data[0:num_ghost_cells, :] = data[2 * num_ghost_cells - 1:num_ghost_cells - 1:-1, :]\n",
" data[-num_ghost_cells:, :] = data[-num_ghost_cells - 1:-2 * num_ghost_cells - 1:-1, :]\n",
" data[:, 0:num_ghost_cells] = data[:, 2 * num_ghost_cells - 1:num_ghost_cells - 1:-1]\n",
" data[:, -num_ghost_cells:] = data[:, -num_ghost_cells - 1:-2 * num_ghost_cells - 1:-1]\n",
" return data"
]
},
@@ -232,49 +233,49 @@
],
"source": [
"def gen_test_data(nx, ref_nx, g, num_ghost_cells):\n",
" assert(num_ghost_cells >= 1)\n",
" assert (num_ghost_cells >= 1)\n",
" ref_ny = num_ghost_cells\n",
" ny = ref_ny\n",
" \n",
"\n",
" dx = 10.0 / float(nx)\n",
" dy = 1.0\n",
"\n",
" h_highres = np.zeros((ref_ny, ref_nx))\n",
" h = np.zeros((ny+2*num_ghost_cells, nx+2*num_ghost_cells), dtype=np.float32); \n",
" hu = np.zeros((ny+2*num_ghost_cells, nx+2*num_ghost_cells), dtype=np.float32);\n",
" hv = np.zeros((ny+2*num_ghost_cells, nx+2*num_ghost_cells), dtype=np.float32);\n",
" h = np.zeros((ny + 2 * num_ghost_cells, nx + 2 * num_ghost_cells), dtype=np.float32)\n",
" hu = np.zeros((ny + 2 * num_ghost_cells, nx + 2 * num_ghost_cells), dtype=np.float32)\n",
" hv = np.zeros((ny + 2 * num_ghost_cells, nx + 2 * num_ghost_cells), dtype=np.float32)\n",
"\n",
" #Create dam break\n",
" x = np.linspace(0, 10, ref_nx)\n",
" \n",
" h_highres = np.where(x < 5, 0.005, 0.001); \n",
" h[num_ghost_cells:-num_ghost_cells, num_ghost_cells:-num_ghost_cells] = downsample(h_highres, ref_nx/nx)\n",
" h = wall_boundary_conditions(h, num_ghost_cells);\n",
" \n",
"\n",
" h_highres = np.where(x < 5, 0.005, 0.001)\n",
" h[num_ghost_cells:-num_ghost_cells, num_ghost_cells:-num_ghost_cells] = downsample(h_highres, ref_nx / nx)\n",
" h = wall_boundary_conditions(h, num_ghost_cells)\n",
"\n",
" max_h_estimate = 0.005\n",
" max_u_estimate = 0.3\n",
" dt = min(dx, dy) / (max_u_estimate + np.sqrt(g*max_h_estimate))\n",
" \n",
" dt = min(dx, dy) / (max_u_estimate + np.sqrt(g * max_h_estimate))\n",
"\n",
" return h, hu, hv, dx, dy, dt\n",
"\n",
"\n",
"h, hu, hv, dx, dy, dt = gen_test_data(nx=100, ref_nx=1000, g=9.81, num_ghost_cells=2)\n",
"\n",
"plt.figure()\n",
"ax=plt.subplot(1,3,1)\n",
"im=plt.imshow(h, extent=[0, 100, 0, 100])\n",
"ax = plt.subplot(1, 3, 1)\n",
"im = plt.imshow(h, extent=[0, 100, 0, 100])\n",
"divider = make_axes_locatable(ax)\n",
"cax = divider.append_axes(\"right\", size=\"5%\", pad=0.05)\n",
"plt.colorbar(im, cax=cax)\n",
"\n",
"ax=plt.subplot(1,3,2)\n",
"im=plt.imshow(hu, extent=[0, 100, 0, 100])\n",
"ax = plt.subplot(1, 3, 2)\n",
"im = plt.imshow(hu, extent=[0, 100, 0, 100])\n",
"divider = make_axes_locatable(ax)\n",
"cax = divider.append_axes(\"right\", size=\"5%\", pad=0.05)\n",
"plt.colorbar(im, cax=cax)\n",
"\n",
"ax=plt.subplot(1,3,3)\n",
"im=plt.imshow(hv, extent=[0, 100, 0, 100])\n",
"ax = plt.subplot(1, 3, 3)\n",
"im = plt.imshow(hv, extent=[0, 100, 0, 100])\n",
"divider = make_axes_locatable(ax)\n",
"cax = divider.append_axes(\"right\", size=\"5%\", pad=0.05)\n",
"plt.colorbar(im, cax=cax)\n",
@@ -305,10 +306,12 @@
" minval = np.amin(variable)\n",
" assert np.isnan(maxval) == False and np.isnan(minval) == False, variable_name + \" contains NaN values\"\n",
"\n",
"def run_benchmark(datafilename, simulator, courant_number, nx, ref_nx, ghost_cells, g=9.81, force_rerun=False, transpose=False):\n",
"\n",
"def run_benchmark(datafilename, simulator, courant_number, nx, ref_nx, ghost_cells, g=9.81, force_rerun=False,\n",
" transpose=False):\n",
" datafilename = gen_filename(simulator, nx)\n",
"\n",
" if (datafilename and os.path.isfile(datafilename) and force_rerun == False):\n",
" if datafilename and os.path.isfile(datafilename) and force_rerun == False:\n",
" print(\"WARNING: Existing simulation exists, skipping simulation\")\n",
" return [0, 0, 0]\n",
" else:\n",
@@ -317,23 +320,23 @@
" dt = dt * courant_number\n",
"\n",
" #Initialize simulator\n",
" with Common.Timer(simulator.__name__ + \"_\" + str(nx)) as timer:\n",
" if (transpose):\n",
" with Timer(simulator.__name__ + \"_\" + str(nx)) as timer:\n",
" if transpose:\n",
" h0 = np.ascontiguousarray(np.transpose(h0))\n",
" hu0, hv0 = np.ascontiguousarray(np.transpose(hv0)), np.ascontiguousarray(np.transpose(hu0))\n",
" dx, dy = dy, dx\n",
" nx, ny = ny, nx\n",
" \n",
" sim = simulator(cuda_context, \\\n",
" h0, hu0, hv0, \\\n",
" nx, ny, \\\n",
" dx, dy, dt, \\\n",
"\n",
" sim = simulator(cuda_context,\n",
" h0, hu0, hv0,\n",
" nx, ny,\n",
" dx, dy, dt,\n",
" g)\n",
"\n",
" t, nt = sim.simulate(6.0)\n",
" h, hu, hv = sim.download()\n",
"\n",
" if (transpose):\n",
" if transpose:\n",
" h = np.transpose(h)\n",
" hu, hv = np.ascontiguousarray(np.transpose(hv)), np.ascontiguousarray(np.transpose(hu))\n",
"\n",
@@ -341,14 +344,14 @@
" sanity_check(hu, \"hu\")\n",
" sanity_check(hv, \"hv\")\n",
"\n",
" if (datafilename):\n",
" if datafilename:\n",
" dirname = os.path.dirname(datafilename)\n",
" if (dirname and not os.path.isdir(dirname)):\n",
" if dirname and not os.path.isdir(dirname):\n",
" os.makedirs(dirname)\n",
" np.savez_compressed(datafilename, h=h, hu=hu, hv=hv)\n",
" \n",
" gc.collect() # Force run garbage collection to free up memory\n",
" \n",
"\n",
" gc.collect() # Force run garbage collection to free up memory\n",
"\n",
" return [timer.secs, nt, dt]"
]
},
@@ -463,8 +466,8 @@
"simulator_ghost_cells = [1, 1, 1, 2, 2, 2, 2]\n",
"simulator_reference = [LxF.LxF, FORCE.FORCE, HLL.HLL, HLL2.HLL2, KP07.KP07, KP07_dimsplit.KP07_dimsplit, WAF.WAF]\n",
"simulator_courant_numbers = [0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95]\n",
"force_rerun=True\n",
"transpose=True\n",
"force_rerun = True\n",
"transpose = True\n",
"\n",
"#simulators = [HLL2.HLL2, WAF.WAF]\n",
"#simulator_ghost_cells = [2, 2]\n",
@@ -473,40 +476,40 @@
"#force_rerun=True\n",
"\n",
"#Warmup to compile kernels\n",
"with Common.Timer(\"warmup\") as timer:\n",
"with Timer(\"warmup\") as timer:\n",
" for i in range(len(simulators)):\n",
" run_benchmark(\"warmup\", \n",
" simulators[i], \n",
" simulator_courant_numbers[i], \n",
" 16, 16, \n",
" simulator_ghost_cells[i], \n",
" run_benchmark(\"warmup\",\n",
" simulators[i],\n",
" simulator_courant_numbers[i],\n",
" 16, 16,\n",
" simulator_ghost_cells[i],\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
"print(\"Completed warmup in {:} s\".format(timer.secs))\n",
" \n",
"with Common.Timer(\"simulation\") as timer:\n",
" for i in range(len(simulators)): \n",
"print(f\"Completed warmup in {timer.secs} s\")\n",
"\n",
"with Timer(\"simulation\") as timer:\n",
" for i in range(len(simulators)):\n",
" for j, nx in enumerate(domain_sizes):\n",
" datafilename = gen_filename(simulators[i], nx)\n",
" run_benchmark(datafilename, \\\n",
" simulators[i], \\\n",
" simulator_courant_numbers[i], \\\n",
" nx, reference_nx, \\\n",
" simulator_ghost_cells[i], \\\n",
" run_benchmark(datafilename,\n",
" simulators[i],\n",
" simulator_courant_numbers[i],\n",
" nx, reference_nx,\n",
" simulator_ghost_cells[i],\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
" gc.collect()\n",
" \n",
"\n",
" # Run reference with a low CFL-number\n",
" datafilename = gen_filename(simulators[i], reference_nx)\n",
" run_benchmark(datafilename, \n",
" simulators[i], \n",
" simulator_courant_numbers[i] * 0.5, \n",
" reference_nx, reference_nx, \n",
" simulator_ghost_cells[i], \n",
" run_benchmark(datafilename,\n",
" simulators[i],\n",
" simulator_courant_numbers[i] * 0.5,\n",
" reference_nx, reference_nx,\n",
" simulator_ghost_cells[i],\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
"print(\"Completed simulation in {:} s\".format(timer.secs))"
"print(f\"Completed simulation in {timer.secs} s\")"
]
},
{
@@ -541,7 +544,7 @@
"error = np.zeros((len(simulators), len(domain_sizes)))\n",
"\n",
"for i, simulator in enumerate(simulators):\n",
" if (use_analytic):\n",
" if use_analytic:\n",
" #Read the reference solution\n",
" _, h_ref, _ = gen_reference(reference_nx)\n",
" else:\n",
@@ -549,37 +552,36 @@
" datafilename = gen_filename(simulator, reference_nx)\n",
" with np.load(datafilename) as data:\n",
" h_ref = data['h']\n",
" \n",
" \n",
"\n",
" #Compute conservation and error for subsequent resolutions\n",
" for j, nx in enumerate(domain_sizes):\n",
" #Read current resolution\n",
" datafilename = gen_filename(simulator, nx)\n",
" with np.load(datafilename) as data:\n",
" h = data['h']\n",
" \n",
"\n",
" #Downsample reference to current resolution\n",
" h_ref_downsampled = downsample(h_ref, h_ref.size/h.size)\n",
" \n",
" h_ref_downsampled = downsample(h_ref, h_ref.size / h.size)\n",
"\n",
" #Compute error per cell\n",
" error[i, j] = np.linalg.norm((h_ref_downsampled - h).flatten(), ord=1) / nx\n",
" conservation[i, j] = (np.sum(h_ref_downsampled) - np.sum(h))\n",
" \n",
"\n",
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"\n",
"x = np.linspace(domain_sizes[0], domain_sizes[-1], 100);\n",
"x = np.linspace(domain_sizes[0], domain_sizes[-1], 100)\n",
"\n",
"scaling = np.max(error[:,0]) * domain_sizes[0] * 1.5\n",
"plt.loglog(x, scaling/x, '-', color='gray', label='Order 1')\n",
"scaling = np.max(error[:, 0]) * domain_sizes[0] * 1.5\n",
"plt.loglog(x, scaling / x, '-', color='gray', label='Order 1')\n",
"\n",
"scaling = np.min(error[:,0]) * domain_sizes[0]**2 * 0.5\n",
"plt.loglog(x, scaling/(x*x), '-', color='gray', label='Order 2')\n",
"scaling = np.min(error[:, 0]) * domain_sizes[0] ** 2 * 0.5\n",
"plt.loglog(x, scaling / (x * x), '-', color='gray', label='Order 2')\n",
"\n",
"for i, simulator in enumerate(simulators):\n",
" plt.loglog(domain_sizes, error[i,:], label=simulators[i].__name__)\n",
" \n",
"plt.axis([10**1, 10**4, 10**-8, 10**-2])\n",
" plt.loglog(domain_sizes, error[i, :], label=simulators[i].__name__)\n",
"\n",
"plt.axis([10 ** 1, 10 ** 4, 10 ** -8, 10 ** -2])\n",
"plt.xlabel('Number of cells')\n",
"plt.ylabel('L2 error')\n",
"plt.legend()\n",
@@ -616,14 +618,14 @@
],
"source": [
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"\n",
"for i, simulator in enumerate(simulators):\n",
" a = np.arange(1, len(domain_sizes))\n",
" order = np.log2(error[i, :-1] / error[i, 1:])\n",
" #order = np.log2(error[i, 0] / (error[i, 1:]) )/a\n",
" plt.plot(order, label=simulator.__name__)\n",
" \n",
"\n",
"plt.ylim([0, 3])\n",
"plt.xlabel('Number of cells')\n",
"plt.ylabel('Order')\n",
@@ -722,38 +724,39 @@
"source": [
"def plot_solution(simulator, nx, label, **kwargs):\n",
" datafilename = gen_filename(simulator, nx)\n",
" \n",
"\n",
" #Read the solution\n",
" with np.load(datafilename) as data:\n",
" h = data['h']\n",
" \n",
" x = np.linspace(0.5, nx-0.5, nx)*10.0/float(nx)\n",
" y = h[0,:]\n",
" \n",
"\n",
" x = np.linspace(0.5, nx - 0.5, nx) * 10.0 / float(nx)\n",
" y = h[0, :]\n",
"\n",
" plt.plot(x, y, label=label, **kwargs)\n",
" \n",
"\n",
" h = None\n",
" x = None\n",
" gc.collect() # Force run garbage collection to free up memory\n",
" gc.collect() # Force run garbage collection to free up memory\n",
"\n",
"\n",
"for simulator in simulators:\n",
" fig = plt.figure()\n",
" setBwStyles(fig.gca())\n",
" \n",
" set_bw_styles(fig.gca())\n",
"\n",
" #Plot reference solution\n",
" ref_x, ref_h, _ = gen_reference(4096)\n",
" plt.plot(ref_x, ref_h, '-', label='Reference')\n",
" ref_x = None\n",
" ref_h = None\n",
" gc.collect()\n",
" \n",
"\n",
" #Plot simulation results\n",
" for j, nx in enumerate(domain_sizes):\n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ')\n",
" plt.title(str(simulator.__name__))\n",
" plt.axis([0, 10, 0, 0.0055])\n",
" plt.legend()\n",
" \n",
"\n",
" save_figure(fig, simulator.__name__)"
]
},
@@ -774,9 +777,9 @@
}
],
"source": [
"figLegend = plt.figure()#figsize = (1.5,1.3))\n",
"figLegend = plt.figure() #figsize = (1.5,1.3))\n",
"plt.plot(0, 0)\n",
"plt.figlegend(*fig.gca().get_legend_handles_labels(), loc = 'center')\n",
"plt.figlegend(*fig.gca().get_legend_handles_labels(), loc='center')\n",
"plt.axis('off')\n",
"plt.show()"
]
@@ -850,20 +853,20 @@
" for i, simulator in enumerate(simulators):\n",
" nx = 512\n",
" plot_solution(simulator, nx, simulator.__name__, **kwargs)\n",
" \n",
"\n",
"\n",
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"\n",
"boxes = np.array([\\\n",
" [[3.5, 4.1], [0.0044, 0.00505]], \\\n",
" [[4.7, 5.05], [0.0025, 0.0027]], \\\n",
" [[6.1, 6.40], [0.0009, 0.0026]] \\\n",
" ])\n",
"boxes = np.array([\n",
" [[3.5, 4.1], [0.0044, 0.00505]],\n",
" [[4.7, 5.05], [0.0025, 0.0027]],\n",
" [[6.1, 6.40], [0.0009, 0.0026]]\n",
"])\n",
"\n",
"plot_comparison(marker=' ')\n",
"for box_x, box_y in boxes:\n",
" x = np.kron(box_x, np.ones((2)))\n",
" x = np.kron(box_x, np.ones(2))\n",
" y = np.hstack((box_y, box_y[::-1]))\n",
" fig.gca().fill(x, y, fill=True, linestyle='-', color='gray', alpha=0.5)\n",
"plt.xlim([0, 10])\n",
@@ -873,12 +876,12 @@
"\n",
"for i, bbox in enumerate(boxes):\n",
" fig = plt.figure()\n",
" setBwStyles(fig.gca())\n",
" set_bw_styles(fig.gca())\n",
" plot_comparison()\n",
" plt.xlim(bbox[0])\n",
" plt.ylim(bbox[1])\n",
" plt.legend()\n",
" \n",
"\n",
" save_figure(fig, \"solution_zoom_\" + str(i))\n"
]
},

291
ConvergenceSmooth1D.ipynb
View File

@@ -61,7 +61,7 @@
"source": [
"#Finally, import our simulator\n",
"from GPUSimulators import LxF, FORCE, HLL, HLL2, KP07, KP07_dimsplit, WAF\n",
"from GPUSimulators.common import common\n",
"from GPUSimulators.common import Timer\n",
"from GPUSimulators.helpers import InitialConditions"
]
},
@@ -100,13 +100,14 @@
"outputs": [],
"source": [
"#Misc plotting setup\n",
"def setBwStyles(ax):\n",
"def set_bw_styles(ax):\n",
" from cycler import cycler\n",
"\n",
" ax.set_prop_cycle( cycler('marker', ['.', 'x', 4, '+', '*', '1', 5]) +\n",
" cycler('linestyle', ['-.', '--', ':', '-.', '--', ':', '-.']) +\n",
" #cycler('markersize', [5, 5, 5, 5, 5, 5]) +\n",
" cycler('color', ['k', 'k', 'k', 'k', 'k', 'k', 'k']) ) \n",
" ax.set_prop_cycle(cycler('marker', ['.', 'x', 4, '+', '*', '1', 5]) +\n",
" cycler('linestyle', ['-.', '--', ':', '-.', '--', ':', '-.']) +\n",
" #cycler('markersize', [5, 5, 5, 5, 5, 5]) +\n",
" cycler('color', ['k', 'k', 'k', 'k', 'k', 'k', 'k']))\n",
"\n",
"\n",
"#Set large figure sizes\n",
"plt.rcParams['figure.figsize'] = [12, 8]\n",
@@ -125,26 +126,26 @@
"outputs": [],
"source": [
"def save_figure(fig, stem):\n",
" if (not os.path.isdir(\"figures\")):\n",
" if not os.path.isdir(\"figures\"):\n",
" os.mkdir(\"figures\")\n",
" \n",
"\n",
" fig_filename = os.path.join(\"figures\", \"convergence_smooth1d_\" + stem + \".pdf\")\n",
" \n",
"\n",
" metadata = {\n",
" 'created': time.strftime(\"%Y_%m_%d-%H_%M_%S\"),\n",
" 'hostname': socket.gethostname()\n",
" }\n",
" \n",
"\n",
" legend = fig.gca().legend_\n",
" if (legend != None):\n",
" fig.savefig(fig_filename, dpi=300,format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None, \n",
" metadata=metadata, \n",
" bbox_extra_artists=(legend, ), bbox_inches='tight')\n",
" if legend != None:\n",
" fig.savefig(fig_filename, dpi=300, format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None,\n",
" metadata=metadata,\n",
" bbox_extra_artists=(legend,), bbox_inches='tight')\n",
" else:\n",
" fig.savefig(fig_filename, dpi=300,format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None, \n",
" metadata=metadata)"
" fig.savefig(fig_filename, dpi=300, format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None,\n",
" metadata=metadata)"
]
},
{
@@ -168,21 +169,22 @@
"ny = 2\n",
"ref_nx = 1000\n",
"width = 50\n",
"height = width/5\n",
"height = width / 5\n",
"bump_size = 10\n",
"h, hu, hv, dx, dy = InitialConditions.bump(nx=nx, ny=ny, width=width, height=height, bump_size=bump_size, ref_nx=ref_nx, v_amp=0, v_ref=0)\n",
"h, hu, hv, dx, dy = InitialConditions.bump(nx=nx, ny=ny, width=width, height=height, bump_size=bump_size, ref_nx=ref_nx,\n",
" v_amp=0, v_ref=0)\n",
"\n",
"plt.figure()\n",
"ax=plt.subplot(1,3,1)\n",
"im=plt.imshow(h, extent=[0, nx*dx, 0, ny*dy])\n",
"ax = plt.subplot(1, 3, 1)\n",
"im = plt.imshow(h, extent=[0, nx * dx, 0, ny * dy])\n",
"plt.colorbar(orientation='horizontal')\n",
"\n",
"ax=plt.subplot(1,3,2)\n",
"im=plt.imshow(hu, extent=[0, nx*dx, 0, ny*dy])\n",
"ax = plt.subplot(1, 3, 2)\n",
"im = plt.imshow(hu, extent=[0, nx * dx, 0, ny * dy])\n",
"plt.colorbar(orientation='horizontal')\n",
"\n",
"ax=plt.subplot(1,3,3)\n",
"im=plt.imshow(hv, extent=[0, nx*dx, 0, ny*dy])\n",
"ax = plt.subplot(1, 3, 3)\n",
"im = plt.imshow(hv, extent=[0, nx * dx, 0, ny * dy])\n",
"plt.colorbar(orientation='horizontal')\n",
"\n",
"plt.show()\n"
@@ -207,18 +209,17 @@
"def run_benchmark(datafilename, simulator, simulator_args, nx, reference_nx, ny, reference_ny,\n",
" h_ref=0.5, h_amp=0.1, u_ref=0.0, u_amp=0.1, v_ref=0.0, v_amp=0.1,\n",
" dt=None, force_rerun=False, transpose=False):\n",
"\n",
" if (datafilename and os.path.isfile(datafilename) and force_rerun == False):\n",
" if datafilename and os.path.isfile(datafilename) and force_rerun == False:\n",
" print(\"WARNING: Existing simulation exists, skipping simulation\")\n",
" return [0, 0, 0]\n",
" else:\n",
" width=100\n",
" test_data_args= {\n",
" width = 100\n",
" test_data_args = {\n",
" 'nx': nx,\n",
" 'ny': ny,\n",
" 'width': width,\n",
" 'height': 2*width/nx,\n",
" 'bump_size': width/5,\n",
" 'height': 2 * width / nx,\n",
" 'bump_size': width / 5,\n",
" 'ref_nx': reference_nx,\n",
" 'ref_ny': reference_ny,\n",
" 'h_ref': h_ref, 'h_amp': h_amp,\n",
@@ -228,13 +229,13 @@
" h0, hu0, hv0, dx, dy = InitialConditions.bump(**test_data_args)\n",
"\n",
" #Initialize simulator\n",
" with Common.Timer(simulator.__name__ + \"_\" + str(nx)) as timer:\n",
" if (transpose):\n",
" with Timer(simulator.__name__ + \"_\" + str(nx)) as timer:\n",
" if transpose:\n",
" h0 = np.ascontiguousarray(np.transpose(h0))\n",
" hu0, hv0 = np.ascontiguousarray(np.transpose(hv0)), np.ascontiguousarray(np.transpose(hu0))\n",
" dx, dy = dy, dx\n",
" nx, ny = ny, nx\n",
" \n",
"\n",
" sim_args = {\n",
" 'h0': h0, 'hu0': hu0, 'hv0': hv0,\n",
" 'nx': nx, 'ny': ny,\n",
@@ -245,24 +246,23 @@
" sim = simulator(**sim_args)\n",
" sim.simulate(1.0, dt=dt)\n",
" sim.check()\n",
" \n",
"\n",
" nt = sim.sim_steps()\n",
" dt = sim.sim_time() / nt\n",
" h, hu, hv = sim.download()\n",
" \n",
" if (transpose):\n",
"\n",
" if transpose:\n",
" h = np.ascontiguousarray(np.transpose(h))\n",
" hu, hv = np.ascontiguousarray(np.transpose(hv)), np.ascontiguousarray(np.transpose(hu))\n",
" \n",
"\n",
" if (datafilename):\n",
" if datafilename:\n",
" dirname = os.path.dirname(datafilename)\n",
" if (dirname and not os.path.isdir(dirname)):\n",
" if dirname and not os.path.isdir(dirname):\n",
" os.makedirs(dirname)\n",
" np.savez_compressed(datafilename, h=h, hu=hu, hv=hv)\n",
" \n",
" gc.collect() # Force run garbage collection to free up memory\n",
" \n",
"\n",
" gc.collect() # Force run garbage collection to free up memory\n",
"\n",
" return [timer.secs, nt, dt]"
]
},
@@ -287,14 +287,14 @@
"domain_sizes = [8, 16, 32, 64, 128, 256, 512, 1024, 2048, 4096]\n",
"\n",
"simulators = [LxF.LxF, FORCE.FORCE, HLL.HLL, HLL2.HLL2, KP07.KP07, KP07_dimsplit.KP07_dimsplit, WAF.WAF]\n",
"force_rerun=True\n",
"transpose=False\n",
"force_rerun = True\n",
"transpose = False\n",
"g = 9.81\n",
"h_ref = 0.5\n",
"h_amp = 0.1\n",
"u_ref = 0.0\n",
"u_amp = 0.1\n",
"dt = 0.7*(width/reference_nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp)))\n",
"dt = 0.7 * (width / reference_nx) / (u_ref + u_amp + np.sqrt(g * (h_ref + h_amp)))\n",
"dt = None\n",
"sim_args = {\n",
" 'context': cuda_context,\n",
@@ -309,41 +309,40 @@
"sim_dt = np.zeros((len(simulators), len(domain_sizes)))\n",
"sim_nt = np.zeros((len(simulators), len(domain_sizes)))\n",
"\n",
"\n",
"#Warmup to compile kernels\n",
"with Common.Timer(\"warmup\") as timer:\n",
"with Timer(\"warmup\") as timer:\n",
" for i in range(len(simulators)):\n",
" run_benchmark(None, \n",
" simulators[i], \n",
" run_benchmark(None,\n",
" simulators[i],\n",
" sim_args,\n",
" 16, 16, 2, 2,\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
"print(\"Completed warmup in {:} s\".format(timer.secs))\n",
" \n",
"with Common.Timer(\"simulation\") as timer:\n",
" for i in range(len(simulators)): \n",
"print(f\"Completed warmup in {timer.secs} s\")\n",
"\n",
"with Timer(\"simulation\") as timer:\n",
" for i in range(len(simulators)):\n",
" # Run the simulation for all the domain sizes\n",
" for j, nx in enumerate(domain_sizes):\n",
" datafilename = gen_filename(simulators[i], nx)\n",
" [secs, nt, dt] = run_benchmark(datafilename, \n",
" simulators[i], \n",
" sim_args,\n",
" nx, reference_nx, 2, 2,\n",
" h_ref=h_ref, h_amp=h_amp,\n",
" u_ref=u_ref, u_amp=u_amp, \n",
" v_ref=0.0, v_amp=0.0,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
" [secs, nt, dt] = run_benchmark(datafilename,\n",
" simulators[i],\n",
" sim_args,\n",
" nx, reference_nx, 2, 2,\n",
" h_ref=h_ref, h_amp=h_amp,\n",
" u_ref=u_ref, u_amp=u_amp,\n",
" v_ref=0.0, v_amp=0.0,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
" sim_elapsed_time[i, j] = secs\n",
" sim_dt[i, j] = dt\n",
" sim_nt[i, j] = nt\n",
"\n",
" # Run reference with a low CFL-number\n",
" datafilename = gen_filename(simulators[i], reference_nx)\n",
" run_benchmark(datafilename, \n",
" simulators[i], \n",
" run_benchmark(datafilename,\n",
" simulators[i],\n",
" sim_args,\n",
" reference_nx, reference_nx, 2, 2,\n",
" h_ref=h_ref, h_amp=h_amp,\n",
@@ -352,7 +351,7 @@
" #dt=0.25*0.7*(width/reference_nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
"print(\"Completed simulation in {:} s\".format(timer.secs))"
"print(f\"Completed simulation in {timer.secs} s\")"
]
},
{
@@ -405,24 +404,24 @@
],
"source": [
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"plt.title(\"Elapsed time vs number of cells\")\n",
"for i, simulator in enumerate(simulators):\n",
" plt.plot(domain_sizes, sim_elapsed_time[i,:], label=simulator.__name__)\n",
" plt.plot(domain_sizes, sim_elapsed_time[i, :], label=simulator.__name__)\n",
"plt.legend()\n",
" \n",
"\n",
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"plt.title(\"Dt vs number of cells\")\n",
"for i, simulator in enumerate(simulators):\n",
" plt.loglog(domain_sizes, sim_dt[i,:], label=simulator.__name__)\n",
" plt.loglog(domain_sizes, sim_dt[i, :], label=simulator.__name__)\n",
"plt.legend()\n",
" \n",
"\n",
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"plt.title(\"Number of timesteps vs number of cells\")\n",
"for i, simulator in enumerate(simulators):\n",
" plt.plot(domain_sizes, sim_nt[i,:], label=simulator.__name__)\n",
" plt.plot(domain_sizes, sim_nt[i, :], label=simulator.__name__)\n",
"plt.legend()"
]
},
@@ -457,46 +456,46 @@
"source": [
"conservation = np.zeros((len(simulators), len(domain_sizes)))\n",
"error = np.zeros((len(simulators), len(domain_sizes)))\n",
"order = np.zeros((len(simulators), len(domain_sizes)-1))\n",
"order = np.zeros((len(simulators), len(domain_sizes) - 1))\n",
"\n",
"for i, simulator in enumerate(simulators):\n",
" #Read the reference solution\n",
" datafilename = gen_filename(simulator, reference_nx)\n",
" with np.load(datafilename) as data:\n",
" h_ref = data['h'].astype(np.float64)\n",
" \n",
"\n",
" #Compute conservation and error for subsequent resolutions\n",
" for j, nx in enumerate(domain_sizes):\n",
" #Read current resolution\n",
" datafilename = gen_filename(simulator, nx)\n",
" with np.load(datafilename) as data:\n",
" h = data['h'].astype(np.float64)\n",
" \n",
"\n",
" #Downsample reference to current resolution\n",
" h_ref_downsampled = InitialConditions.downsample(h_ref, reference_nx/nx, 1)\n",
" \n",
" h_ref_downsampled = InitialConditions.downsample(h_ref, reference_nx / nx, 1)\n",
"\n",
" #Compute error per cell\n",
" error[i, j] = np.linalg.norm((h_ref_downsampled - h).flatten(), ord=1) / nx\n",
" conservation[i, j] = (np.sum(h_ref) - np.sum(h))\n",
" \n",
"\n",
" #Compute the numerical order\n",
" a = np.arange(1, len(domain_sizes))\n",
" order[i,:] = np.log2(error[i, :-1] / error[i, 1:])\n",
" \n",
" order[i, :] = np.log2(error[i, :-1] / error[i, 1:])\n",
"\n",
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"\n",
"x = np.linspace(domain_sizes[0], domain_sizes[-1], 100);\n",
"x = np.linspace(domain_sizes[0], domain_sizes[-1], 100)\n",
"\n",
"scaling = np.max(error[:,0]) * domain_sizes[0] * 1.5\n",
"plt.loglog(x, scaling/x, '-', color='gray', label='Order 1')\n",
"scaling = np.max(error[:, 0]) * domain_sizes[0] * 1.5\n",
"plt.loglog(x, scaling / x, '-', color='gray', label='Order 1')\n",
"\n",
"scaling = np.min(error[:,0]) * domain_sizes[0]**2 * 0.5\n",
"plt.loglog(x, scaling/(x*x), '-', color='gray', label='Order 2')\n",
"scaling = np.min(error[:, 0]) * domain_sizes[0] ** 2 * 0.5\n",
"plt.loglog(x, scaling / (x * x), '-', color='gray', label='Order 2')\n",
"\n",
"for i, simulator in enumerate(simulators):\n",
" plt.loglog(domain_sizes, error[i,:], label=simulators[i].__name__)\n",
" \n",
" plt.loglog(domain_sizes, error[i, :], label=simulators[i].__name__)\n",
"\n",
"plt.xlabel('Number of cells')\n",
"plt.ylabel('Error')\n",
"plt.title(\"Convergence\")\n",
@@ -748,11 +747,12 @@
"table = pd.DataFrame(order, columns=domain_sizes[:-1], index=[simulator.__name__ for simulator in simulators])\n",
"display(table)\n",
"\n",
"order_summary = pd.DataFrame(np.empty((len(simulators), 3)), columns=['Max', 'Min', 'Mean'], index=[simulator.__name__ for simulator in simulators])\n",
"order_summary = pd.DataFrame(np.empty((len(simulators), 3)), columns=['Max', 'Min', 'Mean'],\n",
" index=[simulator.__name__ for simulator in simulators])\n",
"for i, simulator in enumerate(simulators):\n",
" order_summary.loc[simulator.__name__] = [np.max(order[i,:]), \n",
" np.min(order[i,:]), \n",
" np.mean(order[i,:])]\n",
" order_summary.loc[simulator.__name__] = [np.max(order[i, :]),\n",
" np.min(order[i, :]),\n",
" np.mean(order[i, :])]\n",
"display(order_summary)\n"
]
},
@@ -784,11 +784,11 @@
],
"source": [
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"\n",
"for i, simulator in enumerate(simulators):\n",
" plt.plot(order[i,:], label=simulator.__name__)\n",
" \n",
" plt.plot(order[i, :], label=simulator.__name__)\n",
"\n",
"plt.ylim([0, 3])\n",
"plt.xlabel('Number of cells')\n",
"plt.ylabel('Order')\n",
@@ -805,20 +805,21 @@
"source": [
"def plot_solution(simulator, nx, label, **kwargs):\n",
" datafilename = gen_filename(simulator, nx)\n",
" \n",
"\n",
" #Read the solution\n",
" with np.load(datafilename) as data:\n",
" h = data['h']\n",
" \n",
" x = np.linspace(0.5, nx-0.5, nx)*100.0/float(nx)\n",
" y = h[0,:]\n",
" \n",
"\n",
" x = np.linspace(0.5, nx - 0.5, nx) * 100.0 / float(nx)\n",
" y = h[0, :]\n",
"\n",
" plt.plot(x, y, label=label, **kwargs)\n",
" \n",
"\n",
" h = None\n",
" x = None\n",
" gc.collect() # Force run garbage collection to free up memory\n",
" \n",
" gc.collect() # Force run garbage collection to free up memory\n",
"\n",
"\n",
"def plot_comparison(nx, **kwargs):\n",
" plot_solution(HLL2.HLL2, reference_nx, 'Reference', marker=' ', linestyle='-')\n",
"\n",
@@ -845,41 +846,41 @@
}
],
"source": [
"boxes = np.array([\\\n",
" [[28, 35], [0.499, 0.516]], \\\n",
" [[45, 56], [0.58, 0.60]]\n",
" ])\n",
"boxes = np.array([\n",
" [[28, 35], [0.499, 0.516]],\n",
" [[45, 56], [0.58, 0.60]]\n",
"])\n",
"comparison_nx = 32\n",
"\n",
"fig = plt.figure(figsize=(12, 12))\n",
"gs = gridspec.GridSpec(2, 2)\n",
" \n",
"ax1 = fig.add_subplot(gs[0,:])\n",
"setBwStyles(ax1)\n",
"\n",
"ax1 = fig.add_subplot(gs[0, :])\n",
"set_bw_styles(ax1)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, comparison_nx, simulator.__name__, marker=' ')\n",
"for box_x, box_y in boxes:\n",
" x = np.kron(box_x, np.ones((2)))\n",
" x = np.kron(box_x, np.ones(2))\n",
" y = np.hstack((box_y, box_y[::-1]))\n",
" fig.gca().fill(x, y, fill=True, linestyle='-', color='gray', alpha=0.5)\n",
"plt.title(\"Comparison of solutions at nx=\" + str(comparison_nx))\n",
"plt.legend()\n",
"\n",
"ax2 = fig.add_subplot(gs[1,0])\n",
"setBwStyles(ax2)\n",
"ax2 = fig.add_subplot(gs[1, 0])\n",
"set_bw_styles(ax2)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, comparison_nx, simulator.__name__, markevery=2)\n",
"plt.xlim(boxes[0, 0])\n",
"plt.ylim(boxes[0, 1])\n",
"\n",
"ax3 = fig.add_subplot(gs[1,1])\n",
"setBwStyles(ax3)\n",
"setBwStyles(ax2)\n",
"ax3 = fig.add_subplot(gs[1, 1])\n",
"set_bw_styles(ax3)\n",
"set_bw_styles(ax2)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, comparison_nx, simulator.__name__, markevery=2)\n",
"plt.xlim(boxes[1, 0])\n",
"plt.ylim(boxes[1, 1])\n",
" \n",
"\n",
"save_figure(fig, 'solution')"
]
},
@@ -902,33 +903,33 @@
"source": [
"fig = plt.figure(figsize=(12, 12))\n",
"gs = gridspec.GridSpec(2, 2)\n",
" \n",
"ax1 = fig.add_subplot(gs[0,:])\n",
"setBwStyles(ax1)\n",
"\n",
"ax1 = fig.add_subplot(gs[0, :])\n",
"set_bw_styles(ax1)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, reference_nx, simulator.__name__, marker=' ')\n",
"for box_x, box_y in boxes:\n",
" x = np.kron(box_x, np.ones((2)))\n",
" x = np.kron(box_x, np.ones(2))\n",
" y = np.hstack((box_y, box_y[::-1]))\n",
" fig.gca().fill(x, y, fill=True, linestyle='-', color='gray', alpha=0.5)\n",
"plt.title(\"Comparison of reference solutions\")\n",
"plt.legend()\n",
"\n",
"ax2 = fig.add_subplot(gs[1,0])\n",
"setBwStyles(ax2)\n",
"ax2 = fig.add_subplot(gs[1, 0])\n",
"set_bw_styles(ax2)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, reference_nx, simulator.__name__, marker=' ')\n",
"plt.xlim(boxes[0, 0])\n",
"plt.ylim(boxes[0, 1])\n",
"\n",
"ax3 = fig.add_subplot(gs[1,1])\n",
"setBwStyles(ax3)\n",
"setBwStyles(ax2)\n",
"ax3 = fig.add_subplot(gs[1, 1])\n",
"set_bw_styles(ax3)\n",
"set_bw_styles(ax2)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, reference_nx, simulator.__name__, marker=' ')\n",
"plt.xlim(boxes[1, 0])\n",
"plt.ylim(boxes[1, 1])\n",
" \n",
"\n",
"save_figure(fig, 'reference_solutions')"
]
},
@@ -1013,37 +1014,37 @@
"source": [
"for simulator in simulators:\n",
" fig = plt.figure(figsize=(12, 12))\n",
" \n",
"\n",
" gs = gridspec.GridSpec(2, 2)\n",
" \n",
" ax1 = fig.add_subplot(gs[0,:])\n",
" setBwStyles(ax1)\n",
"\n",
" ax1 = fig.add_subplot(gs[0, :])\n",
" set_bw_styles(ax1)\n",
" plot_solution(simulator, reference_nx, 'Reference', linestyle='-', marker=' ')\n",
" for j, nx in enumerate(domain_sizes):\n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ')\n",
" for box_x, box_y in boxes:\n",
" x = np.kron(box_x, np.ones((2)))\n",
" x = np.kron(box_x, np.ones(2))\n",
" y = np.hstack((box_y, box_y[::-1]))\n",
" fig.gca().fill(x, y, fill=True, linestyle='-', color='gray', alpha=0.5)\n",
" plt.title(str(simulator.__name__))\n",
" plt.legend()\n",
" \n",
" ax2 = fig.add_subplot(gs[1,0])\n",
" setBwStyles(ax2)\n",
"\n",
" ax2 = fig.add_subplot(gs[1, 0])\n",
" set_bw_styles(ax2)\n",
" plot_solution(simulator, reference_nx, 'Reference', linestyle='-', marker=' ')\n",
" for j, nx in enumerate(domain_sizes):\n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ') \n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ')\n",
" plt.xlim(boxes[0, 0])\n",
" plt.ylim(boxes[0, 1])\n",
" \n",
" ax3 = fig.add_subplot(gs[1,1])\n",
" setBwStyles(ax3)\n",
"\n",
" ax3 = fig.add_subplot(gs[1, 1])\n",
" set_bw_styles(ax3)\n",
" plot_solution(simulator, reference_nx, 'Reference', linestyle='-', marker=' ')\n",
" for j, nx in enumerate(domain_sizes):\n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ') \n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ')\n",
" plt.xlim(boxes[1, 0])\n",
" plt.ylim(boxes[1, 1])\n",
" \n",
"\n",
" save_figure(fig, simulator.__name__)"
]
},

320
ConvergenceSmooth2D.ipynb
View File

@@ -61,7 +61,7 @@
"source": [
"#Finally, import our simulator\n",
"from GPUSimulators import LxF, FORCE, HLL, HLL2, KP07, KP07_dimsplit, WAF\n",
"from GPUSimulators.common import common\n",
"from GPUSimulators.common import Timer\n",
"from GPUSimulators.helpers import InitialConditions"
]
},
@@ -100,13 +100,14 @@
"outputs": [],
"source": [
"#Misc plotting setup\n",
"def setBwStyles(ax):\n",
"def set_bw_styles(ax):\n",
" from cycler import cycler\n",
"\n",
" ax.set_prop_cycle( cycler('marker', ['.', 'x', 4, '+', '*', '1', 5]) +\n",
" cycler('linestyle', ['-.', '--', ':', '-.', '--', ':', '-.']) +\n",
" #cycler('markersize', [5, 5, 5, 5, 5, 5]) +\n",
" cycler('color', ['k', 'k', 'k', 'k', 'k', 'k', 'k']) ) \n",
" ax.set_prop_cycle(cycler('marker', ['.', 'x', 4, '+', '*', '1', 5]) +\n",
" cycler('linestyle', ['-.', '--', ':', '-.', '--', ':', '-.']) +\n",
" #cycler('markersize', [5, 5, 5, 5, 5, 5]) +\n",
" cycler('color', ['k', 'k', 'k', 'k', 'k', 'k', 'k']))\n",
"\n",
"\n",
"#Set large figure sizes\n",
"plt.rcParams['figure.figsize'] = [12, 8]\n",
@@ -125,26 +126,26 @@
"outputs": [],
"source": [
"def save_figure(fig, stem):\n",
" if (not os.path.isdir(\"figures\")):\n",
" if not os.path.isdir(\"figures\"):\n",
" os.mkdir(\"figures\")\n",
" \n",
"\n",
" fig_filename = os.path.join(\"figures\", \"convergence_smooth2d_\" + stem + \".pdf\")\n",
" \n",
"\n",
" metadata = {\n",
" 'created': time.strftime(\"%Y_%m_%d-%H_%M_%S\"),\n",
" 'hostname': socket.gethostname()\n",
" }\n",
" \n",
"\n",
" legend = fig.gca().legend_\n",
" if (legend != None):\n",
" fig.savefig(fig_filename, dpi=300,format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None, \n",
" metadata=metadata, \n",
" bbox_extra_artists=(legend, ), bbox_inches='tight')\n",
" if legend != None:\n",
" fig.savefig(fig_filename, dpi=300, format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None,\n",
" metadata=metadata,\n",
" bbox_extra_artists=(legend,), bbox_inches='tight')\n",
" else:\n",
" fig.savefig(fig_filename, dpi=300,format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None, \n",
" metadata=metadata)"
" fig.savefig(fig_filename, dpi=300, format='pdf',\n",
" transparent=True, pad_inches=0.0, frameon=None,\n",
" metadata=metadata)"
]
},
{
@@ -164,26 +165,27 @@
}
],
"source": [
"nx=100\n",
"ref_nx=1000\n",
"ny=500\n",
"ref_ny=500\n",
"nx = 100\n",
"ref_nx = 1000\n",
"ny = 500\n",
"ref_ny = 500\n",
"width = 50\n",
"height = width\n",
"bump_size = 10\n",
"h, hu, hv, dx, dy = InitialConditions.bump(nx=nx, ny=ny, ref_nx=ref_nx, ref_ny=ref_ny, width=width, height=height, bump_size=bump_size)\n",
"h, hu, hv, dx, dy = InitialConditions.bump(nx=nx, ny=ny, ref_nx=ref_nx, ref_ny=ref_ny, width=width, height=height,\n",
" bump_size=bump_size)\n",
"\n",
"plt.figure()\n",
"ax=plt.subplot(1,3,1)\n",
"im=plt.imshow(h, extent=[0, nx*dx, 0, ny*dy])\n",
"ax = plt.subplot(1, 3, 1)\n",
"im = plt.imshow(h, extent=[0, nx * dx, 0, ny * dy])\n",
"plt.colorbar(orientation='horizontal')\n",
"\n",
"ax=plt.subplot(1,3,2)\n",
"im=plt.imshow(hu, extent=[0, nx*dx, 0, ny*dy])\n",
"ax = plt.subplot(1, 3, 2)\n",
"im = plt.imshow(hu, extent=[0, nx * dx, 0, ny * dy])\n",
"plt.colorbar(orientation='horizontal')\n",
"\n",
"ax=plt.subplot(1,3,3)\n",
"im=plt.imshow(hv, extent=[0, nx*dx, 0, ny*dy])\n",
"ax = plt.subplot(1, 3, 3)\n",
"im = plt.imshow(hv, extent=[0, nx * dx, 0, ny * dy])\n",
"plt.colorbar(orientation='horizontal')\n",
"\n",
"plt.show()"
@@ -208,12 +210,11 @@
"def run_benchmark(datafilename, simulator, simulator_args, nx, reference_nx, ny, reference_ny,\n",
" h_ref=0.5, h_amp=0.1, u_ref=0.0, u_amp=0.1, v_ref=0.0, v_amp=0.1,\n",
" dt=None, force_rerun=False, transpose=False):\n",
"\n",
" if (datafilename and os.path.isfile(datafilename) and force_rerun == False):\n",
" if datafilename and os.path.isfile(datafilename) and force_rerun == False:\n",
" print(\"WARNING: Existing simulation exists, skipping simulation\")\n",
" return [0, 0, 0]\n",
" else:\n",
" test_data_args= {\n",
" test_data_args = {\n",
" 'nx': nx,\n",
" 'ny': ny,\n",
" 'width': 100,\n",
@@ -228,13 +229,13 @@
" h0, hu0, hv0, dx, dy = InitialConditions.bump(**test_data_args)\n",
"\n",
" #Initialize simulator\n",
" with Common.Timer(simulator.__name__ + \"_\" + str(nx)) as timer:\n",
" if (transpose):\n",
" with Timer(simulator.__name__ + \"_\" + str(nx)) as timer:\n",
" if transpose:\n",
" h0 = np.ascontiguousarray(np.transpose(h0))\n",
" hu0, hv0 = np.ascontiguousarray(np.transpose(hv0)), np.ascontiguousarray(np.transpose(hu0))\n",
" dx, dy = dy, dx\n",
" nx, ny = ny, nx\n",
" \n",
"\n",
" sim_args = {\n",
" 'h0': h0, 'hu0': hu0, 'hv0': hv0,\n",
" 'nx': nx, 'ny': ny,\n",
@@ -245,24 +246,23 @@
" sim = simulator(**sim_args)\n",
" sim.simulate(1.0, dt=dt)\n",
" sim.check()\n",
" \n",
"\n",
" nt = sim.sim_steps()\n",
" dt = sim.sim_time() / nt\n",
" h, hu, hv = sim.download()\n",
" \n",
" if (transpose):\n",
"\n",
" if transpose:\n",
" h = np.ascontiguousarray(np.transpose(h))\n",
" hu, hv = np.ascontiguousarray(np.transpose(hv)), np.ascontiguousarray(np.transpose(hu))\n",
" \n",
"\n",
" if (datafilename):\n",
" if datafilename:\n",
" dirname = os.path.dirname(datafilename)\n",
" if (dirname and not os.path.isdir(dirname)):\n",
" if dirname and not os.path.isdir(dirname):\n",
" os.makedirs(dirname)\n",
" np.savez_compressed(datafilename, h=h, hu=hu, hv=hv)\n",
" \n",
" gc.collect() # Force run garbage collection to free up memory\n",
" \n",
"\n",
" gc.collect() # Force run garbage collection to free up memory\n",
"\n",
" return [timer.secs, nt, dt]"
]
},
@@ -310,26 +310,25 @@
" 'g': 9.81\n",
"}\n",
"\n",
"\n",
"elapsed, nt, dt = run_benchmark(datafilename, \n",
" KP07.KP07, \n",
" sim_args,\n",
" nx=50, reference_nx=500, ny=100, reference_ny=200,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=True,\n",
" transpose=False)\n",
"elapsed, nt, dt = run_benchmark(datafilename,\n",
" KP07.KP07,\n",
" sim_args,\n",
" nx=50, reference_nx=500, ny=100, reference_ny=200,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=True,\n",
" transpose=False)\n",
"\n",
"with np.load(datafilename) as data:\n",
" h = data['h']\n",
" rh = InitialConditions.downsample(h, 2, 4)\n",
" \n",
"\n",
" print(str(h.shape) + \" =downsample=> \" + str(rh.shape))\n",
"\n",
" plt.figure(figsize=(8, 6))\n",
" plt.imshow(rh, extent=[0, 100, 0, 100])\n",
" plt.colorbar()\n",
" plt.show()\n",
" \n",
"\n",
"gc.collect()"
]
},
@@ -410,55 +409,53 @@
],
"source": [
"reference_nx = 2048\n",
"force_rerun=True\n",
"transpose=True\n",
"force_rerun = True\n",
"transpose = True\n",
"domain_sizes = [8, 16, 32, 64, 128, 256, 512, 1024]\n",
"simulators = [LxF.LxF, FORCE.FORCE, HLL.HLL, HLL2.HLL2, WAF.WAF, KP07.KP07, KP07_dimsplit.KP07_dimsplit]\n",
"\n",
"\n",
"\n",
"sim_elapsed_time = np.zeros((len(simulators), len(domain_sizes)))\n",
"sim_dt = np.zeros((len(simulators), len(domain_sizes)))\n",
"sim_nt = np.zeros((len(simulators), len(domain_sizes)))\n",
"\n",
"#Warmup to compile kernels\n",
"with Common.Timer(\"warmup\") as timer:\n",
"with Timer(\"warmup\") as timer:\n",
" for i in range(len(simulators)):\n",
" run_benchmark(datafilename=None, \n",
" simulator=simulators[i],\n",
" simulator_args=sim_args,\n",
" nx=16, reference_nx=16, ny=16, reference_ny=16,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
" run_benchmark(datafilename=None,\n",
" simulator=simulators[i],\n",
" simulator_args=sim_args,\n",
" nx=16, reference_nx=16, ny=16, reference_ny=16,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
"print(\"Completed warmup in {:.2f} s\".format(timer.secs))\n",
" \n",
"with Common.Timer(\"simulation\") as timer:\n",
"\n",
"with Timer(\"simulation\") as timer:\n",
" for i in range(len(simulators)):\n",
" with Common.Timer(simulators[i].__name__) as sim_timer:\n",
" with Timer(simulators[i].__name__) as sim_timer:\n",
" # Run the simulation for all the domain sizes\n",
" for j, nx in enumerate(domain_sizes):\n",
" datafilename = gen_filename(simulators[i], nx)\n",
" secs, nt, dt = run_benchmark(datafilename=datafilename, \n",
" simulator=simulators[i],\n",
" simulator_args=sim_args,\n",
" nx=nx, ny=nx, reference_nx=reference_nx, reference_ny=reference_nx,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
" secs, nt, dt = run_benchmark(datafilename=datafilename,\n",
" simulator=simulators[i],\n",
" simulator_args=sim_args,\n",
" nx=nx, ny=nx, reference_nx=reference_nx, reference_ny=reference_nx,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
" sim_elapsed_time[i, j] = secs\n",
" sim_dt[i, j] = dt\n",
" sim_nt[i, j] = nt\n",
"\n",
" # Run reference with a low CFL-number\n",
" datafilename = gen_filename(simulators[i], reference_nx)\n",
" run_benchmark(datafilename=datafilename, \n",
" simulator=simulators[i],\n",
" simulator_args=sim_args,\n",
" nx=reference_nx, ny=reference_nx, reference_nx=reference_nx, reference_ny=reference_nx,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
" run_benchmark(datafilename=datafilename,\n",
" simulator=simulators[i],\n",
" simulator_args=sim_args,\n",
" nx=reference_nx, ny=reference_nx, reference_nx=reference_nx, reference_ny=reference_nx,\n",
" #dt=0.25*0.7*(width/nx)/(u_ref+u_amp + np.sqrt(g*(h_ref+h_amp))),\n",
" force_rerun=force_rerun,\n",
" transpose=transpose)\n",
" print(\"Completed simulation ({:}) in {:.2f} s\".format(simulators[i].__name__, sim_timer.secs))\n",
"print(\"Completed simulation in {:.2f} s\".format(timer.secs))"
]
@@ -513,24 +510,24 @@
],
"source": [
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"plt.title(\"Elapsed time vs number of cells\")\n",
"for i, simulator in enumerate(simulators):\n",
" plt.plot(domain_sizes, sim_elapsed_time[i,:], label=simulator.__name__)\n",
" plt.plot(domain_sizes, sim_elapsed_time[i, :], label=simulator.__name__)\n",
"plt.legend()\n",
" \n",
"\n",
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"plt.title(\"Dt vs number of cells\")\n",
"for i, simulator in enumerate(simulators):\n",
" plt.loglog(domain_sizes, sim_dt[i,:], label=simulator.__name__)\n",
" plt.loglog(domain_sizes, sim_dt[i, :], label=simulator.__name__)\n",
"plt.legend()\n",
" \n",
"\n",
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"plt.title(\"Number of timesteps vs number of cells\")\n",
"for i, simulator in enumerate(simulators):\n",
" plt.plot(domain_sizes, sim_nt[i,:], label=simulator.__name__)\n",
" plt.plot(domain_sizes, sim_nt[i, :], label=simulator.__name__)\n",
"plt.legend()"
]
},
@@ -769,42 +766,42 @@
"\n",
"conservation = np.zeros((len(simulators), len(domain_sizes)))\n",
"error = np.zeros((len(simulators), len(domain_sizes)))\n",
"order = np.zeros((len(simulators), len(domain_sizes)-1))\n",
"order = np.zeros((len(simulators), len(domain_sizes) - 1))\n",
"\n",
"for i, simulator in enumerate(simulators):\n",
" #Read the reference solution\n",
" datafilename = gen_filename(simulator, reference_nx)\n",
" with np.load(datafilename) as data:\n",
" h_ref = data['h'].astype(np.float64)\n",
" \n",
"\n",
" #Compute conservation and error for subsequent resolutions\n",
" for j, nx in enumerate(domain_sizes):\n",
" #Read current resolution\n",
" datafilename = gen_filename(simulator, nx)\n",
" with np.load(datafilename) as data:\n",
" h = data['h'].astype(np.float64)\n",
" \n",
"\n",
" #Downsample reference to current resolution\n",
" h_ref_downsampled = InitialConditions.downsample(h_ref, reference_nx/nx)\n",
" \n",
" h_ref_downsampled = InitialConditions.downsample(h_ref, reference_nx / nx)\n",
"\n",
" #Compute error per cell\n",
" error[i, j] = np.linalg.norm((h_ref_downsampled - h).flatten(), ord=1) / (nx*nx)\n",
" error[i, j] = np.linalg.norm((h_ref_downsampled - h).flatten(), ord=1) / (nx * nx)\n",
" conservation[i, j] = (np.sum(h_ref_downsampled) - np.sum(h))\n",
" \n",
"\n",
" #Compute the numerical order\n",
" a = np.arange(1, len(domain_sizes))\n",
" order[i,:] = np.log2(error[i, :-1] / error[i, 1:])\n",
"\n",
" order[i, :] = np.log2(error[i, :-1] / error[i, 1:])\n",
"\n",
"table = pd.DataFrame(order, columns=domain_sizes[:-1], index=[simulator.__name__ for simulator in simulators])\n",
"display(table)\n",
"\n",
"order_summary = pd.DataFrame(np.empty((len(simulators), 4)), columns=['Max', 'Min', 'Mean', 'Stdev'], index=[simulator.__name__ for simulator in simulators])\n",
"order_summary = pd.DataFrame(np.empty((len(simulators), 4)), columns=['Max', 'Min', 'Mean', 'Stdev'],\n",
" index=[simulator.__name__ for simulator in simulators])\n",
"for i, simulator in enumerate(simulators):\n",
" order_summary.loc[simulator.__name__] = [np.max(order[i,:]), \n",
" np.min(order[i,:]), \n",
" np.mean(order[i,:]), \n",
" np.std(order[i,:])]\n",
" order_summary.loc[simulator.__name__] = [np.max(order[i, :]),\n",
" np.min(order[i, :]),\n",
" np.mean(order[i, :]),\n",
" np.std(order[i, :])]\n",
"display(order_summary)"
]
},
@@ -836,18 +833,18 @@
],
"source": [
"fig = plt.figure()\n",
"setBwStyles(fig.gca())\n",
"set_bw_styles(fig.gca())\n",
"\n",
"x = np.linspace(domain_sizes[0], domain_sizes[-1], 100);\n",
"\n",
"scaling = np.max(error[:,0]) * domain_sizes[0] * 1.5\n",
"plt.loglog(x, scaling/x, '-', color='gray', label='Order 1')\n",
"scaling = np.max(error[:, 0]) * domain_sizes[0] * 1.5\n",
"plt.loglog(x, scaling / x, '-', color='gray', label='Order 1')\n",
"\n",
"scaling = np.min(error[:,0]) * domain_sizes[0]**2 * 0.5\n",
"plt.loglog(x, scaling/(x*x), '-', color='gray', label='Order 2')\n",
"scaling = np.min(error[:, 0]) * domain_sizes[0] ** 2 * 0.5\n",
"plt.loglog(x, scaling / (x * x), '-', color='gray', label='Order 2')\n",
"\n",
"for i, simulator in enumerate(simulators):\n",
" plt.loglog(domain_sizes, error[i,:], label=simulators[i].__name__)\n",
" plt.loglog(domain_sizes, error[i, :], label=simulators[i].__name__)\n",
"\n",
"plt.xlabel('Number of cells')\n",
"plt.ylabel('Error')\n",
@@ -864,20 +861,21 @@
"source": [
"def plot_solution(simulator, nx, label, **kwargs):\n",
" datafilename = gen_filename(simulator, nx)\n",
" \n",
"\n",
" #Read the solution\n",
" with np.load(datafilename) as data:\n",
" h = data['h']\n",
" h = h[int(nx/2),:]\n",
" \n",
" x = np.linspace(0.5, nx-0.5, nx)*100.0/float(nx)\n",
" \n",
" h = h[int(nx / 2), :]\n",
"\n",
" x = np.linspace(0.5, nx - 0.5, nx) * 100.0 / float(nx)\n",
"\n",
" plt.plot(x, h, label=label, **kwargs)\n",
" \n",
"\n",
" h = None\n",
" x = None\n",
" gc.collect() # Force run garbage collection to free up memory\n",
" \n",
" gc.collect() # Force run garbage collection to free up memory\n",
"\n",
"\n",
"def plot_comparison(nx, **kwargs):\n",
" plot_solution(HLL2.HLL2, reference_nx, 'Reference', marker=' ', linestyle='-')\n",
"\n",
@@ -902,41 +900,41 @@
}
],
"source": [
"boxes = np.array([\\\n",
" [[28, 35], [0.499, 0.516]], \\\n",
" [[45, 56], [0.58, 0.60]]\n",
" ])\n",
"boxes = np.array([\n",
" [[28, 35], [0.499, 0.516]],\n",
" [[45, 56], [0.58, 0.60]]\n",
"])\n",
"comparison_nx = 32\n",
"\n",
"fig = plt.figure(figsize=(12, 12))\n",
"gs = gridspec.GridSpec(2, 2)\n",
" \n",
"ax1 = fig.add_subplot(gs[0,:])\n",
"setBwStyles(ax1)\n",
"\n",
"ax1 = fig.add_subplot(gs[0, :])\n",
"set_bw_styles(ax1)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, comparison_nx, simulator.__name__, marker=' ')\n",
"for box_x, box_y in boxes:\n",
" x = np.kron(box_x, np.ones((2)))\n",
" x = np.kron(box_x, np.ones(2))\n",
" y = np.hstack((box_y, box_y[::-1]))\n",
" fig.gca().fill(x, y, fill=True, linestyle='-', color='gray', alpha=0.5)\n",
"plt.title(\"Comparison of solutions at nx=\" + str(comparison_nx))\n",
"plt.legend()\n",
"\n",
"ax2 = fig.add_subplot(gs[1,0])\n",
"setBwStyles(ax2)\n",
"ax2 = fig.add_subplot(gs[1, 0])\n",
"set_bw_styles(ax2)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, comparison_nx, simulator.__name__, markevery=2)\n",
"plt.xlim(boxes[0, 0])\n",
"plt.ylim(boxes[0, 1])\n",
"\n",
"ax3 = fig.add_subplot(gs[1,1])\n",
"setBwStyles(ax3)\n",
"setBwStyles(ax2)\n",
"ax3 = fig.add_subplot(gs[1, 1])\n",
"set_bw_styles(ax3)\n",
"set_bw_styles(ax2)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, comparison_nx, simulator.__name__, markevery=2)\n",
"plt.xlim(boxes[1, 0])\n",
"plt.ylim(boxes[1, 1])\n",
" \n",
"\n",
"save_figure(fig, 'solution')"
]
},
@@ -959,33 +957,33 @@
"source": [
"fig = plt.figure(figsize=(12, 12))\n",
"gs = gridspec.GridSpec(2, 2)\n",
" \n",
"ax1 = fig.add_subplot(gs[0,:])\n",
"setBwStyles(ax1)\n",
"\n",
"ax1 = fig.add_subplot(gs[0, :])\n",
"set_bw_styles(ax1)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, reference_nx, simulator.__name__, marker=' ')\n",
"for box_x, box_y in boxes:\n",
" x = np.kron(box_x, np.ones((2)))\n",
" x = np.kron(box_x, np.ones(2))\n",
" y = np.hstack((box_y, box_y[::-1]))\n",
" fig.gca().fill(x, y, fill=True, linestyle='-', color='gray', alpha=0.5)\n",
"plt.title(\"Comparison of reference solutions\")\n",
"plt.legend()\n",
"\n",
"ax2 = fig.add_subplot(gs[1,0])\n",
"setBwStyles(ax2)\n",
"ax2 = fig.add_subplot(gs[1, 0])\n",
"set_bw_styles(ax2)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, reference_nx, simulator.__name__, marker=' ')\n",
"plt.xlim(boxes[0, 0])\n",
"plt.ylim(boxes[0, 1])\n",
"\n",
"ax3 = fig.add_subplot(gs[1,1])\n",
"setBwStyles(ax3)\n",
"setBwStyles(ax2)\n",
"ax3 = fig.add_subplot(gs[1, 1])\n",
"set_bw_styles(ax3)\n",
"set_bw_styles(ax2)\n",
"for simulator in simulators:\n",
" plot_solution(simulator, reference_nx, simulator.__name__, marker=' ')\n",
"plt.xlim(boxes[1, 0])\n",
"plt.ylim(boxes[1, 1])\n",
" \n",
"\n",
"save_figure(fig, 'reference_solutions')"
]
},
@@ -1068,44 +1066,44 @@
}
],
"source": [
"boxes = np.array([\\\n",
" [[28, 35], [0.499, 0.516]], \\\n",
" [[45, 56], [0.58, 0.60]]\n",
" ])\n",
"boxes = np.array([\n",
" [[28, 35], [0.499, 0.516]],\n",
" [[45, 56], [0.58, 0.60]]\n",
"])\n",
"\n",
"for simulator in simulators:\n",
" fig = plt.figure(figsize=(12, 12))\n",
" \n",
"\n",
" gs = gridspec.GridSpec(2, 2)\n",
" \n",
" ax1 = fig.add_subplot(gs[0,:])\n",
" setBwStyles(ax1)\n",
"\n",
" ax1 = fig.add_subplot(gs[0, :])\n",
" set_bw_styles(ax1)\n",
" plot_solution(simulator, reference_nx, 'Reference', linestyle='-', marker=' ')\n",
" for j, nx in enumerate(domain_sizes):\n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ')\n",
" for box_x, box_y in boxes:\n",
" x = np.kron(box_x, np.ones((2)))\n",
" x = np.kron(box_x, np.ones(2))\n",
" y = np.hstack((box_y, box_y[::-1]))\n",
" fig.gca().fill(x, y, fill=True, linestyle='-', color='gray', alpha=0.5)\n",
" plt.title(str(simulator.__name__))\n",
" plt.legend()\n",
" \n",
" ax2 = fig.add_subplot(gs[1,0])\n",
" setBwStyles(ax2)\n",
"\n",
" ax2 = fig.add_subplot(gs[1, 0])\n",
" set_bw_styles(ax2)\n",
" plot_solution(simulator, reference_nx, 'Reference', linestyle='-', marker=' ')\n",
" for j, nx in enumerate(domain_sizes):\n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ') \n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ')\n",
" plt.xlim(boxes[0, 0])\n",
" plt.ylim(boxes[0, 1])\n",
" \n",
" ax3 = fig.add_subplot(gs[1,1])\n",
" setBwStyles(ax3)\n",
"\n",
" ax3 = fig.add_subplot(gs[1, 1])\n",
" set_bw_styles(ax3)\n",
" plot_solution(simulator, reference_nx, 'Reference', linestyle='-', marker=' ')\n",
" for j, nx in enumerate(domain_sizes):\n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ') \n",
" plot_solution(simulator, nx, 'nx=' + str(nx), marker=' ')\n",
" plt.xlim(boxes[1, 0])\n",
" plt.ylim(boxes[1, 1])\n",
" \n",
"\n",
" save_figure(fig, simulator.__name__)"
]
},

52
Desingularization.ipynb
View File

@@ -30,13 +30,13 @@
"metadata": {},
"outputs": [],
"source": [
"def setBwStyles(ax):\n",
"def set_bw_styles(ax):\n",
" from cycler import cycler\n",
"\n",
" ax.set_prop_cycle( #cycler('marker', ['.', 'x', 4, '+', '*', '1']) +\n",
" #cycler('linestyle', ['-.', '--', ':', '-.', '--', ':']) +\n",
" #cycler('markersize', [15, 15, 15, 15, 15, 15]) +\n",
" cycler('color', ['k', 'k', 'k', 'k', 'k', 'k']) )"
" ax.set_prop_cycle( #cycler('marker', ['.', 'x', 4, '+', '*', '1']) +\n",
" #cycler('linestyle', ['-.', '--', ':', '-.', '--', ':']) +\n",
" #cycler('markersize', [15, 15, 15, 15, 15, 15]) +\n",
" cycler('color', ['k', 'k', 'k', 'k', 'k', 'k']))"
]
},
{
@@ -87,46 +87,47 @@
],
"source": [
"def desingularise_our(hu, h, eps):\n",
" return hu / np.maximum(np.minimum(h*h/(2.0*eps)+0.5*eps, eps), np.abs(h))\n",
" return hu / np.maximum(np.minimum(h * h / (2.0 * eps) + 0.5 * eps, eps), np.abs(h))\n",
"\n",
"\n",
"def desingularise_kp1(hu, h, eps):\n",
" return np.sqrt(2.0)*h*hu/np.sqrt(h**4+np.maximum(h**4, eps**4))\n",
" return np.sqrt(2.0) * h * hu / np.sqrt(h ** 4 + np.maximum(h ** 4, eps ** 4))\n",
"\n",
"\n",
"def desingularise_kp2(hu, h, eps):\n",
" return h*hu/(h**2+eps**2)\n",
" return h * hu / (h ** 2 + eps ** 2)\n",
"\n",
"\n",
"def desingularise_kp3(hu, h, eps):\n",
" return 2*h*hu/(h**2+np.maximum(h**2, eps**2))\n",
"\n",
"\n",
" return 2 * h * hu / (h ** 2 + np.maximum(h ** 2, eps ** 2))\n",
"\n",
"\n",
"def desingularised_h_our(hu, h, eps):\n",
" return np.maximum(np.minimum(h*h/(2.0*eps)+0.5*eps, eps), np.abs(h))\n",
" return np.maximum(np.minimum(h * h / (2.0 * eps) + 0.5 * eps, eps), np.abs(h))\n",
"\n",
"\n",
"def desingularised_h_kp1(hu, h, eps):\n",
" return np.sqrt(h**4+np.maximum(h**4, eps**4)) / (np.sqrt(2.0)*h)\n",
" return np.sqrt(h ** 4 + np.maximum(h ** 4, eps ** 4)) / (np.sqrt(2.0) * h)\n",
"\n",
"\n",
"def desingularised_h_kp2(hu, h, eps):\n",
" return (h**2+eps**2) / h\n",
" return (h ** 2 + eps ** 2) / h\n",
"\n",
"\n",
"def desingularised_h_kp3(hu, h, eps):\n",
" return (h**2+np.maximum(h**2, eps**2)) / (2*h)\n",
"\n",
"\n",
" return (h ** 2 + np.maximum(h ** 2, eps ** 2)) / (2 * h)\n",
"\n",
"\n",
"x = np.linspace(0, 3.0e-11, 500).astype(np.float32)\n",
"h = x.copy()\n",
"u = np.ones_like(h)*1.0\n",
"hu = h*u\n",
"u_true = hu/h\n",
"u = np.ones_like(h) * 1.0\n",
"hu = h * u\n",
"u_true = hu / h\n",
"u_our = desingularise_our(hu, h, 1.0e-11)\n",
"u_kp1 = desingularise_kp1(hu, h, 1.0e-11)\n",
"u_kp2 = desingularise_kp2(hu, h, 1.0e-11)\n",
"u_kp3 = desingularise_kp3(hu, h, 1.0e-11)\n",
"\n",
"\n",
"fig = plt.figure(figsize=(5, 4), dpi=150)\n",
"#setBwStyles(fig.gca())\n",
"plt.plot(x, u_kp1, label='$u^*_1$', markevery=50, marker='.', color='k', markersize=10)\n",
@@ -139,13 +140,12 @@
"plt.xlabel('h')\n",
"plt.ylabel('Desingularized $hu/h$')\n",
"plt.legend()\n",
"plt.savefig(\"desing_hu.pdf\", bbox_inches = \"tight\")\n",
"plt.savefig(\"desing_hu.pdf\", bbox_inches=\"tight\")\n",
"plt.show()\n",
"\n",
"\n",
"h = x.copy()\n",
"u = np.ones_like(h)*1.01\n",
"hu = h*u\n",
"u = np.ones_like(h) * 1.01\n",
"hu = h * u\n",
"h_true = h\n",
"h_our = desingularised_h_our(hu, h, 1.0e-11)\n",
"h_kp1 = desingularised_h_kp1(hu, h, 1.0e-11)\n",
@@ -164,7 +164,7 @@
"plt.xlabel('h')\n",
"plt.ylabel('Desingularized $h$')\n",
"plt.legend()\n",
"plt.savefig(\"desing_h.pdf\", bbox_inches = \"tight\")\n",
"plt.savefig(\"desing_h.pdf\", bbox_inches=\"tight\")\n",
"plt.show()"
]
},

265
EulerTesting.ipynb
View File

@@ -1,8 +1,10 @@
{
"cells": [
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"#Lets have matplotlib \"inline\"\n",
"%matplotlib inline\n",
@@ -28,21 +30,17 @@
"except ImportError:\n",
" from io import StringIO\n",
"\n",
"from GPUSimulators.common import common\n",
"from GPUSimulators.common import Timer, DataDumper, ProgressPrinter\n",
"from GPUSimulators.EE2D_KP07_dimsplit import EE2D_KP07_dimsplit\n",
"from GPUSimulators.helpers import InitialConditions, Visualization"
],
"outputs": [],
"execution_count": null
]
},
{
"cell_type": "code",
"metadata": {},
"source": [
"%setup_logging --out test_schemes.log logger"
],
"cell_type": "code",
"outputs": [],
"execution_count": null
"execution_count": null,
"source": "%setup_logging --out test_schemes.log logger"
},
{
"cell_type": "code",
@@ -66,56 +64,62 @@
"execution_count": null
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"def plotVars(rho, rho_u, rho_v, E):\n",
"def plot_vars(rho, rho_u, rho_v, E):\n",
" plt.figure()\n",
" plt.subplot(1,4,1)\n",
" plt.subplot(1, 4, 1)\n",
" plt.imshow(rho, origin='lower')\n",
" plt.colorbar(orientation='horizontal', pad=0.02, shrink=0.9)\n",
"\n",
" plt.subplot(1,4,2)\n",
" plt.subplot(1, 4, 2)\n",
" plt.imshow(rho_u, origin='lower')\n",
" plt.colorbar(orientation='horizontal', pad=0.02, shrink=0.9)\n",
"\n",
" plt.subplot(1,4,3)\n",
" plt.subplot(1, 4, 3)\n",
" plt.imshow(rho_v, origin='lower')\n",
" plt.colorbar(orientation='horizontal', pad=0.02, shrink=0.9)\n",
"\n",
" plt.subplot(1,4,4)\n",
" plt.subplot(1, 4, 4)\n",
" plt.imshow(E, origin='lower')\n",
" plt.colorbar(orientation='horizontal', pad=0.02, shrink=0.9)"
],
"outputs": [],
"execution_count": null
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"def runSimulation(outfile, t_end, sim_args):\n",
" with Common.Timer(\"construct\") as t:\n",
"def run_simulation(outfile, t_end, sim_args):\n",
" with Timer(\"construct\") as t:\n",
" sim = EE2D_KP07_dimsplit(**sim_args)\n",
" print(\"Constructed in \" + str(t.secs) + \" seconds\")\n",
"\n",
" #Create netcdf file and simulate\n",
" with Common.DataDumper(outfile, mode='w', clobber=False) as outdata:\n",
" #Create a netcdf file and simulate\n",
" with DataDumper(outfile, mode='w', clobber=False) as outdata:\n",
" outdata.ncfile.createDimension('time', None)\n",
" outdata.ncfile.createDimension('x', sim.nx)\n",
" outdata.ncfile.createDimension('y', sim.ny)\n",
"\n",
" #Create variables\n",
" outdata.time_var = outdata.ncfile.createVariable('time', np.dtype('float32').char, 'time')\n",
" outdata.x_var = outdata.ncfile.createVariable('x', np.dtype('float32').char, 'x')\n",
" outdata.y_var = outdata.ncfile.createVariable('y', np.dtype('float32').char, 'y')\n",
" outdata.rho_var = outdata.ncfile.createVariable('rho', np.dtype('float32').char, ('time', 'y', 'x'), zlib=True, least_significant_digit=3)\n",
" outdata.rho_u_var = outdata.ncfile.createVariable('rho_u', np.dtype('float32').char, ('time', 'y', 'x'), zlib=True, least_significant_digit=3)\n",
" outdata.rho_v_var = outdata.ncfile.createVariable('rho_v', np.dtype('float32').char, ('time', 'y', 'x'), zlib=True, least_significant_digit=3)\n",
" outdata.E_var = outdata.ncfile.createVariable('E', np.dtype('float32').char, ('time', 'y', 'x'), zlib=True, least_significant_digit=3)\n",
" \n",
" outdata.time_var = outdata.ncfile.createVariable('time', np.dtype('float32').char, 'time')\n",
" outdata.x_var = outdata.ncfile.createVariable('x', np.dtype('float32').char, 'x')\n",
" outdata.y_var = outdata.ncfile.createVariable('y', np.dtype('float32').char, 'y')\n",
" outdata.rho_var = outdata.ncfile.createVariable('rho', np.dtype('float32').char, ('time', 'y', 'x'), zlib=True,\n",
" least_significant_digit=3)\n",
" outdata.rho_u_var = outdata.ncfile.createVariable('rho_u', np.dtype('float32').char, ('time', 'y', 'x'),\n",
" zlib=True, least_significant_digit=3)\n",
" outdata.rho_v_var = outdata.ncfile.createVariable('rho_v', np.dtype('float32').char, ('time', 'y', 'x'),\n",
" zlib=True, least_significant_digit=3)\n",
" outdata.E_var = outdata.ncfile.createVariable('E', np.dtype('float32').char, ('time', 'y', 'x'), zlib=True,\n",
" least_significant_digit=3)\n",
"\n",
" #Create attributes\n",
" def toJson(in_dict):\n",
" def to_json(in_dict):\n",
" out_dict = in_dict.copy()\n",
"\n",
" for key in out_dict:\n",
@@ -128,25 +132,26 @@
" out_dict[key] = str(out_dict[key])\n",
"\n",
" return json.dumps(out_dict)\n",
"\n",
" outdata.ncfile.created = time.ctime(time.time())\n",
" outdata.ncfile.sim_args = toJson(sim_args)\n",
" outdata.ncfile.sim_args = to_json(sim_args)\n",
"\n",
" outdata.x_var[:] = np.linspace(0, sim.nx*sim.dx, sim.nx)\n",
" outdata.y_var[:] = np.linspace(0, sim.ny*sim.dy, sim.ny)\n",
" outdata.x_var[:] = np.linspace(0, sim.nx * sim.dx, sim.nx)\n",
" outdata.y_var[:] = np.linspace(0, sim.ny * sim.dy, sim.ny)\n",
"\n",
" progress_printer = Common.ProgressPrinter(n_save, print_every=10)\n",
" print(\"Simulating to t={:f}. Estimated {:d} timesteps (dt={:f})\".format(t_end, int(t_end / sim.dt), sim.dt))\n",
" for i in range(n_save+1):\n",
" progress_printer = ProgressPrinter(n_save, print_every=10)\n",
" print(f\"Simulating to t={t_end}. Estimated {int(t_end / sim.dt)} timesteps (dt={sim.dt})\")\n",
" for i in range(n_save + 1):\n",
" #Sanity check simulator\n",
" try:\n",
" sim.check()\n",
" except AssertionError as e:\n",
" print(\"Error after {:d} steps (t={:f}: {:s}\".format(sim.sim_steps(), sim.sim_time(), str(e)))\n",
" print(f\"Error after {sim.sim_steps()} steps (t={sim.sim_time()}: {str(e)}\")\n",
" return outdata.filename\n",
"\n",
" #Simulate\n",
" if (i > 0):\n",
" sim.simulate(t_end/n_save)\n",
" if i > 0:\n",
" sim.simulate(t_end / n_save)\n",
"\n",
" #Save to file\n",
" #rho = sim.u0[0].download(sim.stream)\n",
@@ -159,49 +164,51 @@
"\n",
" #Write progress to screen\n",
" print_string = progress_printer.getPrintString(i)\n",
" if (print_string):\n",
" if print_string:\n",
" print(print_string)\n",
"\n",
" return outdata.filename "
],
"outputs": [],
"execution_count": null
" return outdata.filename"
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"class VisType(Enum):\n",
" Schlieren = 0\n",
" Density = 1\n",
"\n",
"def createAnimation(infile, vis_type, vmin, vmax, save_anim=False, cmap=plt.cm.coolwarm, fig_args=None):\n",
"\n",
"def create_animation(infile, vis_type, vmin, vmax, save_anim=False, cmap=plt.cm.coolwarm, fig_args=None):\n",
" fig = plt.figure(**fig_args)\n",
" \n",
" with Common.DataDumper(infile, 'r') as indata:\n",
"\n",
" with DataDumper(infile, 'r') as indata:\n",
" time = indata.ncfile.variables['time'][:]\n",
" x = indata.ncfile.variables['x'][:]\n",
" y = indata.ncfile.variables['y'][:]\n",
" rho = indata.ncfile.variables['rho'][0]\n",
" rho_u = indata.ncfile.variables['rho_u'][0]\n",
" rho_v = indata.ncfile.variables['rho_v'][0]\n",
" \n",
"\n",
" created = indata.ncfile.created\n",
" sim_args = json.loads(indata.ncfile.sim_args)\n",
" for key in sim_args:\n",
" if isinstance(sim_args[key], list):\n",
" sim_args[key] = \"[...]\"\n",
" num_frames = len(time)\n",
" print(\"{:s} created {:s} contains {:d} timesteps\".format(infile, created, num_frames))\n",
" print(f\"{infile} created {created} contains {num_frames} timesteps\")\n",
" print(\"Simulator arguments: \\n\", sim_args)\n",
"\n",
" ax1 = fig.gca()\n",
" extents = [0, x.max(), 0, y.max()]\n",
" \n",
" if (vis_type == VisType.Schlieren):\n",
" im = ax1.imshow(Visualization.genColors(rho, rho_u, rho_v, cmap, vmax, vmin), origin='lower', extent=extents, cmap='gray', vmin=0.0, vmax=1.0)\n",
"\n",
" if vis_type == VisType.Schlieren:\n",
" im = ax1.imshow(Visualization.genColors(rho, rho_u, rho_v, cmap, vmax, vmin), origin='lower',\n",
" extent=extents, cmap='gray', vmin=0.0, vmax=1.0)\n",
" fig.suptitle(\"Schlieren / vorticity at t={:.2f}\".format(time[0]))\n",
" elif (vis_type == VisType.Density):\n",
" elif vis_type == VisType.Density:\n",
" im = ax1.imshow(rho, origin='lower', extent=extents, cmap=cmap, vmin=vmin, vmax=vmax)\n",
" fig.suptitle(\"Density at t={:.2f}\".format(time[0]))\n",
" else:\n",
@@ -213,32 +220,32 @@
" divider = make_axes_locatable(ax1)\n",
" ax2 = divider.append_axes(\"right\", size=0.1, pad=0.05)\n",
" cb1 = ColorbarBase(ax2, cmap=cmap,\n",
" norm=Normalize(vmin=vmin, vmax=vmax),\n",
" orientation='vertical')\n",
" \n",
" norm=Normalize(vmin=vmin, vmax=vmax),\n",
" orientation='vertical')\n",
"\n",
" #Label colorbar\n",
" if (vis_type == VisType.Schlieren):\n",
" if vis_type == VisType.Schlieren:\n",
" cb1.set_label('Vorticity')\n",
" elif (vis_type == VisType.Density):\n",
" elif vis_type == VisType.Density:\n",
" cb1.set_label('Density')\n",
"\n",
" progress_printer = Common.ProgressPrinter(num_frames, print_every=5)\n",
" \n",
" progress_printer = ProgressPrinter(num_frames, print_every=5)\n",
"\n",
" def animate(i):\n",
" rho = indata.ncfile.variables['rho'][i]\n",
" rho_u = indata.ncfile.variables['rho_u'][i]\n",
" rho_v = indata.ncfile.variables['rho_v'][i]\n",
" \n",
" if (vis_type == VisType.Schlieren):\n",
"\n",
" if vis_type == VisType.Schlieren:\n",
" im.set_data(Visualization.genColors(rho, rho_u, rho_v, cmap, vmax, vmin))\n",
" fig.suptitle(\"Schlieren / vorticity at t={:.2f}\".format(time[i]))\n",
" elif (vis_type == VisType.Density):\n",
" im.set_data(rho) \n",
" elif vis_type == VisType.Density:\n",
" im.set_data(rho)\n",
" fig.suptitle(\"Density at t={:.2f}\".format(time[i]))\n",
" \n",
"\n",
" #Write progress to screen\n",
" print_string = progress_printer.getPrintString(i)\n",
" if (print_string):\n",
" if print_string:\n",
" print(print_string)\n",
"\n",
" anim = animation.FuncAnimation(fig, animate, interval=50, frames=range(num_frames))\n",
@@ -267,9 +274,7 @@
" #plt.rcParams[\"animation.html\"] = \"html5\"\n",
" plt.rcParams[\"animation.html\"] = \"jshtml\"\n",
" display(anim)"
],
"outputs": [],
"execution_count": null
]
},
{
"cell_type": "markdown",
@@ -279,47 +284,48 @@
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"nx = 400\n",
"arguments = InitialConditions.genShockBubble(nx, nx//4, 1.4)\n",
"arguments = InitialConditions.genShockBubble(nx, nx // 4, 1.4)\n",
"plt.figure()\n",
"plt.imshow(Visualization.genSchlieren(arguments['rho']), cmap='gray')\n",
"plt.colorbar(orientation='vertical', aspect=20, pad=0.02, shrink=0.3)"
],
"outputs": [],
"execution_count": null
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"nx = 800#1600\n",
"ny = nx//4\n",
"nx = 800 #1600\n",
"ny = nx // 4\n",
"g = 0.0\n",
"gamma = 1.4\n",
"t_end = 0.5#3.0\n",
"n_save = 20#500\n",
"t_end = 0.5 #3.0\n",
"n_save = 20 #500\n",
"outfile = \"data/euler_shock-bubble.nc\"\n",
"outdata = None\n",
"\n",
"arguments = InitialConditions.genShockBubble(nx, ny, gamma)\n",
"arguments['context'] = my_context\n",
"outfile = runSimulation(outfile, t_end, arguments)"
],
"outputs": [],
"execution_count": null
"outfile = run_simulation(outfile, t_end, arguments)"
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"#outfile = 'data/euler_shock-bubble_0008.nc'\n",
"createAnimation(outfile, vis_type=VisType.Schlieren, vmin=-0.2, vmax=0.2, cmap=plt.cm.RdBu, save_anim=False, fig_args={'figsize':(16, 5)})"
],
"outputs": [],
"execution_count": null
"create_animation(outfile, vis_type=VisType.Schlieren, vmin=-0.2, vmax=0.2, cmap=plt.cm.RdBu, save_anim=False,\n",
" fig_args={'figsize': (16, 5)})"
]
},
{
"cell_type": "markdown",
@@ -329,47 +335,47 @@
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"nx = 400\n",
"arguments = InitialConditions.genKelvinHelmholtz(nx, nx//4, 1.4)\n",
"arguments = InitialConditions.genKelvinHelmholtz(nx, nx // 4, 1.4)\n",
"\n",
"plt.figure()\n",
"plt.imshow(arguments['rho'])\n",
"plt.colorbar(orientation='vertical', aspect=20, pad=0.02, shrink=0.3)"
],
"outputs": [],
"execution_count": null
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"nx = 400\n",
"ny = nx//2\n",
"ny = nx // 2\n",
"roughness = 0.125\n",
"t_end = 10.0\n",
"n_save = 100#1000\n",
"n_save = 100 #1000\n",
"outfile = \"data/euler_kelvin_helmholtz.nc\"\n",
"outdata = None\n",
"\n",
"arguments = InitialConditions.genKelvinHelmholtz(nx, ny, gamma, roughness)\n",
"arguments['context'] = my_context\n",
"outfile = runSimulation(outfile, t_end, arguments)"
],
"outputs": [],
"execution_count": null
"outfile = run_simulation(outfile, t_end, arguments)"
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"#outfile='data/euler_kelvin_helmholtz_0012.nc'\n",
"createAnimation(outfile, vis_type=VisType.Density, vmin=1.0, vmax=2.0, save_anim=False, fig_args={'figsize':(16, 9)})"
],
"outputs": [],
"execution_count": null
"create_animation(outfile, vis_type=VisType.Density, vmin=1.0, vmax=2.0, save_anim=False, fig_args={'figsize': (16, 9)})"
]
},
{
"cell_type": "markdown",
@@ -383,22 +389,24 @@
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"nx = 400\n",
"arguments = InitialConditions.genRayleighTaylor(nx, nx*3, 1.4, version=0)\n",
"plotVars(arguments['rho'], arguments['rho_u'], arguments['rho_v'], arguments['E'])"
],
"outputs": [],
"execution_count": null
"arguments = InitialConditions.genRayleighTaylor(nx, nx * 3, 1.4, version=0)\n",
"plot_vars(arguments['rho'], arguments['rho_u'], arguments['rho_v'], arguments['E'])"
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"nx = 151\n",
"ny = nx*3\n",
"ny = nx * 3\n",
"g = 0.1\n",
"gamma = 1.4\n",
"t_end = 30\n",
@@ -408,27 +416,26 @@
"\n",
"arguments = InitialConditions.genRayleighTaylor(nx, ny, gamma)\n",
"arguments['context'] = my_context\n",
"outfile = runSimulation(outfile, t_end, arguments)"
],
"outputs": [],
"execution_count": null
"outfile = run_simulation(outfile, t_end, arguments)"
]
},
{
"cell_type": "code",
"metadata": {},
"cell_type": "code",
"outputs": [],
"execution_count": null,
"source": [
"#outfile = 'data/euler_rayleigh-taylor_0007.nc'\n",
"createAnimation(outfile, vis_type=VisType.Density, vmin=1, vmax=2, cmap=plt.cm.coolwarm, save_anim=False, fig_args={'figsize':(3.4, 8)})"
],
"outputs": [],
"execution_count": null
"create_animation(outfile, vis_type=VisType.Density, vmin=1, vmax=2, cmap=plt.cm.coolwarm, save_anim=False,\n",
" fig_args={'figsize': (3.4, 8)})"
]
},
{
"cell_type": "code",
"metadata": {},
"source": [],
"cell_type": "code",
"outputs": [],
"execution_count": null
"execution_count": null,
"source": ""
}
],
"metadata": {

129
Figures_Saga.ipynb
View File

@@ -31,7 +31,7 @@
"outputs": [],
"source": [
"def speedup(t_serial, t):\n",
" return t_serial/t"
" return t_serial / t"
]
},
{
@@ -45,10 +45,10 @@
"\n",
" json_filenames = [file for file in os.listdir(profile_dir_path) if file.endswith(\"_profiling.json\")]\n",
"\n",
" if(drop_singlenode):\n",
" if drop_singlenode:\n",
" json_filenames = [file for file in json_filenames if \"1_nodes\" not in file]\n",
"\n",
" if(drop_multinode):\n",
" if drop_multinode:\n",
" json_filenames = [file for file in json_filenames if \"1_nodes\" in file]\n",
"\n",
" for json_filename in json_filenames:\n",
@@ -79,11 +79,11 @@
"evalue": "[Errno 2] No such file or directory: 'output_saga/weak_scaling/2022-06-21T235600/'",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mFileNotFoundError\u001b[0m Traceback (most recent call last)",
"Cell \u001b[0;32mIn[5], line 12\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[38;5;66;03m# DGX-2\u001b[39;00m\n\u001b[1;32m 2\u001b[0m \u001b[38;5;66;03m#weak_scaling_profiling_data = read_profiling_files(\"output_dgx-2/weak_scaling/2022-06-09T134809/\")\u001b[39;00m\n\u001b[1;32m 3\u001b[0m \u001b[38;5;66;03m#weak_scaling_profiling_data = read_profiling_files(\"output_dgx-2/weak_scaling/2022-06-23T154025/\")\u001b[39;00m\n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 10\u001b[0m \u001b[38;5;66;03m#singlenode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_multinode=True)\u001b[39;00m\n\u001b[1;32m 11\u001b[0m \u001b[38;5;66;03m#multinode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_singlenode=True)\u001b[39;00m\n\u001b[0;32m---> 12\u001b[0m singlenode_weak_scaling_profiling_data \u001b[38;5;241m=\u001b[39m \u001b[43mread_profiling_files\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43moutput_saga/weak_scaling/2022-06-21T235600/\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mdrop_multinode\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;28;43;01mTrue\u001b[39;49;00m\u001b[43m)\u001b[49m\n\u001b[1;32m 13\u001b[0m multinode_weak_scaling_profiling_data \u001b[38;5;241m=\u001b[39m read_profiling_files(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124moutput_saga/weak_scaling/2022-06-21T235600/\u001b[39m\u001b[38;5;124m\"\u001b[39m, drop_singlenode\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mTrue\u001b[39;00m)\n\u001b[1;32m 15\u001b[0m \u001b[38;5;28mprint\u001b[39m(singlenode_weak_scaling_profiling_data)\n",
"Cell \u001b[0;32mIn[4], line 4\u001b[0m, in \u001b[0;36mread_profiling_files\u001b[0;34m(profile_dir_path, drop_multinode, drop_singlenode)\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mread_profiling_files\u001b[39m(profile_dir_path\u001b[38;5;241m=\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m.\u001b[39m\u001b[38;5;124m\"\u001b[39m, drop_multinode\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m, drop_singlenode\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m):\n\u001b[1;32m 2\u001b[0m profiling_data \u001b[38;5;241m=\u001b[39m pd\u001b[38;5;241m.\u001b[39mDataFrame()\n\u001b[0;32m----> 4\u001b[0m json_filenames \u001b[38;5;241m=\u001b[39m [file \u001b[38;5;28;01mfor\u001b[39;00m file \u001b[38;5;129;01min\u001b[39;00m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mlistdir\u001b[49m\u001b[43m(\u001b[49m\u001b[43mprofile_dir_path\u001b[49m\u001b[43m)\u001b[49m \u001b[38;5;28;01mif\u001b[39;00m file\u001b[38;5;241m.\u001b[39mendswith(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m_profiling.json\u001b[39m\u001b[38;5;124m\"\u001b[39m)]\n\u001b[1;32m 6\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m(drop_singlenode):\n\u001b[1;32m 7\u001b[0m json_filenames \u001b[38;5;241m=\u001b[39m [file \u001b[38;5;28;01mfor\u001b[39;00m file \u001b[38;5;129;01min\u001b[39;00m json_filenames \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m1_nodes\u001b[39m\u001b[38;5;124m\"\u001b[39m \u001b[38;5;129;01mnot\u001b[39;00m \u001b[38;5;129;01min\u001b[39;00m file]\n",
"\u001b[0;31mFileNotFoundError\u001b[0m: [Errno 2] No such file or directory: 'output_saga/weak_scaling/2022-06-21T235600/'"
"\u001B[0;31m---------------------------------------------------------------------------\u001B[0m",
"\u001B[0;31mFileNotFoundError\u001B[0m Traceback (most recent call last)",
"Cell \u001B[0;32mIn[5], line 12\u001B[0m\n\u001B[1;32m 1\u001B[0m \u001B[38;5;66;03m# DGX-2\u001B[39;00m\n\u001B[1;32m 2\u001B[0m \u001B[38;5;66;03m#weak_scaling_profiling_data = read_profiling_files(\"output_dgx-2/weak_scaling/2022-06-09T134809/\")\u001B[39;00m\n\u001B[1;32m 3\u001B[0m \u001B[38;5;66;03m#weak_scaling_profiling_data = read_profiling_files(\"output_dgx-2/weak_scaling/2022-06-23T154025/\")\u001B[39;00m\n\u001B[0;32m (...)\u001B[0m\n\u001B[1;32m 10\u001B[0m \u001B[38;5;66;03m#singlenode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_multinode=True)\u001B[39;00m\n\u001B[1;32m 11\u001B[0m \u001B[38;5;66;03m#multinode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_singlenode=True)\u001B[39;00m\n\u001B[0;32m---> 12\u001B[0m singlenode_weak_scaling_profiling_data \u001B[38;5;241m=\u001B[39m \u001B[43mread_profiling_files\u001B[49m\u001B[43m(\u001B[49m\u001B[38;5;124;43m\"\u001B[39;49m\u001B[38;5;124;43moutput_saga/weak_scaling/2022-06-21T235600/\u001B[39;49m\u001B[38;5;124;43m\"\u001B[39;49m\u001B[43m,\u001B[49m\u001B[43m \u001B[49m\u001B[43mdrop_multinode\u001B[49m\u001B[38;5;241;43m=\u001B[39;49m\u001B[38;5;28;43;01mTrue\u001B[39;49;00m\u001B[43m)\u001B[49m\n\u001B[1;32m 13\u001B[0m multinode_weak_scaling_profiling_data \u001B[38;5;241m=\u001B[39m read_profiling_files(\u001B[38;5;124m\"\u001B[39m\u001B[38;5;124moutput_saga/weak_scaling/2022-06-21T235600/\u001B[39m\u001B[38;5;124m\"\u001B[39m, drop_singlenode\u001B[38;5;241m=\u001B[39m\u001B[38;5;28;01mTrue\u001B[39;00m)\n\u001B[1;32m 15\u001B[0m \u001B[38;5;28mprint\u001B[39m(singlenode_weak_scaling_profiling_data)\n",
"Cell \u001B[0;32mIn[4], line 4\u001B[0m, in \u001B[0;36mread_profiling_files\u001B[0;34m(profile_dir_path, drop_multinode, drop_singlenode)\u001B[0m\n\u001B[1;32m 1\u001B[0m \u001B[38;5;28;01mdef\u001B[39;00m \u001B[38;5;21mread_profiling_files\u001B[39m(profile_dir_path\u001B[38;5;241m=\u001B[39m\u001B[38;5;124m\"\u001B[39m\u001B[38;5;124m.\u001B[39m\u001B[38;5;124m\"\u001B[39m, drop_multinode\u001B[38;5;241m=\u001B[39m\u001B[38;5;28;01mFalse\u001B[39;00m, drop_singlenode\u001B[38;5;241m=\u001B[39m\u001B[38;5;28;01mFalse\u001B[39;00m):\n\u001B[1;32m 2\u001B[0m profiling_data \u001B[38;5;241m=\u001B[39m pd\u001B[38;5;241m.\u001B[39mDataFrame()\n\u001B[0;32m----> 4\u001B[0m json_filenames \u001B[38;5;241m=\u001B[39m [file \u001B[38;5;28;01mfor\u001B[39;00m file \u001B[38;5;129;01min\u001B[39;00m \u001B[43mos\u001B[49m\u001B[38;5;241;43m.\u001B[39;49m\u001B[43mlistdir\u001B[49m\u001B[43m(\u001B[49m\u001B[43mprofile_dir_path\u001B[49m\u001B[43m)\u001B[49m \u001B[38;5;28;01mif\u001B[39;00m file\u001B[38;5;241m.\u001B[39mendswith(\u001B[38;5;124m\"\u001B[39m\u001B[38;5;124m_profiling.json\u001B[39m\u001B[38;5;124m\"\u001B[39m)]\n\u001B[1;32m 6\u001B[0m \u001B[38;5;28;01mif\u001B[39;00m(drop_singlenode):\n\u001B[1;32m 7\u001B[0m json_filenames \u001B[38;5;241m=\u001B[39m [file \u001B[38;5;28;01mfor\u001B[39;00m file \u001B[38;5;129;01min\u001B[39;00m json_filenames \u001B[38;5;28;01mif\u001B[39;00m \u001B[38;5;124m\"\u001B[39m\u001B[38;5;124m1_nodes\u001B[39m\u001B[38;5;124m\"\u001B[39m \u001B[38;5;129;01mnot\u001B[39;00m \u001B[38;5;129;01min\u001B[39;00m file]\n",
"\u001B[0;31mFileNotFoundError\u001B[0m: [Errno 2] No such file or directory: 'output_saga/weak_scaling/2022-06-21T235600/'"
]
}
],
@@ -99,8 +99,10 @@
"# Saga\n",
"#singlenode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_multinode=True)\n",
"#multinode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_singlenode=True)\n",
"singlenode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-21T235600/\", drop_multinode=True)\n",
"multinode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-21T235600/\", drop_singlenode=True)\n",
"singlenode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-21T235600/\",\n",
" drop_multinode=True)\n",
"multinode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-21T235600/\",\n",
" drop_singlenode=True)\n",
"\n",
"print(singlenode_weak_scaling_profiling_data)\n",
"print(multinode_weak_scaling_profiling_data)"
@@ -217,8 +219,10 @@
"#strong_scaling_profiling_data = read_profiling_files(\"output_hgx/strong_scaling/2022-06-16T152945/\")\n",
"\n",
"# Saga\n",
"singlenode_strong_scaling_profiling_data = read_profiling_files(\"output_saga/strong_scaling/2022-06-16T190721/\", drop_multinode=True)\n",
"multinode_strong_scaling_profiling_data = read_profiling_files(\"output_saga/strong_scaling/2022-06-16T190721/\", drop_singlenode=True)\n",
"singlenode_strong_scaling_profiling_data = read_profiling_files(\"output_saga/strong_scaling/2022-06-16T190721/\",\n",
" drop_multinode=True)\n",
"multinode_strong_scaling_profiling_data = read_profiling_files(\"output_saga/strong_scaling/2022-06-16T190721/\",\n",
" drop_singlenode=True)\n",
"\n",
"print(singlenode_strong_scaling_profiling_data)\n",
"print(multinode_strong_scaling_profiling_data)"
@@ -255,22 +259,23 @@
"plt.rcParams['axes.linewidth'] = 2\n",
"plt.rcParams['lines.linewidth'] = 2\n",
"\n",
"fig, (ax_weak, ax_strong) = plt.subplots(1, 2, figsize=(16,6))\n",
"fig, (ax_weak, ax_strong) = plt.subplots(1, 2, figsize=(16, 6))\n",
"\n",
"t_total_no_init_or_file_io = weak_scaling_profiling_data[\"t_total\"] \\\n",
" -weak_scaling_profiling_data[\"t_init\"] \\\n",
" -weak_scaling_profiling_data[\"t_nc_write\"] \\\n",
" -weak_scaling_profiling_data[\"t_sim_init\"]\n",
" - weak_scaling_profiling_data[\"t_init\"] \\\n",
" - weak_scaling_profiling_data[\"t_nc_write\"] \\\n",
" - weak_scaling_profiling_data[\"t_sim_init\"]\n",
"\n",
"t_total_halo_exchange = weak_scaling_profiling_data[\"t_mpi_halo_exchange_download\"] \\\n",
" +weak_scaling_profiling_data[\"t_mpi_halo_exchange_sendreceive\"] \\\n",
" +weak_scaling_profiling_data[\"t_mpi_halo_exchange_upload\"]\n",
" + weak_scaling_profiling_data[\"t_mpi_halo_exchange_sendreceive\"] \\\n",
" + weak_scaling_profiling_data[\"t_mpi_halo_exchange_upload\"]\n",
"\n",
"#ax_weak.plot(weak_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
"# speedup(t_total_no_init_or_file_io[0], t_total_no_init_or_file_io), label=\"Total (no init or file I/O)\")\n",
"\n",
"ax_weak.plot(weak_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
" speedup(weak_scaling_profiling_data[\"t_full_step\"][0], weak_scaling_profiling_data[\"t_full_step\"]), label=\"Runtime (except init and file I/O)\")\n",
"ax_weak.plot(weak_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"),\n",
" speedup(weak_scaling_profiling_data[\"t_full_step\"][0], weak_scaling_profiling_data[\"t_full_step\"]),\n",
" label=\"Runtime (except init and file I/O)\")\n",
"ax_weak.locator_params(axis=\"x\", nbins=16)\n",
"\n",
"\"\"\"\n",
@@ -308,19 +313,20 @@
"#fig, ax = plt.subplots(figsize=(8,6))\n",
"\n",
"t_total_no_init_or_file_io = strong_scaling_profiling_data[\"t_total\"] \\\n",
" -strong_scaling_profiling_data[\"t_init\"] \\\n",
" -strong_scaling_profiling_data[\"t_nc_write\"] \\\n",
" -strong_scaling_profiling_data[\"t_sim_init\"]\n",
" - strong_scaling_profiling_data[\"t_init\"] \\\n",
" - strong_scaling_profiling_data[\"t_nc_write\"] \\\n",
" - strong_scaling_profiling_data[\"t_sim_init\"]\n",
"\n",
"t_total_halo_exchange = strong_scaling_profiling_data[\"t_mpi_halo_exchange_download\"] \\\n",
" +strong_scaling_profiling_data[\"t_mpi_halo_exchange_sendreceive\"] \\\n",
" +strong_scaling_profiling_data[\"t_mpi_halo_exchange_upload\"]\n",
" + strong_scaling_profiling_data[\"t_mpi_halo_exchange_sendreceive\"] \\\n",
" + strong_scaling_profiling_data[\"t_mpi_halo_exchange_upload\"]\n",
"\n",
"#ax_strong.plot(strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
"# speedup(t_total_no_init_or_file_io[0], t_total_no_init_or_file_io)*4, label=\"Total (no init or file I/O)\")\n",
"\n",
"ax_strong.plot(strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
" speedup(strong_scaling_profiling_data[\"t_full_step\"][0], strong_scaling_profiling_data[\"t_full_step\"])*4, label=\"Runtime (except init and file I/O)\")\n",
"ax_strong.plot(strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"),\n",
" speedup(strong_scaling_profiling_data[\"t_full_step\"][0],\n",
" strong_scaling_profiling_data[\"t_full_step\"]) * 4, label=\"Runtime (except init and file I/O)\")\n",
"\n",
"\"\"\"\n",
"ax_strong.plot(strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
@@ -398,21 +404,24 @@
"plt.rcParams['axes.linewidth'] = 2\n",
"plt.rcParams['lines.linewidth'] = 2\n",
"\n",
"fig, (ax_weak, ax_strong) = plt.subplots(1, 2, figsize=(16,6))\n",
"fig, (ax_weak, ax_strong) = plt.subplots(1, 2, figsize=(16, 6))\n",
"\n",
"ax_weak.plot(singlenode_weak_scaling_profiling_data[\"n_processes\"][0:4].to_numpy(dtype=\"int\"), \n",
" speedup(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0], singlenode_weak_scaling_profiling_data[\"t_full_step\"][0:4]), \n",
" label=\"Shared-memory run time\", marker=\"x\")\n",
"ax_weak.plot(singlenode_weak_scaling_profiling_data[\"n_processes\"][0:4].to_numpy(dtype=\"int\"),\n",
" speedup(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0],\n",
" singlenode_weak_scaling_profiling_data[\"t_full_step\"][0:4]),\n",
" label=\"Shared-memory run time\", marker=\"x\")\n",
"\n",
"ax_weak.plot(singlenode_weak_scaling_profiling_data[\"n_processes\"][0:4].to_numpy(dtype=\"int\"), \n",
" speedup(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0], \n",
" np.append(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0], multinode_weak_scaling_profiling_data[\"t_full_step\"][0:3])), \n",
" label=\"Distributed-memory (1-4 nodes) run time\", marker=\"o\", color=\"green\")\n",
"ax_weak.plot(singlenode_weak_scaling_profiling_data[\"n_processes\"][0:4].to_numpy(dtype=\"int\"),\n",
" speedup(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0],\n",
" np.append(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0],\n",
" multinode_weak_scaling_profiling_data[\"t_full_step\"][0:3])),\n",
" label=\"Distributed-memory (1-4 nodes) run time\", marker=\"o\", color=\"green\")\n",
"\n",
"ax_weak.locator_params(axis=\"x\", nbins=4)\n",
"ax_weak.set_ylim(0.89, 1.01)\n",
"\n",
"ax_weak.plot(nproc[0:4], np.ones(len(nproc[0:4])), label=\"Ideal run time (constant)\", linestyle=\"dotted\", color=\"orange\")\n",
"ax_weak.plot(nproc[0:4], np.ones(len(nproc[0:4])), label=\"Ideal run time (constant)\", linestyle=\"dotted\",\n",
" color=\"orange\")\n",
"\n",
"ax_weak.set_xlabel(\"Number of ranks/GPUs\")\n",
"ax_weak.set_ylabel(\"Efficiency\")\n",
@@ -425,9 +434,10 @@
"# speedup(singlenode_strong_scaling_profiling_data[\"t_full_step\"][0], singlenode_strong_scaling_profiling_data[\"t_full_step\"])*4, \n",
"# label=\"Single-node (no init or file I/O)\", marker=\"x\")\n",
"\n",
"ax_strong.plot(multinode_strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
" speedup(multinode_strong_scaling_profiling_data[\"t_full_step\"][0], multinode_strong_scaling_profiling_data[\"t_full_step\"])*4, \n",
" label=\"Distributed-memory (4 nodes) run time\", marker=\"o\")\n",
"ax_strong.plot(multinode_strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"),\n",
" speedup(multinode_strong_scaling_profiling_data[\"t_full_step\"][0],\n",
" multinode_strong_scaling_profiling_data[\"t_full_step\"]) * 4,\n",
" label=\"Distributed-memory (4 nodes) run time\", marker=\"o\")\n",
"\n",
"ax_strong.locator_params(axis=\"x\", nbins=16)\n",
"\n",
@@ -464,18 +474,19 @@
"source": [
"bar_width = 0.8\n",
"\n",
"share_init = (profiling_data[\"t_init\"]+profiling_data[\"t_sim_init\"])/profiling_data[\"t_total\"]\n",
"share_mpi = profiling_data[\"t_step_mpi_halo_exchange\"]/profiling_data[\"t_total\"]\n",
"share_simulate = profiling_data[\"t_step_mpi\"]/profiling_data[\"t_total\"]\n",
"share_init = (profiling_data[\"t_init\"] + profiling_data[\"t_sim_init\"]) / profiling_data[\"t_total\"]\n",
"share_mpi = profiling_data[\"t_step_mpi_halo_exchange\"] / profiling_data[\"t_total\"]\n",
"share_simulate = profiling_data[\"t_step_mpi\"] / profiling_data[\"t_total\"]\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_init, color='red', edgecolor='white', width=bar_width, label=\"Init\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_init, color='red', edgecolor='white', width=bar_width, label=\"Init\")\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_mpi, bottom=share_init, color='orange', edgecolor='white', width=bar_width, label=\"MPI\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_mpi, bottom=share_init, color='orange', edgecolor='white', width=bar_width, label=\"MPI\")\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_simulate, bottom=share_init+share_mpi , color='green', edgecolor='white', width=bar_width, label=\"Simulate\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_simulate, bottom=share_init + share_mpi, color='green', edgecolor='white', width=bar_width,\n",
" label=\"Simulate\")\n",
"\n",
"plt.xlabel(\"Number of cores/GPUs\")\n",
"plt.ylabel(\"Share of normalized runtime\")\n",
@@ -505,18 +516,19 @@
"source": [
"bar_width = 0.8\n",
"\n",
"share_init = (profiling_data[\"t_init\"]+profiling_data[\"t_sim_init\"])\n",
"share_init = (profiling_data[\"t_init\"] + profiling_data[\"t_sim_init\"])\n",
"share_mpi = profiling_data[\"t_step_mpi_halo_exchange\"]\n",
"share_simulate = profiling_data[\"t_step_mpi\"]\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_init, color='red', edgecolor='white', width=bar_width, label=\"Init\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_init, color='red', edgecolor='white', width=bar_width, label=\"Init\")\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_mpi, bottom=share_init, color='orange', edgecolor='white', width=bar_width, label=\"MPI\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_mpi, bottom=share_init, color='orange', edgecolor='white', width=bar_width, label=\"MPI\")\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_simulate, bottom=share_init+share_mpi , color='green', edgecolor='white', width=bar_width, label=\"Simulate\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_simulate, bottom=share_init + share_mpi, color='green', edgecolor='white', width=bar_width,\n",
" label=\"Simulate\")\n",
"\n",
"plt.xlabel(\"Number of cores/GPUs\")\n",
"plt.ylabel(\"Runtime (s)\")\n",
@@ -546,13 +558,12 @@
"source": [
"bar_width = 0.8\n",
"\n",
"share_init = (profiling_data[\"t_init\"]+profiling_data[\"t_sim_init\"])\n",
"share_init = (profiling_data[\"t_init\"] + profiling_data[\"t_sim_init\"])\n",
"share_mpi = profiling_data[\"t_step_mpi_halo_exchange\"]\n",
"share_simulate = profiling_data[\"t_step\"]\n",
"\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_simulate, color='green', edgecolor='white', width=bar_width, label=\"Simulate\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_simulate, color='green', edgecolor='white', width=bar_width, label=\"Simulate\")\n",
"\n",
"plt.xlabel(\"Number of cores/GPUs\")\n",
"plt.ylabel(\"Runtime (s)\")\n",

105
Figures_Simula.ipynb
View File

@@ -31,7 +31,7 @@
"outputs": [],
"source": [
"def speedup(t_serial, t):\n",
" return t_serial/t"
" return t_serial / t"
]
},
{
@@ -45,10 +45,10 @@
"\n",
" json_filenames = [file for file in os.listdir(profile_dir_path) if file.endswith(\"_profiling.json\")]\n",
"\n",
" if(drop_singlenode):\n",
" if drop_singlenode:\n",
" json_filenames = [file for file in json_filenames if \"1_nodes\" not in file]\n",
"\n",
" if(drop_multinode):\n",
" if drop_multinode:\n",
" json_filenames = [file for file in json_filenames if \"1_nodes\" in file]\n",
"\n",
" for json_filename in json_filenames:\n",
@@ -79,11 +79,11 @@
"evalue": "[Errno 2] No such file or directory: 'output_dgx-2/weak_scaling/2022-06-23T154025/'",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mFileNotFoundError\u001b[0m Traceback (most recent call last)",
"Cell \u001b[0;32mIn[11], line 3\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[38;5;66;03m# DGX-2\u001b[39;00m\n\u001b[1;32m 2\u001b[0m \u001b[38;5;66;03m#weak_scaling_profiling_data = read_profiling_files(\"output_dgx-2/weak_scaling/2022-06-09T134809/\")\u001b[39;00m\n\u001b[0;32m----> 3\u001b[0m weak_scaling_profiling_data \u001b[38;5;241m=\u001b[39m \u001b[43mread_profiling_files\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43moutput_dgx-2/weak_scaling/2022-06-23T154025/\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\n\u001b[1;32m 5\u001b[0m \u001b[38;5;66;03m# HGX\u001b[39;00m\n\u001b[1;32m 6\u001b[0m \u001b[38;5;66;03m#weak_scaling_profiling_data = read_profiling_files(\"output_hgx/weak_scaling/2022-06-16T162931/\")\u001b[39;00m\n\u001b[1;32m 7\u001b[0m \u001b[38;5;66;03m##weak_scaling_profiling_data = read_profiling_files(\"output_hgx/weak_scaling/2022-06-16T170630/\")\u001b[39;00m\n\u001b[0;32m (...)\u001b[0m\n\u001b[1;32m 10\u001b[0m \u001b[38;5;66;03m#singlenode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_multinode=True)\u001b[39;00m\n\u001b[1;32m 11\u001b[0m \u001b[38;5;66;03m#multinode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_singlenode=True)\u001b[39;00m\n\u001b[1;32m 13\u001b[0m \u001b[38;5;28mprint\u001b[39m(weak_scaling_profiling_data)\n",
"Cell \u001b[0;32mIn[10], line 4\u001b[0m, in \u001b[0;36mread_profiling_files\u001b[0;34m(profile_dir_path, drop_multinode, drop_singlenode)\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mread_profiling_files\u001b[39m(profile_dir_path\u001b[38;5;241m=\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m.\u001b[39m\u001b[38;5;124m\"\u001b[39m, drop_multinode\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m, drop_singlenode\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mFalse\u001b[39;00m):\n\u001b[1;32m 2\u001b[0m profiling_data \u001b[38;5;241m=\u001b[39m pd\u001b[38;5;241m.\u001b[39mDataFrame()\n\u001b[0;32m----> 4\u001b[0m json_filenames \u001b[38;5;241m=\u001b[39m [file \u001b[38;5;28;01mfor\u001b[39;00m file \u001b[38;5;129;01min\u001b[39;00m \u001b[43mos\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mlistdir\u001b[49m\u001b[43m(\u001b[49m\u001b[43mprofile_dir_path\u001b[49m\u001b[43m)\u001b[49m \u001b[38;5;28;01mif\u001b[39;00m file\u001b[38;5;241m.\u001b[39mendswith(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m_profiling.json\u001b[39m\u001b[38;5;124m\"\u001b[39m)]\n\u001b[1;32m 6\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m(drop_singlenode):\n\u001b[1;32m 7\u001b[0m json_filenames \u001b[38;5;241m=\u001b[39m [file \u001b[38;5;28;01mfor\u001b[39;00m file \u001b[38;5;129;01min\u001b[39;00m json_filenames \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m1_nodes\u001b[39m\u001b[38;5;124m\"\u001b[39m \u001b[38;5;129;01mnot\u001b[39;00m \u001b[38;5;129;01min\u001b[39;00m file]\n",
"\u001b[0;31mFileNotFoundError\u001b[0m: [Errno 2] No such file or directory: 'output_dgx-2/weak_scaling/2022-06-23T154025/'"
"\u001B[0;31m---------------------------------------------------------------------------\u001B[0m",
"\u001B[0;31mFileNotFoundError\u001B[0m Traceback (most recent call last)",
"Cell \u001B[0;32mIn[11], line 3\u001B[0m\n\u001B[1;32m 1\u001B[0m \u001B[38;5;66;03m# DGX-2\u001B[39;00m\n\u001B[1;32m 2\u001B[0m \u001B[38;5;66;03m#weak_scaling_profiling_data = read_profiling_files(\"output_dgx-2/weak_scaling/2022-06-09T134809/\")\u001B[39;00m\n\u001B[0;32m----> 3\u001B[0m weak_scaling_profiling_data \u001B[38;5;241m=\u001B[39m \u001B[43mread_profiling_files\u001B[49m\u001B[43m(\u001B[49m\u001B[38;5;124;43m\"\u001B[39;49m\u001B[38;5;124;43moutput_dgx-2/weak_scaling/2022-06-23T154025/\u001B[39;49m\u001B[38;5;124;43m\"\u001B[39;49m\u001B[43m)\u001B[49m\n\u001B[1;32m 5\u001B[0m \u001B[38;5;66;03m# HGX\u001B[39;00m\n\u001B[1;32m 6\u001B[0m \u001B[38;5;66;03m#weak_scaling_profiling_data = read_profiling_files(\"output_hgx/weak_scaling/2022-06-16T162931/\")\u001B[39;00m\n\u001B[1;32m 7\u001B[0m \u001B[38;5;66;03m##weak_scaling_profiling_data = read_profiling_files(\"output_hgx/weak_scaling/2022-06-16T170630/\")\u001B[39;00m\n\u001B[0;32m (...)\u001B[0m\n\u001B[1;32m 10\u001B[0m \u001B[38;5;66;03m#singlenode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_multinode=True)\u001B[39;00m\n\u001B[1;32m 11\u001B[0m \u001B[38;5;66;03m#multinode_weak_scaling_profiling_data = read_profiling_files(\"output_saga/weak_scaling/2022-06-16T151516/\", drop_singlenode=True)\u001B[39;00m\n\u001B[1;32m 13\u001B[0m \u001B[38;5;28mprint\u001B[39m(weak_scaling_profiling_data)\n",
"Cell \u001B[0;32mIn[10], line 4\u001B[0m, in \u001B[0;36mread_profiling_files\u001B[0;34m(profile_dir_path, drop_multinode, drop_singlenode)\u001B[0m\n\u001B[1;32m 1\u001B[0m \u001B[38;5;28;01mdef\u001B[39;00m \u001B[38;5;21mread_profiling_files\u001B[39m(profile_dir_path\u001B[38;5;241m=\u001B[39m\u001B[38;5;124m\"\u001B[39m\u001B[38;5;124m.\u001B[39m\u001B[38;5;124m\"\u001B[39m, drop_multinode\u001B[38;5;241m=\u001B[39m\u001B[38;5;28;01mFalse\u001B[39;00m, drop_singlenode\u001B[38;5;241m=\u001B[39m\u001B[38;5;28;01mFalse\u001B[39;00m):\n\u001B[1;32m 2\u001B[0m profiling_data \u001B[38;5;241m=\u001B[39m pd\u001B[38;5;241m.\u001B[39mDataFrame()\n\u001B[0;32m----> 4\u001B[0m json_filenames \u001B[38;5;241m=\u001B[39m [file \u001B[38;5;28;01mfor\u001B[39;00m file \u001B[38;5;129;01min\u001B[39;00m \u001B[43mos\u001B[49m\u001B[38;5;241;43m.\u001B[39;49m\u001B[43mlistdir\u001B[49m\u001B[43m(\u001B[49m\u001B[43mprofile_dir_path\u001B[49m\u001B[43m)\u001B[49m \u001B[38;5;28;01mif\u001B[39;00m file\u001B[38;5;241m.\u001B[39mendswith(\u001B[38;5;124m\"\u001B[39m\u001B[38;5;124m_profiling.json\u001B[39m\u001B[38;5;124m\"\u001B[39m)]\n\u001B[1;32m 6\u001B[0m \u001B[38;5;28;01mif\u001B[39;00m(drop_singlenode):\n\u001B[1;32m 7\u001B[0m json_filenames \u001B[38;5;241m=\u001B[39m [file \u001B[38;5;28;01mfor\u001B[39;00m file \u001B[38;5;129;01min\u001B[39;00m json_filenames \u001B[38;5;28;01mif\u001B[39;00m \u001B[38;5;124m\"\u001B[39m\u001B[38;5;124m1_nodes\u001B[39m\u001B[38;5;124m\"\u001B[39m \u001B[38;5;129;01mnot\u001B[39;00m \u001B[38;5;129;01min\u001B[39;00m file]\n",
"\u001B[0;31mFileNotFoundError\u001B[0m: [Errno 2] No such file or directory: 'output_dgx-2/weak_scaling/2022-06-23T154025/'"
]
}
],
@@ -298,22 +298,23 @@
"plt.rcParams['axes.linewidth'] = 2\n",
"plt.rcParams['lines.linewidth'] = 2\n",
"\n",
"fig, (ax_weak, ax_strong) = plt.subplots(1, 2, figsize=(16,6))\n",
"fig, (ax_weak, ax_strong) = plt.subplots(1, 2, figsize=(16, 6))\n",
"\n",
"t_total_no_init_or_file_io = weak_scaling_profiling_data[\"t_total\"] \\\n",
" -weak_scaling_profiling_data[\"t_init\"] \\\n",
" -weak_scaling_profiling_data[\"t_nc_write\"] \\\n",
" -weak_scaling_profiling_data[\"t_sim_init\"]\n",
" - weak_scaling_profiling_data[\"t_init\"] \\\n",
" - weak_scaling_profiling_data[\"t_nc_write\"] \\\n",
" - weak_scaling_profiling_data[\"t_sim_init\"]\n",
"\n",
"t_total_halo_exchange = weak_scaling_profiling_data[\"t_mpi_halo_exchange_download\"] \\\n",
" +weak_scaling_profiling_data[\"t_mpi_halo_exchange_sendreceive\"] \\\n",
" +weak_scaling_profiling_data[\"t_mpi_halo_exchange_upload\"]\n",
" + weak_scaling_profiling_data[\"t_mpi_halo_exchange_sendreceive\"] \\\n",
" + weak_scaling_profiling_data[\"t_mpi_halo_exchange_upload\"]\n",
"\n",
"#ax_weak.plot(weak_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
"# speedup(t_total_no_init_or_file_io[0], t_total_no_init_or_file_io), label=\"Total (no init or file I/O)\")\n",
"\n",
"ax_weak.plot(weak_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
" speedup(weak_scaling_profiling_data[\"t_full_step\"][0], weak_scaling_profiling_data[\"t_full_step\"]), label=\"Shared-memory run time\", marker=\"o\")\n",
"ax_weak.plot(weak_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"),\n",
" speedup(weak_scaling_profiling_data[\"t_full_step\"][0], weak_scaling_profiling_data[\"t_full_step\"]),\n",
" label=\"Shared-memory run time\", marker=\"o\")\n",
"ax_weak.locator_params(axis=\"x\", nbins=16)\n",
"ax_weak.set_ylim(0.89, 1.01)\n",
"\n",
@@ -353,19 +354,20 @@
"#fig, ax = plt.subplots(figsize=(8,6))\n",
"\n",
"t_total_no_init_or_file_io = strong_scaling_profiling_data[\"t_total\"] \\\n",
" -strong_scaling_profiling_data[\"t_init\"] \\\n",
" -strong_scaling_profiling_data[\"t_nc_write\"] \\\n",
" -strong_scaling_profiling_data[\"t_sim_init\"]\n",
" - strong_scaling_profiling_data[\"t_init\"] \\\n",
" - strong_scaling_profiling_data[\"t_nc_write\"] \\\n",
" - strong_scaling_profiling_data[\"t_sim_init\"]\n",
"\n",
"t_total_halo_exchange = strong_scaling_profiling_data[\"t_mpi_halo_exchange_download\"] \\\n",
" +strong_scaling_profiling_data[\"t_mpi_halo_exchange_sendreceive\"] \\\n",
" +strong_scaling_profiling_data[\"t_mpi_halo_exchange_upload\"]\n",
" + strong_scaling_profiling_data[\"t_mpi_halo_exchange_sendreceive\"] \\\n",
" + strong_scaling_profiling_data[\"t_mpi_halo_exchange_upload\"]\n",
"\n",
"#ax_strong.plot(strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
"# speedup(t_total_no_init_or_file_io[0], t_total_no_init_or_file_io)*4, label=\"Total (no init or file I/O)\")\n",
"\n",
"ax_strong.plot(strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
" speedup(strong_scaling_profiling_data[\"t_full_step\"][0], strong_scaling_profiling_data[\"t_full_step\"])*4, label=\"Shared-memory run time\", marker=\"o\")\n",
"ax_strong.plot(strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"),\n",
" speedup(strong_scaling_profiling_data[\"t_full_step\"][0],\n",
" strong_scaling_profiling_data[\"t_full_step\"]) * 4, label=\"Shared-memory run time\", marker=\"o\")\n",
"\n",
"\"\"\"\n",
"ax_strong.plot(strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
@@ -445,15 +447,17 @@
"plt.rcParams['axes.linewidth'] = 2\n",
"plt.rcParams['lines.linewidth'] = 2\n",
"\n",
"fig, (ax_weak, ax_strong) = plt.subplots(1, 2, figsize=(16,6))\n",
"fig, (ax_weak, ax_strong) = plt.subplots(1, 2, figsize=(16, 6))\n",
"\n",
"ax_weak.plot(singlenode_weak_scaling_profiling_data[\"n_processes\"][0:4].to_numpy(dtype=\"int\"), \n",
" speedup(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0], singlenode_weak_scaling_profiling_data[\"t_full_step\"][0:4]), \n",
" label=\"Single-node runtime (no init or file I/O)\", marker=\"x\")\n",
"ax_weak.plot(singlenode_weak_scaling_profiling_data[\"n_processes\"][0:4].to_numpy(dtype=\"int\"),\n",
" speedup(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0],\n",
" singlenode_weak_scaling_profiling_data[\"t_full_step\"][0:4]),\n",
" label=\"Single-node runtime (no init or file I/O)\", marker=\"x\")\n",
"\n",
"ax_weak.plot(multinode_weak_scaling_profiling_data[\"n_processes\"][0:3].to_numpy(dtype=\"int\"), \n",
" speedup(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0], multinode_weak_scaling_profiling_data[\"t_full_step\"][0:3]), \n",
" label=\"2–4 nodes runtime (no init or file I/O)\", marker=\"o\", color=\"green\")\n",
"ax_weak.plot(multinode_weak_scaling_profiling_data[\"n_processes\"][0:3].to_numpy(dtype=\"int\"),\n",
" speedup(singlenode_weak_scaling_profiling_data[\"t_full_step\"][0],\n",
" multinode_weak_scaling_profiling_data[\"t_full_step\"][0:3]),\n",
" label=\"2–4 nodes runtime (no init or file I/O)\", marker=\"o\", color=\"green\")\n",
"\n",
"ax_weak.locator_params(axis=\"x\", nbins=4)\n",
"\n",
@@ -469,9 +473,10 @@
"# speedup(singlenode_strong_scaling_profiling_data[\"t_full_step\"][0], singlenode_strong_scaling_profiling_data[\"t_full_step\"])*4, \n",
"# label=\"Single-node (no init or file I/O)\", marker=\"x\")\n",
"\n",
"ax_strong.plot(multinode_strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"), \n",
" speedup(multinode_strong_scaling_profiling_data[\"t_full_step\"][0], multinode_strong_scaling_profiling_data[\"t_full_step\"])*4, \n",
" label=\"Four nodes runtime (no init or file I/O)\", marker=\"o\")\n",
"ax_strong.plot(multinode_strong_scaling_profiling_data[\"n_processes\"].to_numpy(dtype=\"int\"),\n",
" speedup(multinode_strong_scaling_profiling_data[\"t_full_step\"][0],\n",
" multinode_strong_scaling_profiling_data[\"t_full_step\"]) * 4,\n",
" label=\"Four nodes runtime (no init or file I/O)\", marker=\"o\")\n",
"\n",
"ax_strong.locator_params(axis=\"x\", nbins=16)\n",
"\n",
@@ -506,18 +511,19 @@
"source": [
"bar_width = 0.8\n",
"\n",
"share_init = (profiling_data[\"t_init\"]+profiling_data[\"t_sim_init\"])/profiling_data[\"t_total\"]\n",
"share_mpi = profiling_data[\"t_step_mpi_halo_exchange\"]/profiling_data[\"t_total\"]\n",
"share_simulate = profiling_data[\"t_step_mpi\"]/profiling_data[\"t_total\"]\n",
"share_init = (profiling_data[\"t_init\"] + profiling_data[\"t_sim_init\"]) / profiling_data[\"t_total\"]\n",
"share_mpi = profiling_data[\"t_step_mpi_halo_exchange\"] / profiling_data[\"t_total\"]\n",
"share_simulate = profiling_data[\"t_step_mpi\"] / profiling_data[\"t_total\"]\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_init, color='red', edgecolor='white', width=bar_width, label=\"Init\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_init, color='red', edgecolor='white', width=bar_width, label=\"Init\")\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_mpi, bottom=share_init, color='orange', edgecolor='white', width=bar_width, label=\"MPI\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_mpi, bottom=share_init, color='orange', edgecolor='white', width=bar_width, label=\"MPI\")\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_simulate, bottom=share_init+share_mpi , color='green', edgecolor='white', width=bar_width, label=\"Simulate\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_simulate, bottom=share_init + share_mpi, color='green', edgecolor='white', width=bar_width,\n",
" label=\"Simulate\")\n",
"\n",
"plt.xlabel(\"Number of cores/GPUs\")\n",
"plt.ylabel(\"Share of normalized runtime\")\n",
@@ -547,18 +553,19 @@
"source": [
"bar_width = 0.8\n",
"\n",
"share_init = (profiling_data[\"t_init\"]+profiling_data[\"t_sim_init\"])\n",
"share_init = (profiling_data[\"t_init\"] + profiling_data[\"t_sim_init\"])\n",
"share_mpi = profiling_data[\"t_step_mpi_halo_exchange\"]\n",
"share_simulate = profiling_data[\"t_step_mpi\"]\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_init, color='red', edgecolor='white', width=bar_width, label=\"Init\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_init, color='red', edgecolor='white', width=bar_width, label=\"Init\")\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_mpi, bottom=share_init, color='orange', edgecolor='white', width=bar_width, label=\"MPI\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_mpi, bottom=share_init, color='orange', edgecolor='white', width=bar_width, label=\"MPI\")\n",
"\n",
"plt.bar(profiling_data[\"n_processes\"], \n",
" share_simulate, bottom=share_init+share_mpi , color='green', edgecolor='white', width=bar_width, label=\"Simulate\")\n",
"plt.bar(profiling_data[\"n_processes\"],\n",
" share_simulate, bottom=share_init + share_mpi, color='green', edgecolor='white', width=bar_width,\n",
" label=\"Simulate\")\n",
"\n",
"plt.xlabel(\"Number of cores/GPUs\")\n",
"plt.ylabel(\"Runtime (s)\")\n",

650
MPITest.ipynb
View File

File diff suppressed because one or more lines are too long

88
TestSchemes.ipynb
View File

File diff suppressed because one or more lines are too long

58
mpiTesting.py
View File

@@ -19,7 +19,6 @@ You should have received a copy of the GNU General Public License
along with this program. If not, see <http://www.gnu.org/licenses/>.
"""
import numpy as np
import gc
import time
@@ -35,31 +34,29 @@ import pycuda.driver as cuda
# Simulator engine etc
from GPUSimulators import MPISimulator
from GPUSimulators.common import common
from GPUSimulators.gpu import cuda_context
from GPUSimulators.common import run_simulation, get_git_hash, get_git_status
from GPUSimulators.gpu import CudaContext
from GPUSimulators import EE2D_KP07_dimsplit
from GPUSimulators.helpers import InitialConditions as IC
import argparse
parser = argparse.ArgumentParser(description='Strong and weak scaling experiments.')
parser.add_argument('-nx', type=int, default=128)
parser.add_argument('-ny', type=int, default=128)
parser.add_argument('--profile', action='store_true') # default: False
parser.add_argument('--profile', action='store_true') # default: False
args = parser.parse_args()
if(args.profile):
if args.profile:
profiling_data = {}
# profiling: total run time
t_total_start = time.time()
t_init_start = time.time()
# Get MPI COMM to use
comm = MPI.COMM_WORLD
####
# Initialize logging
####
@@ -72,8 +69,7 @@ logger.setLevel(min(log_level_console, log_level_file))
ch = logging.StreamHandler()
ch.setLevel(log_level_console)
logger.addHandler(ch)
logger.info("Console logger using level %s",
logging.getLevelName(log_level_console))
logger.info(f"Console logger using level {logging.getLevelName(log_level_console)}")
fh = logging.FileHandler(log_filename)
formatter = logging.Formatter(
@@ -81,9 +77,7 @@ formatter = logging.Formatter(
fh.setFormatter(formatter)
fh.setLevel(log_level_file)
logger.addHandler(fh)
logger.info("File logger using level %s to %s",
logging.getLevelName(log_level_file), log_filename)
logger.info(f"File logger using level {logging.getLevelName(log_level_file)} to {log_filename}")
####
# Initialize MPI grid etc
@@ -91,7 +85,6 @@ logger.info("File logger using level %s to %s",
logger.info("Creating MPI grid")
grid = MPISimulator.MPIGrid(MPI.COMM_WORLD)
####
# Initialize CUDA
####
@@ -100,9 +93,8 @@ logger.info("Initializing CUDA")
local_rank = grid.get_local_rank()
num_cuda_devices = cuda.Device.count()
cuda_device = local_rank % num_cuda_devices
logger.info("Process %s using CUDA device %s", str(local_rank), str(cuda_device))
cuda_context = CudaContext.CudaContext(device=cuda_device, autotuning=False)
logger.info(f"Process {str(local_rank)} using CUDA device {str(cuda_device)}")
cuda_context = CudaContext(device=cuda_device, autotuning=False)
####
# Set initial conditions
@@ -118,9 +110,9 @@ ny = args.ny
dt = 0.000001
gamma = 1.4
#save_times = np.linspace(0, 0.000009, 2)
#save_times = np.linspace(0, 0.000099, 11)
#save_times = np.linspace(0, 0.000099, 2)
# save_times = np.linspace(0, 0.000009, 2)
# save_times = np.linspace(0, 0.000099, 11)
# save_times = np.linspace(0, 0.000099, 2)
save_times = np.linspace(0, 0.0000999, 2)
outfile = "mpi_out_" + str(MPI.COMM_WORLD.rank) + ".nc"
save_var_names = ['rho', 'rho_u', 'rho_v', 'E']
@@ -130,7 +122,7 @@ arguments['context'] = cuda_context
arguments['theta'] = 1.2
arguments['grid'] = grid
if(args.profile):
if args.profile:
t_init_end = time.time()
t_init = t_init_end - t_init_start
profiling_data["t_init"] = t_init
@@ -139,6 +131,8 @@ if(args.profile):
# Run simulation
####
logger.info("Running simulation")
# Helper function to create MPI simulator
@@ -148,27 +142,29 @@ def genSim(grid, **kwargs):
return sim
outfile, sim_runner_profiling_data, sim_profiling_data = Common.run_simulation(
outfile, sim_runner_profiling_data, sim_profiling_data = run_simulation(
genSim, arguments, outfile, save_times, save_var_names, dt)
if(args.profile):
if args.profile:
t_total_end = time.time()
t_total = t_total_end - t_total_start
profiling_data["t_total"] = t_total
print("Total run time on rank " + str(MPI.COMM_WORLD.rank) + " is " + str(t_total) + " s")
print(f"Total run time on rank {str(MPI.COMM_WORLD.rank)} is {str(t_total)} s")
# write profiling to json file
if(args.profile and MPI.COMM_WORLD.rank == 0):
# write profiling to JSON file
if args.profile and MPI.COMM_WORLD.rank == 0:
job_id = ""
if "SLURM_JOB_ID" in os.environ:
job_id = int(os.environ["SLURM_JOB_ID"])
allocated_nodes = int(os.environ["SLURM_JOB_NUM_NODES"])
allocated_gpus = int(os.environ["CUDA_VISIBLE_DEVICES"].count(",") + 1)
profiling_file = "MPI_jobid_" + \
str(job_id) + "_" + str(allocated_nodes) + "_nodes_and_" + str(allocated_gpus) + "_GPUs_profiling.json"
str(job_id) + "_" + str(allocated_nodes) + "_nodes_and_" + str(
allocated_gpus) + "_GPUs_profiling.json"
profiling_data["outfile"] = outfile
else:
profiling_file = "MPI_" + str(MPI.COMM_WORLD.size) + "_procs_and_" + str(num_cuda_devices) + "_GPUs_profiling.json"
profiling_file = "MPI_" + str(MPI.COMM_WORLD.size) + "_procs_and_" + str(
num_cuda_devices) + "_GPUs_profiling.json"
for stage in sim_runner_profiling_data["start"].keys():
profiling_data[stage] = sim_runner_profiling_data["end"][stage] - sim_runner_profiling_data["start"][stage]
@@ -184,8 +180,8 @@ if(args.profile and MPI.COMM_WORLD.rank == 0):
profiling_data["slurm_job_id"] = job_id
profiling_data["n_cuda_devices"] = str(num_cuda_devices)
profiling_data["n_processes"] = str(MPI.COMM_WORLD.size)
profiling_data["git_hash"] = Common.get_git_hash()
profiling_data["git_status"] = Common.get_git_status()
profiling_data["git_hash"] = get_git_hash()
profiling_data["git_status"] = get_git_status()
with open(profiling_file, "w") as write_file:
json.dump(profiling_data, write_file)
@@ -200,8 +196,6 @@ arguments = None
logging.shutdown()
gc.collect()
####
# Print completion and exit
####

37
shmemTesting.py
View File

@@ -19,22 +19,21 @@ You should have received a copy of the GNU General Public License
along with this program. If not, see <http://www.gnu.org/licenses/>.
"""
import numpy as np
import gc
import logging
#Simulator engine etc
# Simulator engine etc
from GPUSimulators import SHMEMSimulatorGroup
from GPUSimulators.common import common
from GPUSimulators.common import common, run_simulation
from GPUSimulators import EE2D_KP07_dimsplit
from GPUSimulators.helpers import InitialConditions as IC
####
#Initialize logging
# Initialize logging
####
log_level_console = 20
log_level_file = 10
log_level_file = 10
log_filename = 'shmem.log'
logger = logging.getLogger('GPUSimulators')
logger.setLevel(min(log_level_console, log_level_file))
@@ -42,16 +41,14 @@ logger.setLevel(min(log_level_console, log_level_file))
ch = logging.StreamHandler()
ch.setLevel(log_level_console)
logger.addHandler(ch)
logger.info("Console logger using level %s", logging.getLevelName(log_level_console))
logger.info(f"Console logger using level {logging.getLevelName(log_level_console)}")
fh = logging.FileHandler(log_filename)
formatter = logging.Formatter('%(asctime)s:%(name)s:%(levelname)s: %(message)s')
fh.setFormatter(formatter)
fh.setLevel(log_level_file)
logger.addHandler(fh)
logger.info("File logger using level %s to %s", logging.getLevelName(log_level_file), log_filename)
logger.info(f"File logger using level {logging.getLevelName(log_level_file)} to {log_filename}")
####
# Initialize SHMEM grid etc
@@ -59,10 +56,8 @@ logger.info("File logger using level %s to %s", logging.getLevelName(log_level_f
logger.info("Creating SHMEM grid")
# XXX: need to explicitly set ngpus when testing on single-GPU system
grid = SHMEMSimulatorGroup.SHMEMGrid(ngpus=4)
# XXX: need to explicitly set ngpus when testing on a single-GPU system
grid = SHMEMSimulatorGroup.SHMEMGrid(ngpus=4)
####
# Set initial conditions
@@ -71,7 +66,7 @@ logger.info("Generating initial conditions")
nx = 128
ny = 128
gamma = 1.4
#save_times = np.linspace(0, 0.01, 10)
# save_times = np.linspace(0, 0.01, 10)
save_times = np.linspace(0, 10, 10)
save_var_names = ['rho', 'rho_u', 'rho_v', 'E']
@@ -89,20 +84,20 @@ for i in range(grid.ngpus):
arguments['theta'] = 1.2
sims.append(EE2D_KP07_dimsplit.EE2D_KP07_dimsplit(**arguments))
#sims[i] = SHMEMSimulator(i, local_sim, grid) # 1st attempt: no wrapper (per sim)
# sims[i] = SHMEMSimulator(i, local_sim, grid) # 1st attempt: no wrapper (per sim)
arguments['sims'] = sims
arguments['grid'] = grid
#Helper function to create SHMEM simulator
def genSim(sims, grid, **kwargs):
# Helper function to create SHMEM simulator
def gen_sim(sims, grid, **kwargs):
# XXX: kwargs not used, since the simulators are already instantiated in the for-loop above
sim = SHMEMSimulatorGroup.SHMEMSimulatorGroup(sims, grid)
return sim
outfile = Common.run_simulation(genSim, arguments, outfile, save_times, save_var_names)
outfile = run_simulation(gen_sim, arguments, outfile, save_times, save_var_names)
####
# Clean shutdown
@@ -114,10 +109,8 @@ arguments = None
logging.shutdown()
gc.collect()
####
# Print completion and exit
####
print("Completed!")
exit(0)
exit(0)

28
singleGPUTesting.py
View File

@@ -19,7 +19,6 @@ You should have received a copy of the GNU General Public License
along with this program. If not, see <http://www.gnu.org/licenses/>.
"""
import numpy as np
import gc
import logging
@@ -28,17 +27,17 @@ import logging
import pycuda.driver as cuda
# Simulator engine etc
from GPUSimulators.common import common
from GPUSimulators.gpu import cuda_context
from GPUSimulators.common import run_simulation
from GPUSimulators.gpu import CudaContext
from GPUSimulators import EE2D_KP07_dimsplit
from GPUSimulators.helpers import InitialConditions as IC
import argparse
parser = argparse.ArgumentParser(description='Single GPU testing.')
parser.add_argument('-nx', type=int, default=128)
parser.add_argument('-ny', type=int, default=128)
args = parser.parse_args()
####
@@ -53,8 +52,7 @@ logger.setLevel(min(log_level_console, log_level_file))
ch = logging.StreamHandler()
ch.setLevel(log_level_console)
logger.addHandler(ch)
logger.info("Console logger using level %s",
logging.getLevelName(log_level_console))
logger.info(f"Console logger using level {logging.getLevelName(log_level_console)}")
fh = logging.FileHandler(log_filename)
formatter = logging.Formatter(
@@ -62,17 +60,14 @@ formatter = logging.Formatter(
fh.setFormatter(formatter)
fh.setLevel(log_level_file)
logger.addHandler(fh)
logger.info("File logger using level %s to %s",
logging.getLevelName(log_level_file), log_filename)
logger.info(f"File logger using level {logging.getLevelName(log_level_file)} to {log_filename}")
####
# Initialize CUDA
####
cuda.init(flags=0)
logger.info("Initializing CUDA")
cuda_context = CudaContext.CudaContext(autotuning=False)
cuda_context = CudaContext(autotuning=False)
####
# Set initial conditions
@@ -91,21 +86,22 @@ arguments = IC.genKelvinHelmholtz(nx, ny, gamma)
arguments['context'] = cuda_context
arguments['theta'] = 1.2
####
# Run simulation
####
logger.info("Running simulation")
# Helper function to create MPI simulator
def genSim(**kwargs):
def gen_sim(**kwargs):
local_sim = EE2D_KP07_dimsplit.EE2D_KP07_dimsplit(**kwargs)
return local_sim
outfile = Common.run_simulation(
genSim, arguments, outfile, save_times, save_var_names)
outfile = run_simulation(
gen_sim, arguments, outfile, save_times, save_var_names)
####
# Clean shutdown
@@ -116,8 +112,6 @@ arguments = None
logging.shutdown()
gc.collect()
####
# Print completion and exit
####