Mandelbrot
Mandelbrot
#!/usr/bin/python # -*- coding: utf-8 -*- ## Baustelle: MUSS 1. oder 2. Zeile sein # name mb_seriell # MRies, 09/2015 __author__ = 'ries' import numpy as np import matplotlib.pyplot as plt def mb_iteration(x, y): c = complex(x,y) # x+yj z = 0+0j for iteration in range(255): z = z*z + c if abs( z ) >= 2: # Divergenz return iteration return 255 def mb_matrix(x_min, x_max, y_min, y_max, bild): zeilen_zahl = bild.shape[0] spalten_zahl = bild.shape[1] delta_x = (x_max - x_min) / spalten_zahl delta_y = (y_max - y_min) / zeilen_zahl for spalte_j in range(spalten_zahl): x_pix = x_min + spalte_j * delta_x for zeile_i in range(zeilen_zahl): y_pix = y_min + zeile_i * delta_y bild[(zeilen_zahl-1)-zeile_i, spalte_j] = \ mb_iteration(x_pix, y_pix) return bild bild = np.zeros((400, 400), dtype=np.uint8) mb_matrix(-0.088, -0.085, 0.855, 0.858, bild) plt.imshow(bild) plt.gray() plt.show()
#!/usr/bin/python # -*- coding: utf-8 -*- ## Baustelle: MUSS 1. oder 2. Zeile sein # name mb_thread # MRies, 09/2015 __author__ = 'ries' import numpy as np import matplotlib.pyplot as plt import threading def mb_iteration(x, y): c = complex(x,y) # x+yj z = 0+0j for iteration in range(255): z = z*z + c if abs( z ) >= 2: return iteration return 255 def mb_matrix_t(x_min, x_max, y_min, y_max, bild,threads,i): zeilen_zahl = bild.shape[0] spalten_zahl = bild.shape[1] delta_x = (x_max - x_min) / spalten_zahl delta_y = (y_max - y_min) / zeilen_zahl spalten_start = i * spalten_zahl/threads spalten_ende = (i+1) * spalten_zahl/threads for spalte_j in range(spalten_start, spalten_ende): x_pix = x_min + spalte_j * delta_x for zeile_i in range(zeilen_zahl): y_pix = y_min + zeile_i * delta_y bild[(zeilen_zahl-1)-zeile_i, spalte_j] = \ mb_iteration(x_pix, y_pix) return bild bild = np.zeros((400, 400), dtype=np.uint8) threads = 8 jobs = [] for i in range(0, threads): thread = threading.Thread(target=mb_matrix_t, \ args = (-0.088, -0.085, 0.855, 0.858, bild, \ threads,i) ) thread.start() for i in jobs: i.start() for i in jobs: i.join() plt.imshow(bild) plt.jet() # hsv() # hot() # gray() plt.show()
#!/usr/bin/python # -*- coding: utf-8 -*- ## Baustelle: MUSS 1. oder 2. Zeile sein # name mb_multiprocessing # MRies, 09/2015 __author__ = 'ries' import numpy as np import matplotlib.pyplot as plt import multiprocessing def mb_iteration(x, y): c = complex(x,y) # x+yj z = 0+0j for iteration in range(255): z = z*z + c if abs( z ) >= 2: return iteration return 255 def mb_matrix_t(x_min, x_max, y_min, y_max, bild, procs, p_nr): zeilen_zahl = bild.shape[0] spalten_zahl = bild.shape[1] delta_x = (x_max - x_min) / spalten_zahl delta_y = (y_max - y_min) / zeilen_zahl zeilen_start = p_nr * zeilen_zahl/procs zeilen_ende = zeilen_start + zeilen_zahl/procs for zeile_i in range(zeilen_start, zeilen_ende): y_pix = y_min + zeile_i * delta_y for spalte_j in range(spalten_zahl): # Teilen ueber die zeilen wg python x_pix = x_min + spalte_j * delta_x # fortran besser spaltenbloecke bild[(zeilen_zahl-1)-zeile_i, spalte_j] = \ mb_iteration(x_pix, y_pix) # bild auf dem Kopf return bild if __name__ == '__main__': bild = np.zeros((400, 400), dtype=np.uint8) # 400 400 procs = 4 pool = multiprocessing.Pool(processes = procs) results = [ pool.apply_async(mb_matrix_t, \ #args = (-2.0, 1.0, -2.0,1.5, bild, \ args = (-0.088, -0.085, 0.855, 0.858, bild, \ procs,p_nr) ) for p_nr in range(procs) ] pool.close() pool.join() bb = [ p.get() for p in results ] bild = bb[0] # np.zeros((400, 400), dtype=np.uint8) for i in range(1, procs): bild = bild + bb[i] # besser: einzelne Zeilen-Scheiben plt.imshow(bild) plt.jet() # jet() # hsv() # hot() # gray() plt.show()
ries@ceti:~/mympi/p> cat hostfile ceti slots=4 pallas.uni-trier.de slots=4 puppis.uni-trier.de slots=4 omega.uni-trier.de slots=24
#mr:>FW_SERVICES_EXT_TCP="" FW_SERVICES_EXT_TCP="1024:65535" #mr:>FW_TRUSTED_NETS="" FW_TRUSTED_NETS="136.199.234.1"
# MR , Version g6, 12.06.2015
from mpi4py import MPI
import bmp
import numpy as np
import matplotlib.pyplot as plt
welt = MPI.COMM_WORLD
nummer = MPI.COMM_WORLD.Get_rank()
name = MPI.Get_processor_name()
status = MPI.Status()
#
spalten = 1024 # x, Realteil
zeilen = 1024 # y, Imaginaerteil
MAX_ITER= 1025 # 256
RE_MIN = -0.088 #-2.0
RE_MAX = -0.085 # 2.0
IM_MIN = 0.856 # -2.0
IM_MAX = 0.859 # 2.0
job_anzahl = np.zeros( welt.Get_size(), dtype=np.int16 )
#interactive:>plt.ion()
if welt.Get_rank() == 0: # Master
# init Bild
bild = np.zeros( (zeilen,spalten) )
#sehr langsam:>plt.imshow( bild , vmin=10, vmax=555 , cmap='gray')
#plt.show()
aktive_sklaven = 0
naechste_zeilen_nr = 0
# los geht's: Auftraege an alle Sklaven
for sklave in range( 1, welt.Get_size() ):
welt.send( naechste_zeilen_nr, dest = sklave )
naechste_zeilen_nr = naechste_zeilen_nr + 1 # naechste Zeile
aktive_sklaven = aktive_sklaven + 1
# Master-Schleife
while aktive_sklaven > 0:
puffer = welt.recv( source=MPI.ANY_SOURCE, status=status )
zeile, bildzeile = puffer
sklaven_id = status.Get_source() # so erfahre ich woher das Ergebnis kommt
#print(sklaven_id)
job_anzahl[sklaven_id] += 1
# speicher die Bildzeile im Bild
bild[zeile] = bildzeile
#sehr langsam:>plt.imshow( bild , vmin=10, vmax=555 , cmap='gray')
#plt.draw()
if naechste_zeilen_nr < zeilen:
# neuer Auftrag an diesen Skalven
welt.send( naechste_zeilen_nr, sklaven_id)
naechste_zeilen_nr += 1
else:
# schicke negative Zeilennummer (-1), um Ende zu symbolisieren
welt.send( -1, sklaven_id )
aktive_sklaven = aktive_sklaven - 1
# Statistik
for sklaven_id in range( 1, welt.Get_size() ):
print("{0:5d} Arbeiten durch Sklave {1:3d} erledigt".format( job_anzahl[sklaven_id], sklaven_id ) )
# plot nur am Ende, da sonst SEHR langsam
plt.imshow( bild , vmin=10, vmax=888 , cmap='gray')
plt.show()
else:
# Sklave
while True:
# welche Zeile
zeilennummer = welt.recv( source=0 )
# noch weiter, oder fertig
if zeilennummer < 0:
break
# berechne Bildzeile
bildzeile = np.zeros( spalten) #
# 1j = imaginaere Einheit
c_im = ( IM_MAX - IM_MIN ) *1j * zeilennummer / zeilen + IM_MIN * 1j
for x in range( spalten ):
# und hier der Realteil
c = c_im + ( RE_MAX - RE_MIN ) * x / spalten + RE_MIN
z = 0
for iter in range(MAX_ITER):
z = z * z + c
if abs(z) >= 2.0:
break
bildzeile[x] = iter - 1
# sende Ergebnis an Master
welt.send( [zeilennummer, bildzeile], 0 )
Graph
#!/usr/bin/python # -*- coding: utf-8 -*- ## MUSS in 1. oder 2. Zeile # name plot_1.py import numpy as np import matplotlib.pyplot as plt t = np.linspace(-1.0, 2.0, 10) y = 2*t**3 - 3*t**2 + 1 plt.plot(t, y,'*') plt.show()
#!/usr/bin/python
# -*- coding: utf-8 -*- ## MUSS in 1. oder 2. Zeile
# name plot_2.py
from sympy import * # also matplotlib drin
var('t')
plot( 2*t**3 - 3*t**2 + 1 , (t,-1, 2) )
Räuber-Beute
#!/usr/bin/python
# -*- coding: utf-8 -*- ## Baustelle: MUSS 1. oder 2. Zeile sein
# name myrb2_2.py # MRies, 06/2015
import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate
from matplotlib.widgets import Slider
def mydu(u, t, w,fminus,s, fplus ): # t, da odeint das braucht
# Paramater wie in args=(.. ,)
return np.array( [ w * u[0] - ( fminus * u[1] ) * u[0] ,
-s * u[1] + ( fplus * u[0] ) * u[1] ])
w = 3; # Wachstumsrate
fminus = 2; # Gefressen-Werden-Rate
s = 1; # Sterberate
fplus = 0.5; # Fressrate
t = np.linspace(0,10,300) # 100 Punkte
u0 = np.array([1.1, 0.5] ) # Startwerte
u = integrate.odeint(mydu, u0, t, args=(w,fminus,s, fplus,) )
# args, am Ende ein Komma
fig = plt.figure(figsize=(7,10)) # width, heigt in inch
ax_1 = plt.subplot(3,1,1)
ax_1.set_ylim([0.0,10.0])
u1g, = plt.plot(t, u[...,0],'g+') # erste Spalte
u2g, = plt.plot(t, u[...,1],'rx') # zweite Spalte
ax_2 = plt.subplot(3,1,2)
ax_2.set_ylim([0.0,5.0])
u1u2g, = plt.plot(u[...,0], u[...,1],'bx') # Phasenbild
#plt.subplot(3,1,3)
w_axis = plt.axes([0.15,0.3,0.75,0.03]) # Rechteck
w_slider = Slider(w_axis, 'w', 0,10, valinit=3)
fminus_axis = plt.axes([0.15,0.25,0.75,0.03]) # Rechteck
fminus_slider = Slider(fminus_axis, 'fminus', 0,10, valinit=2)
s_axis = plt.axes([0.15,0.2,0.75,0.03]) # Rechteck
s_slider = Slider(s_axis, 's', 0,10, valinit=1)
fplus_axis = plt.axes([0.15,0.15,0.75,0.03]) # Rechteck
fplus_slider = Slider(fplus_axis, 'fplus', 0,10, valinit=0.5)
def update(val): # muss val heißen
#global w, fminus, s, fplus # sonst in mydu nicht geändert
w = w_slider.val # alternativ als args=(.. ,)
fminus = fminus_slider.val # in odeint-Aufruf
s = s_slider.val
fplus = fplus_slider.val
u = integrate.odeint(mydu, u0, t, args=(w,fminus,s, fplus,) )
u1g.set_ydata(u[...,0])
u2g.set_ydata(u[...,1])
u1u2g.set_data(u[...,0],u[...,1])
plt.draw()
w_slider.on_changed(update)
fminus_slider.on_changed(update)
s_slider.on_changed(update)
fplus_slider.on_changed(update)
plt.show()