{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Ball and stick 4: Parallel vs serial mode"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This is the final part of the series where we build a multicompartment cell and evolve it into a network of cells running on a parallel machine (which is basically all computers made within the last decade). On this page, we translate the classes we have previously constructed so that they operate in either a parallel or serial mode."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"**Note:** If you followed the [installation instructions](../install/install_instructions.html) on the NEURON website, you should have everything you need to run parallel simulations.\n",
"If not, if you do not already have an MPI installation, go to that link and follow the instructions related to MPI.\n",
"If you compiled NEURON yourself instead of using an installer (this is rarely necessary), this part of the tutorial requires you to have used the `-DNRN_ENABLE_MPI=ON` flag at configure time."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Parallel communication in NEURON"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Parallel communication takes place via logical events in network connection objects known as [NetCon](../python/modelspec/programmatic/network/netcon.html).\n",
"NetCon sources are threshold detectors. They monitor some variable, say the membrane potential of a cell, and when the variable reaches some threshold, it triggers an event sent to the targets. Targets are frequently synapses on other cells. When they receive the event, they activate.\n",
"\n",
"In a parallel context across several machines, communication between hosts can be computationally inefficient when the frequency of events is high and when the message being sent is large. NEURON uses an efficient priority queueing mechanism to deliver events to targets after the delay specified by the NetCon. The message passed is succinct. It is an integer, the unique global identifier (gid) of the source.\n",
"The following two figures illustrate these ideas and come from [Hines M.L. and Carnevale N.T, Translating network models to parallel hardware in NEURON, Journal of Neuroscience Methods 169 (2008) 425–455](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2430920).\n",
"Users should also consult the [ParallelContext](../python/modelspec/programmatic/network/parcon.html) reference."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"
"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The main step involved in making a parallel implementation is to assign the global identifiers across the various hosts. Care should also be taken to assign cells to the various hosts such that the system is load balanced. For example, in a network with computationally complex and simple cells, several simple cells may be assigned to a host while few complex cells may be assigned to another host."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Test MPI and Parallel NEURON"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Once again, parallel NEURON requires MPI support. If this is your first time using it in a while, you should test your computer setup."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"To do this, first we will write out a small test script.\n",
"Executing the following cell will create a file called `testmpi.py`:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:47.920863Z",
"iopub.status.busy": "2025-08-18T03:35:47.920714Z",
"iopub.status.idle": "2025-08-18T03:35:47.929759Z",
"shell.execute_reply": "2025-08-18T03:35:47.929406Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Writing testmpi.py\n"
]
}
],
"source": [
"%%writefile testmpi.py\n",
"from neuron import n\n",
"n.nrnmpi_init() # initialize MPI\n",
"pc = n.ParallelContext()\n",
"print('I am {} of {}'.format(pc.id(), pc.nhost()))\n",
"n.quit() # necessary to avoid a warning message on parallel exit on some systems"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"And we will test running this script using `mpiexec`.\n",
"Normally we can just write `python` in place of `$python_exe`, but using `sys.executable` in this way can be necessary on systems with multiple Python versions."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:47.994627Z",
"iopub.status.busy": "2025-08-18T03:35:47.994416Z",
"iopub.status.idle": "2025-08-18T03:35:47.998256Z",
"shell.execute_reply": "2025-08-18T03:35:47.997857Z"
}
},
"outputs": [],
"source": [
"import os, sys\n",
"\n",
"os.environ[\"python_exe\"] = sys.executable"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now we can run the script in parallel with `mpiexec -n 4 python testmpi.py`.\n",
"\n",
"You can also run this directly from the command line in a terminal; most likely this is what you will want to do when running larger simulations or using a shared compute cluster."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:48.000754Z",
"iopub.status.busy": "2025-08-18T03:35:48.000602Z",
"iopub.status.idle": "2025-08-18T03:35:51.724517Z",
"shell.execute_reply": "2025-08-18T03:35:51.722083Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"numprocs=4\r\n",
"I am 3 of 4\r\n",
"[WARNING] yaksa: 1 leaked handle pool objects\r\n",
"I am 1 of 4\r\n",
"[WARNING] yaksa: 1 leaked handle pool objects\r\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"I am 0 of 4\r\n",
"[WARNING] yaksa: 1 leaked handle pool objects\r\n",
"I am 2 of 4\r\n",
"[WARNING] yaksa: 1 leaked handle pool objects\r\n"
]
}
],
"source": [
"!mpiexec -n 4 $python_exe testmpi.py"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"You should see something similar to:\n",
"```bash\n",
"numprocs=4\n",
"I am 1 of 4\n",
"I am 2 of 4\n",
"I am 3 of 4\n",
"I am 0 of 4\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"These could appear in any order since in theory they are running simultaneously but must print out in some order.\n",
"If instead you see four processes claiming to be 0 of 1, then your copy of NEURON was not compiled with support for parallel simulation.\n",
"Reconfigure with the [-DNRN_ENABLE_MPI=ON flag](../cmake_doc/options.html#mpi-options), recompile, and try again.\n",
"\n",
"If you get an error saying that `mpiexec` is an unknown command, then MPI is either not installed or not on your PATH; correct your MPI setup and try again."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Back to the example"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This notebook requires that `ballandstick.py` is in your working directory.\n",
"This is equivalent to the classes we created in the previous part of the tutorial."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We will construct a `ring.py` based on the previous `Ring` class. Changes are indicated with `###`:"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:51.729507Z",
"iopub.status.busy": "2025-08-18T03:35:51.726629Z",
"iopub.status.idle": "2025-08-18T03:35:51.741598Z",
"shell.execute_reply": "2025-08-18T03:35:51.739344Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Writing ring.py\n"
]
}
],
"source": [
"%%writefile ring.py\n",
"from neuron import n\n",
"from ballandstick import BallAndStick\n",
"\n",
"### MPI must be initialized before we create a ParallelContext object\n",
"n.nrnmpi_init()\n",
"pc = n.ParallelContext()\n",
"\n",
"class Ring:\n",
" \"\"\"A network of *N* ball-and-stick cells where cell n makes an\n",
" excitatory synapse onto cell n + 1 and the last, Nth cell in the\n",
" network projects to the first cell.\n",
" \"\"\"\n",
" def __init__(self, N=5, stim_w=0.04, stim_t=9, stim_delay=1, syn_w=0.01, syn_delay=5, r=50):\n",
" \"\"\"\n",
" :param N: Number of cells.\n",
" :param stim_w: Weight of the stimulus\n",
" :param stim_t: time of the stimulus (in ms)\n",
" :param stim_delay: delay of the stimulus (in ms)\n",
" :param syn_w: Synaptic weight\n",
" :param syn_delay: Delay of the synapse\n",
" :param r: radius of the network\n",
" \"\"\" \n",
" self._N = N\n",
" self.set_gids() ### assign gids to processors\n",
" self._syn_w = syn_w\n",
" self._syn_delay = syn_delay\n",
" self._create_cells(r) ### changed to use self._N instead of passing in N\n",
" self._connect_cells()\n",
" ### the 0th cell only exists on one process... that's the only one that gets a netstim\n",
" if pc.gid_exists(0):\n",
" self._netstim = n.NetStim()\n",
" self._netstim.number = 1\n",
" self._netstim.start = stim_t\n",
" self._nc = n.NetCon(self._netstim, pc.gid2cell(0).syn) ### grab cell with gid==0 wherever it exists\n",
" self._nc.delay = stim_delay\n",
" self._nc.weight[0] = stim_w\n",
" \n",
" def set_gids(self):\n",
" \"\"\"Set the gidlist on this host.\"\"\"\n",
" #### Round-robin counting.\n",
" #### Each host has an id from 0 to pc.nhost() - 1.\n",
" self.gidlist = list(range(pc.id(), self._N, pc.nhost()))\n",
" for gid in self.gidlist:\n",
" pc.set_gid2node(gid, pc.id())\n",
" \n",
" def _create_cells(self, r):\n",
" self.cells = []\n",
" for i in self.gidlist: ### only create the cells that exist on this host\n",
" theta = i * 2 * n.PI / self._N\n",
" self.cells.append(BallAndStick(i, n.cos(theta) * r, n.sin(theta) * r, 0, theta))\n",
" ### associate the cell with this host and gid\n",
" for cell in self.cells:\n",
" pc.cell(cell._gid, cell._spike_detector)\n",
"\n",
" def _connect_cells(self):\n",
" ### this method is different because we now must use ids instead of objects\n",
" for target in self.cells:\n",
" source_gid = (target._gid - 1 + self._N) % self._N\n",
" nc = pc.gid_connect(source_gid, target.syn)\n",
" nc.weight[0] = self._syn_w\n",
" nc.delay = self._syn_delay\n",
" target._ncs.append(nc)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The call to `n.nrnmpi_init()` must happen before any use of the [ParallelContext](../python/modelspec/programmatic/network/parcon.html) class -- which forms a key part of any NEURON parallel simulation.\n",
"\n",
"The only conceptually new method here is the `set_gids` method where each process specifies which cells it will simulate. Here we use what is known as a round-robin approach, where the `pc.id()`th process starts at `pc.id()` and skips by however many processes are running (`pc.nhost`)."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In `_create_cells`, we now associate the cell `gid` with the NetCon `_spike_detector`. This allows the `_connect_cells` to make connections based on gids instead of objects, using `pc.gid_connect`."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Test the model"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Here's a basic test, `test_ring1.py` that loads the `Ring` class and plots cell 0's membrane potential timeseries:"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:51.744567Z",
"iopub.status.busy": "2025-08-18T03:35:51.743076Z",
"iopub.status.idle": "2025-08-18T03:35:51.752581Z",
"shell.execute_reply": "2025-08-18T03:35:51.750547Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Writing test_ring1.py\n"
]
}
],
"source": [
"%%writefile test_ring1.py\n",
"from neuron import n\n",
"from neuron.units import ms, mV\n",
"import matplotlib.pyplot as plt\n",
"from ring import Ring\n",
"\n",
"cell_to_plot = 0\n",
"\n",
"ring = Ring()\n",
"\n",
"pc = n.ParallelContext()\n",
"pc.set_maxstep(10 * ms)\n",
"\n",
"t = n.Vector().record(n._ref_t)\n",
"n.finitialize(-65 * mV)\n",
"pc.psolve(100 * ms)\n",
"\n",
"if pc.gid_exists(cell_to_plot):\n",
" plt.figure()\n",
" plt.title(\"Cell {}\".format(cell_to_plot))\n",
" plt.plot(t, pc.gid2cell(cell_to_plot).soma_v)\n",
" plt.xlabel(\"Simulation time [ms]\")\n",
" plt.ylabel(\"Soma voltage [mV]\")\n",
" plt.savefig(\"test_ring1_{}ranks.svgz\".format(pc.nhost()))\n",
"\n",
"pc.barrier()\n",
"pc.done()\n",
"n.quit()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The code above should look very familiar."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The conceptually new pieces are:\n",
"\n",
"* [pc.set_maxstep(10 * ms)](../python/modelspec/programmatic/network/parcon.html#ParallelContext.set_maxstep) -- sets an upper bound on how far MPI can simulate without communicating, here a simulated 10 ms. This *must* be called before attempting a parallel simulation.\n",
"* [pc.psolve(100 * ms)](../python/modelspec/programmatic/network/parcon.html#ParallelContext.psolve) -- a parallel version of [n.continuerun](../python/simctrl/stdrun.html), but does not support updating NEURON graphics during the simulation.\n",
"* [pc.gid_exists](../python/modelspec/programmatic/network/parcon.html#ParallelContext.gid_exists) -- only the process that owns the specified cell should make the plot.\n",
"* [pc.gid2cell](../python/modelspec/programmatic/network/parcon.html#ParallelContext.gid2cell) -- lookup a cell by gid.\n",
"* [pc.barrier()](../python/modelspec/programmatic/network/parcon.html#ParallelContext.barrier) -- wait until all processes reach this point; used to make sure processes don't shut down before the graph is closed."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Start by testing this without using MPI:"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:51.755938Z",
"iopub.status.busy": "2025-08-18T03:35:51.755640Z",
"iopub.status.idle": "2025-08-18T03:35:53.723317Z",
"shell.execute_reply": "2025-08-18T03:35:53.722875Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"numprocs=1\r\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"[WARNING] yaksa: 1 leaked handle pool objects\r\n"
]
},
{
"data": {
"image/svg+xml": [
""
],
"text/plain": [
""
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"!$python_exe test_ring1.py\n",
"\n",
"\n",
"def SVGZ(name):\n",
" from gzip import GzipFile\n",
" from IPython.display import SVG\n",
"\n",
" return SVG(data=GzipFile(name + \".svgz\").read())\n",
"\n",
"\n",
"SVGZ(\"test_ring1_1ranks\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This should look similar to the following reference image:"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:53.728110Z",
"iopub.status.busy": "2025-08-18T03:35:53.727940Z",
"iopub.status.idle": "2025-08-18T03:35:53.743103Z",
"shell.execute_reply": "2025-08-18T03:35:53.742748Z"
}
},
"outputs": [
{
"data": {
"image/svg+xml": [
""
],
"text/plain": [
""
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"SVGZ(\"test_ring1_ref\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"And once you are satisfied that works, you can try MPI, e.g."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:53.744947Z",
"iopub.status.busy": "2025-08-18T03:35:53.744804Z",
"iopub.status.idle": "2025-08-18T03:35:57.206644Z",
"shell.execute_reply": "2025-08-18T03:35:57.206138Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"numprocs=2\r\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"[WARNING] yaksa: 1 leaked handle pool objects\r\n",
"[WARNING] yaksa: 1 leaked handle pool objects\r\n"
]
},
{
"data": {
"image/svg+xml": [
""
],
"text/plain": [
""
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"!mpiexec -n 2 $python_exe test_ring1.py\n",
"SVGZ(\"test_ring1_2ranks\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Gathering spikes"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Our above test runs the simulation successfully, but in the end, no single process knows when all the spikes occurred.\n",
"There are a number of ways to deal with this: one solution is to have each process write its data to a file.\n",
"Instead, we will use [pc.py_alltoall](../python/modelspec/programmatic/network/parcon.html#ParallelContext.py_alltoall) to send all the data to MPI rank 0 (`pc.id() == 0`), at which point it can plot the raster, save data, or whatever."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The following code does this and saves itself to a file called `test_ring2.py`:"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:57.211108Z",
"iopub.status.busy": "2025-08-18T03:35:57.210815Z",
"iopub.status.idle": "2025-08-18T03:35:57.219556Z",
"shell.execute_reply": "2025-08-18T03:35:57.218287Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Writing test_ring2.py\n"
]
}
],
"source": [
"%%writefile test_ring2.py\n",
"from neuron import n\n",
"from neuron.units import ms, mV\n",
"import matplotlib.pyplot as plt\n",
"from ring import Ring\n",
"\n",
"ring = Ring()\n",
"\n",
"pc = n.ParallelContext()\n",
"pc.set_maxstep(10 * ms)\n",
"\n",
"t = n.Vector().record(n._ref_t)\n",
"n.finitialize(-65 * mV)\n",
"pc.psolve(100 * ms)\n",
"\n",
"# send all spike time data to rank 0\n",
"local_data = {cell._gid: list(cell.spike_times) for cell in ring.cells}\n",
"all_data = pc.py_alltoall([local_data] + [None] * (pc.nhost() - 1))\n",
"\n",
"if pc.id() == 0:\n",
" # combine the data from the various processes\n",
" data = {}\n",
" for process_data in all_data:\n",
" data.update(process_data)\n",
" # plot it\n",
" plt.figure()\n",
" plt.title(\"Spike raster\")\n",
" plt.xlabel(\"Simulation time [ms]\")\n",
" plt.ylabel(\"Cell\")\n",
" for i, spike_times in data.items():\n",
" plt.vlines(spike_times, i + 0.5, i + 1.5)\n",
" plt.savefig(\"test_ring2.svgz\")\n",
"\n",
"pc.barrier()\n",
"pc.done()\n",
"n.quit()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We can run this using `mpiexec` as before, here with two processes:"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:35:57.223502Z",
"iopub.status.busy": "2025-08-18T03:35:57.222466Z",
"iopub.status.idle": "2025-08-18T03:36:00.654668Z",
"shell.execute_reply": "2025-08-18T03:36:00.652342Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"numprocs=2\r\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"[WARNING] yaksa: 1 leaked handle pool objects\r\n",
"[WARNING] yaksa: 1 leaked handle pool objects\r\n"
]
},
{
"data": {
"image/svg+xml": [
""
],
"text/plain": [
""
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"!mpiexec -n 2 $python_exe test_ring2.py\n",
"SVGZ(\"test_ring2\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This should display the familiar raster plot.\n",
"If you are wondering why node 0 was the one chosen to make the plot, it is because that is the only node that is guaranteed to exist (nothing else exists if there is only one process being used for the simulation).\n",
"A reference version of the raster plot is shown below; the two plots should look the same:"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"execution": {
"iopub.execute_input": "2025-08-18T03:36:00.659205Z",
"iopub.status.busy": "2025-08-18T03:36:00.658863Z",
"iopub.status.idle": "2025-08-18T03:36:00.673967Z",
"shell.execute_reply": "2025-08-18T03:36:00.673617Z"
}
},
"outputs": [
{
"data": {
"image/svg+xml": [
""
],
"text/plain": [
""
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"SVGZ(\"test_ring2_ref\")"
]
}
],
"metadata": {
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.12.10"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
![](\"https://neuron.yale.edu/neuron/static/docs/neuronpython/_images/parfig2.1.png\"/)
![](\"https://neuron.yale.edu/neuron/static/docs/neuronpython/_images/parfig2.2.png\"/)