Skip to content
Snippets Groups Projects
Commit ca6e659e authored by Paul McCarthy's avatar Paul McCarthy :mountain_bicyclist:
Browse files

added section on shared memory in threading/multiprocessing practical

parent ffaad040
No related branches found
No related tags found
No related merge requests found
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
# Threading and parallel processing # Threading and parallel processing
The Python language has built-in support for multi-threading in the The Python language has built-in support for multi-threading in the
[`threading`](https://docs.python.org/3.5/library/threading.html) module, and [`threading`](https://docs.python.org/3.5/library/threading.html) module, and
true parallelism in the true parallelism in the
[`multiprocessing`](https://docs.python.org/3.5/library/multiprocessing.html) [`multiprocessing`](https://docs.python.org/3.5/library/multiprocessing.html)
module. If you want to be impressed, skip straight to the section on module. If you want to be impressed, skip straight to the section on
[`multiprocessing`](todo). [`multiprocessing`](todo).
## Threading ## Threading
The [`threading`](https://docs.python.org/3.5/library/threading.html) module The [`threading`](https://docs.python.org/3.5/library/threading.html) module
provides a traditional multi-threading API that should be familiar to you if provides a traditional multi-threading API that should be familiar to you if
you have worked with threads in other languages. you have worked with threads in other languages.
Running a task in a separate thread in Python is easy - simply create a Running a task in a separate thread in Python is easy - simply create a
`Thread` object, and pass it the function or method that you want it to `Thread` object, and pass it the function or method that you want it to
run. Then call its `start` method: run. Then call its `start` method:
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
import time import time
import threading import threading
def longRunningTask(niters): def longRunningTask(niters):
for i in range(niters): for i in range(niters):
if i % 2 == 0: print('Tick') if i % 2 == 0: print('Tick')
else: print('Tock') else: print('Tock')
time.sleep(0.5) time.sleep(0.5)
t = threading.Thread(target=longRunningTask, args=(8,)) t = threading.Thread(target=longRunningTask, args=(8,))
t.start() t.start()
while t.is_alive(): while t.is_alive():
time.sleep(0.4) time.sleep(0.4)
print('Waiting for thread to finish...') print('Waiting for thread to finish...')
print('Finished!') print('Finished!')
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
You can also `join` a thread, which will block execution in the current thread You can also `join` a thread, which will block execution in the current thread
until the thread that has been `join`ed has finished: until the thread that has been `join`ed has finished:
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
t = threading.Thread(target=longRunningTask, args=(6, )) t = threading.Thread(target=longRunningTask, args=(6, ))
t.start() t.start()
print('Joining thread ...') print('Joining thread ...')
t.join() t.join()
print('Finished!') print('Finished!')
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### Subclassing `Thread` ### Subclassing `Thread`
It is also possible to sub-class the `Thread` class, and override its `run` It is also possible to sub-class the `Thread` class, and override its `run`
method: method:
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
class LongRunningThread(threading.Thread): class LongRunningThread(threading.Thread):
def __init__(self, niters, *args, **kwargs): def __init__(self, niters, *args, **kwargs):
super().__init__(*args, **kwargs) super().__init__(*args, **kwargs)
self.niters = niters self.niters = niters
def run(self): def run(self):
for i in range(self.niters): for i in range(self.niters):
if i % 2 == 0: print('Tick') if i % 2 == 0: print('Tick')
else: print('Tock') else: print('Tock')
time.sleep(0.5) time.sleep(0.5)
t = LongRunningThread(6) t = LongRunningThread(6)
t.start() t.start()
t.join() t.join()
print('Done') print('Done')
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### Daemon threads ### Daemon threads
By default, a Python application will not exit until _all_ active threads have By default, a Python application will not exit until _all_ active threads have
finished execution. If you want to run a task in the background for the finished execution. If you want to run a task in the background for the
duration of your application, you can mark it as a `daemon` thread - when all duration of your application, you can mark it as a `daemon` thread - when all
non-daemon threads in a Python application have finished, all daemon threads non-daemon threads in a Python application have finished, all daemon threads
will be halted, and the application will exit. will be halted, and the application will exit.
You can mark a thread as being a daemon by setting an attribute on it after You can mark a thread as being a daemon by setting an attribute on it after
creation: creation:
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
t = threading.Thread(target=longRunningTask) t = threading.Thread(target=longRunningTask)
t.daemon = True t.daemon = True
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
See the [`Thread` See the [`Thread`
documentation](https://docs.python.org/3.5/library/threading.html#thread-objects) documentation](https://docs.python.org/3.5/library/threading.html#thread-objects)
for more details. for more details.
### Thread synchronisation ### Thread synchronisation
The `threading` module provides some useful thread-synchronisation primitives The `threading` module provides some useful thread-synchronisation primitives
- the `Lock`, `RLock` (re-entrant `Lock`), and `Event` classes. The - the `Lock`, `RLock` (re-entrant `Lock`), and `Event` classes. The
`threading` module also provides `Condition` and `Semaphore` classes - refer `threading` module also provides `Condition` and `Semaphore` classes - refer
to the [documentation](https://docs.python.org/3.5/library/threading.html) for to the [documentation](https://docs.python.org/3.5/library/threading.html) for
more details. more details.
#### `Lock` #### `Lock`
The [`Lock`](https://docs.python.org/3.5/library/threading.html#lock-objects) The [`Lock`](https://docs.python.org/3.5/library/threading.html#lock-objects)
class (and its re-entrant version, the class (and its re-entrant version, the
[`RLock`](https://docs.python.org/3.5/library/threading.html#rlock-objects)) [`RLock`](https://docs.python.org/3.5/library/threading.html#rlock-objects))
prevents a block of code from being accessed by more than one thread at a prevents a block of code from being accessed by more than one thread at a
time. For example, if we have multiple threads running this `task` function, time. For example, if we have multiple threads running this `task` function,
their [outputs](https://www.youtube.com/watch?v=F5fUFnfPpYU) will inevitably their [outputs](https://www.youtube.com/watch?v=F5fUFnfPpYU) will inevitably
become intertwined: become intertwined:
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
def task(): def task():
for i in range(5): for i in range(5):
print('{} Woozle '.format(i), end='') print('{} Woozle '.format(i), end='')
time.sleep(0.1) time.sleep(0.1)
print('Wuzzle') print('Wuzzle')
threads = [threading.Thread(target=task) for i in range(5)] threads = [threading.Thread(target=task) for i in range(5)]
for t in threads: for t in threads:
t.start() t.start()
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
But if we protect the critical section with a `Lock` object, the output will But if we protect the critical section with a `Lock` object, the output will
look more sensible: look more sensible:
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
lock = threading.Lock() lock = threading.Lock()
def task(): def task():
for i in range(5): for i in range(5):
with lock: with lock:
print('{} Woozle '.format(i), end='') print('{} Woozle '.format(i), end='')
time.sleep(0.1) time.sleep(0.1)
print('Wuzzle') print('Wuzzle')
threads = [threading.Thread(target=task) for i in range(5)] threads = [threading.Thread(target=task) for i in range(5)]
for t in threads: for t in threads:
t.start() t.start()
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
> Instead of using a `Lock` object in a `with` statement, it is also possible > Instead of using a `Lock` object in a `with` statement, it is also possible
> to manually call its `acquire` and `release` methods: > to manually call its `acquire` and `release` methods:
> >
> def task(): > def task():
> for i in range(5): > for i in range(5):
> lock.acquire() > lock.acquire()
> print('{} Woozle '.format(i), end='') > print('{} Woozle '.format(i), end='')
> time.sleep(0.1) > time.sleep(0.1)
> print('Wuzzle') > print('Wuzzle')
> lock.release() > lock.release()
Python does not have any built-in constructs to implement `Lock`-based mutual Python does not have any built-in constructs to implement `Lock`-based mutual
exclusion across several functions or methods - each function/method must exclusion across several functions or methods - each function/method must
explicitly acquire/release a shared `Lock` instance. However, it is relatively explicitly acquire/release a shared `Lock` instance. However, it is relatively
straightforward to implement a decorator which does this for you: straightforward to implement a decorator which does this for you:
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
def mutex(func, lock): def mutex(func, lock):
def wrapper(*args): def wrapper(*args):
with lock: with lock:
func(*args) func(*args)
return wrapper return wrapper
class MyClass(object): class MyClass(object):
def __init__(self): def __init__(self):
lock = threading.Lock() lock = threading.Lock()
self.safeFunc1 = mutex(self.safeFunc1, lock) self.safeFunc1 = mutex(self.safeFunc1, lock)
self.safeFunc2 = mutex(self.safeFunc2, lock) self.safeFunc2 = mutex(self.safeFunc2, lock)
def safeFunc1(self): def safeFunc1(self):
time.sleep(0.1) time.sleep(0.1)
print('safeFunc1 start') print('safeFunc1 start')
time.sleep(0.2) time.sleep(0.2)
print('safeFunc1 end') print('safeFunc1 end')
def safeFunc2(self): def safeFunc2(self):
time.sleep(0.1) time.sleep(0.1)
print('safeFunc2 start') print('safeFunc2 start')
time.sleep(0.2) time.sleep(0.2)
print('safeFunc2 end') print('safeFunc2 end')
mc = MyClass() mc = MyClass()
f1threads = [threading.Thread(target=mc.safeFunc1) for i in range(4)] f1threads = [threading.Thread(target=mc.safeFunc1) for i in range(4)]
f2threads = [threading.Thread(target=mc.safeFunc2) for i in range(4)] f2threads = [threading.Thread(target=mc.safeFunc2) for i in range(4)]
for t in f1threads + f2threads: for t in f1threads + f2threads:
t.start() t.start()
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
Try removing the `mutex` lock from the two methods in the above code, and see Try removing the `mutex` lock from the two methods in the above code, and see
what it does to the output. what it does to the output.
#### `Event` #### `Event`
The The
[`Event`](https://docs.python.org/3.5/library/threading.html#event-objects) [`Event`](https://docs.python.org/3.5/library/threading.html#event-objects)
class is essentially a boolean [semaphore][semaphore-wiki]. It can be used to class is essentially a boolean [semaphore][semaphore-wiki]. It can be used to
signal events between threads. Threads can `wait` on the event, and be awoken signal events between threads. Threads can `wait` on the event, and be awoken
when the event is `set` by another thread: when the event is `set` by another thread:
[semaphore-wiki]: https://en.wikipedia.org/wiki/Semaphore_(programming) [semaphore-wiki]: https://en.wikipedia.org/wiki/Semaphore_(programming)
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
import numpy as np import numpy as np
processingFinished = threading.Event() processingFinished = threading.Event()
def processData(data): def processData(data):
print('Processing data ...') print('Processing data ...')
time.sleep(2) time.sleep(2)
print('Result: {}'.format(data.mean())) print('Result: {}'.format(data.mean()))
processingFinished.set() processingFinished.set()
data = np.random.randint(1, 100, 100) data = np.random.randint(1, 100, 100)
t = threading.Thread(target=processData, args=(data,)) t = threading.Thread(target=processData, args=(data,))
t.start() t.start()
processingFinished.wait() processingFinished.wait()
print('Processing finished!') print('Processing finished!')
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
### The Global Interpreter Lock (GIL) ### The Global Interpreter Lock (GIL)
The [_Global Interpreter The [_Global Interpreter
Lock_](https://docs.python.org/3/c-api/init.html#thread-state-and-the-global-interpreter-lock) Lock_](https://docs.python.org/3/c-api/init.html#thread-state-and-the-global-interpreter-lock)
is an implementation detail of [CPython](https://github.com/python/cpython) is an implementation detail of [CPython](https://github.com/python/cpython)
(the official Python interpreter). The GIL means that a multi-threaded (the official Python interpreter). The GIL means that a multi-threaded
program written in pure Python is not able to take advantage of multiple program written in pure Python is not able to take advantage of multiple
cores - this essentially means that only one thread may be executing at any cores - this essentially means that only one thread may be executing at any
point in time. point in time.
The `threading` module does still have its uses though, as this GIL problem The `threading` module does still have its uses though, as this GIL problem
does not affect tasks which involve calls to system or natively compiled does not affect tasks which involve calls to system or natively compiled
libraries (e.g. file and network I/O, Numpy operations, etc.). So you can, libraries (e.g. file and network I/O, Numpy operations, etc.). So you can,
for example, perform some expensive processing on a Numpy array in a thread for example, perform some expensive processing on a Numpy array in a thread
running on one core, whilst having another thread (e.g. user interaction) running on one core, whilst having another thread (e.g. user interaction)
running on another core. running on another core.
## Multiprocessing ## Multiprocessing
For true parallelism, you should check out the For true parallelism, you should check out the
[`multiprocessing`](https://docs.python.org/3.5/library/multiprocessing.html) [`multiprocessing`](https://docs.python.org/3.5/library/multiprocessing.html)
module. module.
The `multiprocessing` module spawns sub-processes, rather than threads, and so The `multiprocessing` module spawns sub-processes, rather than threads, and so
is not subject to the GIL constraints that the `threading` module suffers is not subject to the GIL constraints that the `threading` module suffers
from. It provides two APIs - a "traditional" equivalent to that provided by from. It provides two APIs - a "traditional" equivalent to that provided by
the `threading` module, and a powerful higher-level API. the `threading` module, and a powerful higher-level API.
### `threading`-equivalent API ### `threading`-equivalent API
The The
[`Process`](https://docs.python.org/3.5/library/multiprocessing.html#the-process-class) [`Process`](https://docs.python.org/3.5/library/multiprocessing.html#the-process-class)
class is the `multiprocessing` equivalent of the class is the `multiprocessing` equivalent of the
[`threading.Thread`](https://docs.python.org/3.5/library/threading.html#thread-objects) [`threading.Thread`](https://docs.python.org/3.5/library/threading.html#thread-objects)
class. `multprocessing` also has equivalents of the [`Lock` and `Event` class. `multprocessing` also has equivalents of the [`Lock` and `Event`
classes](https://docs.python.org/3.5/library/multiprocessing.html#synchronization-between-processes), classes](https://docs.python.org/3.5/library/multiprocessing.html#synchronization-between-processes),
and the other synchronisation primitives provided by `threading`. and the other synchronisation primitives provided by `threading`.
So you can simply replace `threading.Thread` with `multiprocessing.Process`, So you can simply replace `threading.Thread` with `multiprocessing.Process`,
and you will have true parallelism. and you will have true parallelism.
Because your "threads" are now independent processes, you need to be a little Because your "threads" are now independent processes, you need to be a little
careful about how to share information across them. Fortunately, the careful about how to share information across them. Fortunately, the
`multiprocessing` module provides [`Queue` and `Pipe` `multiprocessing` module provides [`Queue` and `Pipe`
classes](https://docs.python.org/3.5/library/multiprocessing.html#exchanging-objects-between-processes) classes](https://docs.python.org/3.5/library/multiprocessing.html#exchanging-objects-between-processes)
which make it easy to share data across processes. which make it easy to share data across processes.
### Higher-level API - the `multiprocessing.Pool` ### Higher-level API - the `multiprocessing.Pool`
The real advantages of `multiprocessing` lie in its higher level API, centered The real advantages of `multiprocessing` lie in its higher level API, centered
around the [`Pool` around the [`Pool`
class](https://docs.python.org/3.5/library/multiprocessing.html#using-a-pool-of-workers). class](https://docs.python.org/3.5/library/multiprocessing.html#using-a-pool-of-workers).
Essentially, you create a `Pool` of worker processes - you specify the number Essentially, you create a `Pool` of worker processes - you specify the number
of processes when you create the pool. of processes when you create the pool.
> The best number of processes to use for a `Pool` will depend on the system > The best number of processes to use for a `Pool` will depend on the system
> you are running on (number of cores), and the tasks you are running (e.g. > you are running on (number of cores), and the tasks you are running (e.g.
> I/O bound or CPU bound). > I/O bound or CPU bound).
Once you have created a `Pool`, you can use its methods to automatically Once you have created a `Pool`, you can use its methods to automatically
parallelise tasks. The most useful are the `map`, `starmap` and parallelise tasks. The most useful are the `map`, `starmap` and
`apply_async` methods. `apply_async` methods.
#### `Pool.map` #### `Pool.map`
The The
[`Pool.map`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.Pool.map) [`Pool.map`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.Pool.map)
method is the multiprocessing equivalent of the built-in method is the multiprocessing equivalent of the built-in
[`map`](https://docs.python.org/3.5/library/functions.html#map) function - it [`map`](https://docs.python.org/3.5/library/functions.html#map) function - it
is given a function, and a sequence, and it applies the function to each is given a function, and a sequence, and it applies the function to each
element in the sequence. element in the sequence.
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
import time import time
import multiprocessing as mp import multiprocessing as mp
import numpy as np import numpy as np
def crunchImage(imgfile): def crunchImage(imgfile):
# Load a nifti image, do stuff # Load a nifti image, do stuff
# to it. Use your imagination # to it. Use your imagination
# to fill in this function. # to fill in this function.
time.sleep(2) time.sleep(2)
# numpy's random number generator # numpy's random number generator
# will be initialised in the same # will be initialised in the same
# way in each process, so let's # way in each process, so let's
# re-seed it. # re-seed it.
np.random.seed() np.random.seed()
result = np.random.randint(1, 100, 1) result = np.random.randint(1, 100, 1)
print(imgfile, ':', result) print(imgfile, ':', result)
return result return result
imgfiles = ['{:02d}.nii.gz'.format(i) for i in range(20)] imgfiles = ['{:02d}.nii.gz'.format(i) for i in range(20)]
p = mp.Pool(processes=16) p = mp.Pool(processes=16)
print('Crunching images...') print('Crunching images...')
start = time.time() start = time.time()
results = p.map(crunchImage, imgfiles) results = p.map(crunchImage, imgfiles)
end = time.time() end = time.time()
print('Total execution time: {:0.2f} seconds'.format(end - start)) print('Total execution time: {:0.2f} seconds'.format(end - start))
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
The `Pool.map` method only works with functions that accept one argument, such The `Pool.map` method only works with functions that accept one argument, such
as our `crunchImage` function above. If you have a function which accepts as our `crunchImage` function above. If you have a function which accepts
multiple arguments, use the multiple arguments, use the
[`Pool.starmap`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.Pool.starmap) [`Pool.starmap`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.Pool.starmap)
method instead: method instead:
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
def crunchImage(imgfile, modality): def crunchImage(imgfile, modality):
time.sleep(2) time.sleep(2)
np.random.seed() np.random.seed()
if modality == 't1': if modality == 't1':
result = np.random.randint(1, 100, 1) result = np.random.randint(1, 100, 1)
elif modality == 't2': elif modality == 't2':
result = np.random.randint(100, 200, 1) result = np.random.randint(100, 200, 1)
print(imgfile, ': ', result) print(imgfile, ': ', result)
return result return result
imgfiles = ['t1_{:02d}.nii.gz'.format(i) for i in range(10)] + \ imgfiles = ['t1_{:02d}.nii.gz'.format(i) for i in range(10)] + \
['t2_{:02d}.nii.gz'.format(i) for i in range(10)] ['t2_{:02d}.nii.gz'.format(i) for i in range(10)]
modalities = ['t1'] * 10 + ['t2'] * 10 modalities = ['t1'] * 10 + ['t2'] * 10
pool = mp.Pool(processes=16) pool = mp.Pool(processes=16)
args = [(f, m) for f, m in zip(imgfiles, modalities)] args = [(f, m) for f, m in zip(imgfiles, modalities)]
print('Crunching images...') print('Crunching images...')
start = time.time() start = time.time()
results = pool.starmap(crunchImage, args) results = pool.starmap(crunchImage, args)
end = time.time() end = time.time()
print('Total execution time: {:0.2f} seconds'.format(end - start)) print('Total execution time: {:0.2f} seconds'.format(end - start))
``` ```
%% Cell type:markdown id: tags: %% Cell type:markdown id: tags:
The `map` and `starmap` methods also have asynchronous equivalents `map_async` The `map` and `starmap` methods also have asynchronous equivalents `map_async`
and `starmap_async`, which return immediately. Refer to the and `starmap_async`, which return immediately. Refer to the
[`Pool`](https://docs.python.org/3.5/library/multiprocessing.html#module-multiprocessing.pool) [`Pool`](https://docs.python.org/3.5/library/multiprocessing.html#module-multiprocessing.pool)
documentation for more details. documentation for more details.
#### `Pool.apply_async` #### `Pool.apply_async`
The The
[`Pool.apply`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.Pool.apply) [`Pool.apply`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.Pool.apply)
method will execute a function on one of the processes, and block until it has method will execute a function on one of the processes, and block until it has
finished. The finished. The
[`Pool.apply_async`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.Pool.apply_async) [`Pool.apply_async`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.Pool.apply_async)
method returns immediately, and is thus more suited to asynchronously method returns immediately, and is thus more suited to asynchronously
scheduling multiple jobs to run in parallel. scheduling multiple jobs to run in parallel.
`apply_async` returns an object of type `apply_async` returns an object of type
[`AsyncResult`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.AsyncResult). [`AsyncResult`](https://docs.python.org/3.5/library/multiprocessing.html#multiprocessing.pool.AsyncResult).
An `AsyncResult` object has `wait` and `get` methods which will block until An `AsyncResult` object has `wait` and `get` methods which will block until
the job has completed. the job has completed.
%% Cell type:code id: tags: %% Cell type:code id: tags:
``` ```
import time import time
import multiprocessing as mp import multiprocessing as mp
import numpy as np import numpy as np
def linear_registration(src, ref): def linear_registration(src, ref):
time.sleep(1) time.sleep(1)
return np.eye(4) return np.eye(4)
def nonlinear_registration(src, ref, affine): def nonlinear_registration(src, ref, affine):
time.sleep(3) time.sleep(3)
# this number represents a non-linear warp # this number represents a non-linear warp
# field - use your imagination people! # field - use your imagination people!
np.random.seed() np.random.seed()
return np.random.randint(1, 100, 1) return np.random.randint(1, 100, 1)
t1s = ['{:02d}_t1.nii.gz'.format(i) for i in range(20)] t1s = ['{:02d}_t1.nii.gz'.format(i) for i in range(20)]
std = 'MNI152_T1_2mm.nii.gz' std = 'MNI152_T1_2mm.nii.gz'
pool = mp.Pool(processes=16) pool = mp.Pool(processes=16)
print('Running structural-to-standard registration ' print('Running structural-to-standard registration '
'on {} subjects...'.format(len(t1s))) 'on {} subjects...'.format(len(t1s)))
# Run linear registration on all the T1s. # Run linear registration on all the T1s.
# #
# We build a list of AsyncResult objects # We build a list of AsyncResult objects
linresults = [pool.apply_async(linear_registration, (t1, std)) linresults = [pool.apply_async(linear_registration, (t1, std))
for t1 in t1s] for t1 in t1s]
# Then we wait for each job to finish, # Then we wait for each job to finish,
# and replace its AsyncResult object # and replace its AsyncResult object
# with the actual result - an affine # with the actual result - an affine
# transformation matrix. # transformation matrix.
start = time.time() start = time.time()
for i, r in enumerate(linresults): for i, r in enumerate(linresults):
linresults[i] = r.get() linresults[i] = r.get()
end = time.time() end = time.time()
print('Linear registrations completed in ' print('Linear registrations completed in '
'{:0.2f} seconds'.format(end - start)) '{:0.2f} seconds'.format(end - start))
# Run non-linear registration on all the T1s, # Run non-linear registration on all the T1s,
# using the linear registrations to initialise. # using the linear registrations to initialise.
nlinresults = [pool.apply_async(nonlinear_registration, (t1, std, aff)) nlinresults = [pool.apply_async(nonlinear_registration, (t1, std, aff))
for (t1, aff) in zip(t1s, linresults)] for (t1, aff) in zip(t1s, linresults)]
# Wait for each non-linear reg to finish, # Wait for each non-linear reg to finish,
# and store the resulting warp field. # and store the resulting warp field.
start = time.time() start = time.time()
for i, r in enumerate(nlinresults): for i, r in enumerate(nlinresults):
nlinresults[i] = r.get() nlinresults[i] = r.get()
end = time.time() end = time.time()
print('Non-linear registrations completed in ' print('Non-linear registrations completed in '
'{:0.2f} seconds'.format(end - start)) '{:0.2f} seconds'.format(end - start))
print('Non linear registrations:') print('Non linear registrations:')
for t1, result in zip(t1s, nlinresults): for t1, result in zip(t1s, nlinresults):
print(t1, ':', result) print(t1, ':', result)
``` ```
%% Cell type:markdown id: tags:
### Sharing data between processes
When you use the `Pool.map` method (or any of the other methods we have shown)
to run a function on a sequence of items, those items must be copied into the
memory of each of the child processes. When the child processes are finished,
the data that they return then has to be copied back to the parent process.
Any items which you wish to pass to a function that is executed by a `Pool`
must be - the built-in
[`pickle`](https://docs.python.org/3.5/library/pickle.html) module is used by
`multiprocessing` to serialise and de-serialise the data passed into and
returned from a child process. The majority of standard Python types (`list`,
`dict`, `str` etc), and Numpy arrays can be pickled and unpickled, so you only
need to worry about this detail if you are passing objects of a custom type
(e.g. instances of classes that you have written, or that are defined in some
third-party library).
There is obviously some overhead in copying data back and forth between the
main process and the worker processes. For most computationally intensive
tasks, this communication overhead is not important - the performance
bottleneck is typically going to be the computation time, rather than I/O
between the parent and child processes. You may need to spend some time
adjusting the way in which you split up your data, and the number of
processes, in order to get the best performance.
However, if you have determined that copying data between processes is having
a substantial impact on your performance, the `multiprocessing` module
provides the [`Value`, `Array`, and `RawArray`
classes](https://docs.python.org/3.5/library/multiprocessing.html#shared-ctypes-objects),
which allow you to share individual values, or arrays of values, respectively.
The `Array` and `RawArray` classes essentially wrap a typed pointer (from the
built-in [`ctypes`](https://docs.python.org/3.5/library/ctypes.html) module)
to a block of memory. We can use the `Array` or `RawArray` class to share a
Numpy array between our worker processes. The difference between an `Array`
and a `RawArray` is that the former offers synchronised (i.e. process-safe)
access to the shared memory. This is necessary if your child processes will be
modifying the same parts of your data.
Due to the way that shared memory works, in order to share a Numpy array
between different processes you need to structure your code so that the
array(s) you want to share are accessible at the _module level_. Furthermore,
we need to make sure that our input and output arrays are located in shared
memory - we can do this via the `Array` or `RawArray`.
As an example, let's say we want to parallelise processing of an image by
having each worker process perform calculations on a chunk of the image.
First, let's define a function which does the calculation on a specified set
of image coordinates:
%% Cell type:code id: tags:
```
import multiprocessing as mp
import ctypes
import numpy as np
np.set_printoptions(suppress=True)
def process_chunk(shape, idxs):
# Get references to our
# input/output data, and
# create Numpy array views
# into them.
sindata = process_chunk.input_data
soutdata = process_chunk.output_data
indata = np.ctypeslib.as_array(sindata) .reshape(shape)
outdata = np.ctypeslib.as_array(soutdata).reshape(shape)
# Do the calculation on
# the specified voxels
outdata[idxs] = indata[idxs] ** 2
```
%% Cell type:markdown id: tags:
Rather than passing the input and output data arrays in as arguments to the
`process_chunk` function, we set them as attributes of the `process_chunk`
function. This makes the input/output data accessible at the module level,
which is required in order to share the data between the main process and the
child processes.
Now let's define a second function which process an entire image. It does the
following:
1. Initialises shared memory areas to store the input and output data.
2. Copies the input data into shared memory.
3. Sets the input and output data as attributes of the `process_chunk` function.
4. Creates sets of indices into the input data which, for each worker process,
specifies the portion of the data that it is responsible for.
5. Creates a worker pool, and runs the `process_chunk` function for each set
of indices.
%% Cell type:code id: tags:
```
def process_dataset(data):
nprocs = 8
origData = data
# Create arrays to store the
# input and output data
sindata = mp.RawArray(ctypes.c_double, data.size)
soutdata = mp.RawArray(ctypes.c_double, data.size)
data = np.ctypeslib.as_array(sindata).reshape(data.shape)
outdata = np.ctypeslib.as_array(soutdata).reshape(data.shape)
# Copy the input data
# into shared memory
data[:] = origData
# Make the input/output data
# accessible to the process_chunk
# function. This must be done
# *before* the worker pool is created.
process_chunk.input_data = sindata
process_chunk.output_data = soutdata
# number of boxels to be computed
# by each worker process.
nvox = int(data.size / nprocs)
# Generate coordinates for
# every voxel in the image
xlen, ylen, zlen = data.shape
xs, ys, zs = np.meshgrid(np.arange(xlen),
np.arange(ylen),
np.arange(zlen))
xs = xs.flatten()
ys = ys.flatten()
zs = zs.flatten()
# We're going to pass each worker
# process a list of indices, which
# specify the data items which that
# worker process needs to compute.
xs = [xs[nvox * i:nvox * i + nvox] for i in range(nprocs)] + \
[xs[nvox * nprocs:]]
ys = [ys[nvox * i:nvox * i + nvox] for i in range(nprocs)] + \
[ys[nvox * nprocs:]]
zs = [zs[nvox * i:nvox * i + nvox] for i in range(nprocs)] + \
[zs[nvox * nprocs:]]
# Build the argument lists for
# each worker process.
args = [(data.shape, (x, y, z)) for x, y, z in zip(xs, ys, zs)]
# Create a pool of worker
# processes and run the jobs.
pool = mp.Pool(processes=nprocs)
pool.starmap(process_chunk, args)
return outdata
```
%% Cell type:markdown id: tags:
Now we can call our `process_data` function just like any other function:
%% Cell type:code id: tags:
```
data = np.array(np.arange(64).reshape((4, 4, 4)), dtype=np.float64)
outdata = process_dataset(data)
print('Input')
print(data)
print('Output')
print(outdata)
```
......
...@@ -514,3 +514,185 @@ print('Non linear registrations:') ...@@ -514,3 +514,185 @@ print('Non linear registrations:')
for t1, result in zip(t1s, nlinresults): for t1, result in zip(t1s, nlinresults):
print(t1, ':', result) print(t1, ':', result)
``` ```
### Sharing data between processes
When you use the `Pool.map` method (or any of the other methods we have shown)
to run a function on a sequence of items, those items must be copied into the
memory of each of the child processes. When the child processes are finished,
the data that they return then has to be copied back to the parent process.
Any items which you wish to pass to a function that is executed by a `Pool`
must be - the built-in
[`pickle`](https://docs.python.org/3.5/library/pickle.html) module is used by
`multiprocessing` to serialise and de-serialise the data passed into and
returned from a child process. The majority of standard Python types (`list`,
`dict`, `str` etc), and Numpy arrays can be pickled and unpickled, so you only
need to worry about this detail if you are passing objects of a custom type
(e.g. instances of classes that you have written, or that are defined in some
third-party library).
There is obviously some overhead in copying data back and forth between the
main process and the worker processes. For most computationally intensive
tasks, this communication overhead is not important - the performance
bottleneck is typically going to be the computation time, rather than I/O
between the parent and child processes. You may need to spend some time
adjusting the way in which you split up your data, and the number of
processes, in order to get the best performance.
However, if you have determined that copying data between processes is having
a substantial impact on your performance, the `multiprocessing` module
provides the [`Value`, `Array`, and `RawArray`
classes](https://docs.python.org/3.5/library/multiprocessing.html#shared-ctypes-objects),
which allow you to share individual values, or arrays of values, respectively.
The `Array` and `RawArray` classes essentially wrap a typed pointer (from the
built-in [`ctypes`](https://docs.python.org/3.5/library/ctypes.html) module)
to a block of memory. We can use the `Array` or `RawArray` class to share a
Numpy array between our worker processes. The difference between an `Array`
and a `RawArray` is that the former offers synchronised (i.e. process-safe)
access to the shared memory. This is necessary if your child processes will be
modifying the same parts of your data.
Due to the way that shared memory works, in order to share a Numpy array
between different processes you need to structure your code so that the
array(s) you want to share are accessible at the _module level_. Furthermore,
we need to make sure that our input and output arrays are located in shared
memory - we can do this via the `Array` or `RawArray`.
As an example, let's say we want to parallelise processing of an image by
having each worker process perform calculations on a chunk of the image.
First, let's define a function which does the calculation on a specified set
of image coordinates:
```
import multiprocessing as mp
import ctypes
import numpy as np
np.set_printoptions(suppress=True)
def process_chunk(shape, idxs):
# Get references to our
# input/output data, and
# create Numpy array views
# into them.
sindata = process_chunk.input_data
soutdata = process_chunk.output_data
indata = np.ctypeslib.as_array(sindata) .reshape(shape)
outdata = np.ctypeslib.as_array(soutdata).reshape(shape)
# Do the calculation on
# the specified voxels
outdata[idxs] = indata[idxs] ** 2
```
Rather than passing the input and output data arrays in as arguments to the
`process_chunk` function, we set them as attributes of the `process_chunk`
function. This makes the input/output data accessible at the module level,
which is required in order to share the data between the main process and the
child processes.
Now let's define a second function which process an entire image. It does the
following:
1. Initialises shared memory areas to store the input and output data.
2. Copies the input data into shared memory.
3. Sets the input and output data as attributes of the `process_chunk` function.
4. Creates sets of indices into the input data which, for each worker process,
specifies the portion of the data that it is responsible for.
5. Creates a worker pool, and runs the `process_chunk` function for each set
of indices.
```
def process_dataset(data):
nprocs = 8
origData = data
# Create arrays to store the
# input and output data
sindata = mp.RawArray(ctypes.c_double, data.size)
soutdata = mp.RawArray(ctypes.c_double, data.size)
data = np.ctypeslib.as_array(sindata).reshape(data.shape)
outdata = np.ctypeslib.as_array(soutdata).reshape(data.shape)
# Copy the input data
# into shared memory
data[:] = origData
# Make the input/output data
# accessible to the process_chunk
# function. This must be done
# *before* the worker pool is created.
process_chunk.input_data = sindata
process_chunk.output_data = soutdata
# number of boxels to be computed
# by each worker process.
nvox = int(data.size / nprocs)
# Generate coordinates for
# every voxel in the image
xlen, ylen, zlen = data.shape
xs, ys, zs = np.meshgrid(np.arange(xlen),
np.arange(ylen),
np.arange(zlen))
xs = xs.flatten()
ys = ys.flatten()
zs = zs.flatten()
# We're going to pass each worker
# process a list of indices, which
# specify the data items which that
# worker process needs to compute.
xs = [xs[nvox * i:nvox * i + nvox] for i in range(nprocs)] + \
[xs[nvox * nprocs:]]
ys = [ys[nvox * i:nvox * i + nvox] for i in range(nprocs)] + \
[ys[nvox * nprocs:]]
zs = [zs[nvox * i:nvox * i + nvox] for i in range(nprocs)] + \
[zs[nvox * nprocs:]]
# Build the argument lists for
# each worker process.
args = [(data.shape, (x, y, z)) for x, y, z in zip(xs, ys, zs)]
# Create a pool of worker
# processes and run the jobs.
pool = mp.Pool(processes=nprocs)
pool.starmap(process_chunk, args)
return outdata
```
Now we can call our `process_data` function just like any other function:
```
data = np.array(np.arange(64).reshape((4, 4, 4)), dtype=np.float64)
outdata = process_dataset(data)
print('Input')
print(data)
print('Output')
print(outdata)
```
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment