Commit a8d50e76 authored by Jochen Schulz's avatar Jochen Schulz
Browse files

fix parallel computation

parent 34361997
Pipeline #225291 failed with stage
in 13 minutes and 48 seconds
......@@ -25,6 +25,6 @@ launch_buttons:
execute:
execute_notebooks : cache # Whether to execute notebooks at build time. Must be one of ("auto", "force", "cache", "off")
cache : "_build/.jupyter_cache/" # A path to the jupyter cache that will be used to store execution artifacs. Defaults to `_build/.jupyter_cache/`
exclude_patterns : ["README.md", "*.ipynb", "parallelisierung_lecture.md", "parallelisierung_exercise.md"] # A list of patterns to *skip* in execution (e.g. a notebook that takes a really long time)
exclude_patterns : ["README.md", "*.ipynb"] # A list of patterns to *skip* in execution (e.g. a notebook that takes a really long time)
stderr_output: show
timeout: 240
......@@ -69,7 +69,7 @@ CD
```
```python cell
```{code-cell}
:tags: ['hide-cell']
print('numpy+numba+dask')
......
......@@ -4,8 +4,10 @@ jupytext:
text_representation:
extension: .md
format_name: myst
format_version: 0.13
jupytext_version: 1.10.3
kernelspec:
display_name: Python 3
display_name: Python 3 (ipykernel)
language: python
name: python3
---
......@@ -28,7 +30,7 @@ Im wissenschaftlichen Rechnen beschäftigt man sich eher mit *Parallelisierung*.
Array-Berechnungen in NumPy werden von C- oder Fortran-Backends durchgeführt. Je nach Backend und Operation sind diese automatisch parallelisiert, z.B. die `dot`-Funktion bzw. der `@`-Operator oder diverse `linalg`-Routinen.
```{code-cell}
```{code-cell} ipython3
import numpy as np
A = np.random.random((1000, 1000))
......@@ -40,13 +42,14 @@ for _ in range(10000):
Elementweise Produkte sind allerdings nicht parallelisiert, da hier der limitierende Faktor oft nicht die CPU, sondern der Speicherzugriff ist.
```{code-cell}
```{code-cell} ipython3
x = np.random.random((100000,))
y = np.random.random((100000,))
for _ in range(100000):
z = x * y
```
(numba)=
## Numba vectorize parallel
......@@ -56,7 +59,7 @@ Parallelisierung auf einem Rechner durchführen. Dazu müssen wir nur das keywor
Nehmen wir das Beispiel des Mandelbrots aus der Vorlesung über Performance.
```{code-cell}
```{code-cell} ipython3
import numpy as np
import numba as nu
......@@ -66,14 +69,14 @@ def mandel_vec(c):
maxiter = 40
for n in range(maxiter):
z = z**2 + c
if abs(z) > 2:
if np.abs(z) > 2:
return n
return maxiter
```
Koordinaten wählen; dem Array `C` geben wir hier explizit seinen Typ `complex64` mit:
```{code-cell}
```{code-cell} ipython3
# Koordinaten
N = 2000
X, Y = np.ogrid[-1:1:N*1j, -1:1:N*1j]
......@@ -82,27 +85,28 @@ C = np.array(X + 1j*Y, dtype=np.complex64)
Dann führen wir erst die sequentielle Variante von `mandel_vec` aus und messen die Laufzeit.
```{code-cell}
```{code-cell} ipython3
print("mandel_vectorize:")
%timeit mandel_vec(C)
```
Und dann die parallele:
```{code-cell}
```{code-cell} ipython3
@nu.vectorize(["int8(complex64)"],target='parallel')
def mandel_vec_par(c):
z = 0
maxiter = 40
for n in range(maxiter):
z = z**2 + c
if abs(z) > 2:
if np.abs(z) > 2:
return n
return maxiter
print("mandel_vectorize_parallel:")
%timeit mandel_vec_par(C)
```
(dask)=
## Mächtiges Parallelisierungs-Modul Dask
......@@ -115,7 +119,7 @@ Dazu gibt es ausgiebige debug und Analyse-Möglichkeiten.
Wir starten mit dem `distributed` Client von Dask, der aber auch gut auf einem
Rechner arbeiten kann.
```{code-cell}
```{code-cell} ipython3
from dask.distributed import Client
```
......@@ -123,7 +127,7 @@ Wenn wir Bokeh installiert haben, startet das Ausführen des Clients
ein diagnostisches dashboard unter http://localhost:8787 (auf dem Jupyter-server ist dies leider nicht
sinnvoll möglich anzusehen).
```{code-cell}
```{code-cell} ipython3
client = Client()
client
```
......@@ -142,7 +146,7 @@ dask.array(*args, **kwargs, chunks=None)
Beispiel:
```{code-cell}
```{code-cell} ipython3
import dask.array as da
x = da.random.uniform(low=0, high=10, size=(10000, 10000), # normal numpy code
chunks=(1000, 1000)) # break into chunks of size 1000x1000
......@@ -151,13 +155,12 @@ y = x + x.T - x.mean(axis=0) # Use normal syntax for high level algorithms
y.sum().compute()
```
### Dask Beispiel: Mandelbrot-Menge
Wir führen unser Mandelbrot-Beispiel fort; so haben wir auch die besten Vergleiche.
Wir erzeugen nochmal die Koordinaten.
```{code-cell}
```{code-cell} ipython3
# Koordinaten
N = 2000
X, Y = np.ogrid[-1:1:N*1j, -1:1:N*1j]
......@@ -167,7 +170,7 @@ C = np.array(X + 1j*Y, dtype=np.complex64)
Dann erzeugen wir ein dask array mit gegebenen chunks. Dask zeigt uns dann
auch noch Meta-Informationen über das Array und seine Aufteilung.
```{code-cell}
```{code-cell} ipython3
import dask.array as da
CD = da.from_array(C, chunks=(-1,1000))
CD
......@@ -182,7 +185,7 @@ dask.array.map_blocks(func)
welches die gegebene Funktion `func` mit den verschiedenen chunks parallel füttert.
```python cell
```{code-cell} ipython3
@nu.njit
def mandel_array(C):
Z = np.zeros_like(C)
......@@ -212,24 +215,26 @@ mit Dask verwenden können.
Die Ausführung dieser Funktion macht Schwierigkeiten mit dem distributed scheduler. Daher stoppen wir den leider vorerst einmal,
bis eine Lösung des Problems gefunden ist.
```{code-cell}
```{code-cell} ipython3
client.close()
```
```{code-cell}
```{code-cell} ipython3
@nu.vectorize(["int8(complex64)"])
def mandel_vec_dask(c):
z = 0
maxiter = 40
for n in range(maxiter):
z = z**2 + c
if abs(z) > 2:
if np.abs(z) > 2:
return n
return maxiter
print("mandel_vec dask array:")
%timeit mandel_vec_dask(CD).compute()
```
(tipps_parallel)=
### Tipps zur Optimierung
......@@ -255,14 +260,14 @@ print("mandel_vec dask array:")
## Bonus: guvectorize normal
Alternative generalisierte `ufunc`, die beliebige Eingabe- und Ausgabe-Array-Typen sowie -Größen haben kann.
Alternative generalisierte `ufunc`, die beliebige Eingabe- und Ausgabe-Array-Typen sowie -Größen haben kann. `vectorize` kann nur mit Skalar-Funktionen hantieren.
Hier wird in der Typ-Annotation noch die Größe der Arrays mithilfe der üblichen
shape-Schreibweise angegeben. Beispielhaft: Hat eine Funktion von einem array `(n,m)`
als Rückgabewert ein Skalar `()`, dann schreibt man `(n,m) -> ()`. Aufgrund der
generalisierten Rückgabe gibt es *kein* return mehr, sondern die Rückgabe geschieht
über das letzte Argument der Funktion und wird innerhalb der Funktion gesetzt.
```python cell
```{code-cell} ipython3
@nu.guvectorize('complex64,int8[:]','()->()')
def mandel_guvec(c, out):
z = 0
......@@ -281,7 +286,10 @@ print("mandel_guvec normal:")
Auch dies kann man direkt mit einem Dask-Array füttern:
```python cell
```{code-cell} ipython3
print("mandel_guvec dask array:")
%timeit mandel_guvec(CD).compute()
```
```{code-cell} ipython3
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment