Modellierung und Simulation 1 - fh-swf.de · • Ergänzen Sie das Programm um die Möglichkeit des...

Modellierung und Simulation 1

Prof. Dr. H. Schulze

Prof. Dr. C. Lüders

Modellierung und

Simulation 1

Kapitel 4

Modellierung und Simulation

von Bediensystemen

Prof. Dr. Henrik Schulze, Prof. Dr. Christian Lüders

Fachhochschule Südwestfalen

Standort Meschede

Fachbereich Ingenieur- und Wirtschaftswissenschaften

Tel.: 0291 / 99 10 -4300 | -4261

E-Mail: schulze.henrik | [email protected]




4.1 Allgemeines




Bediensystem

Quellen

Sources Warteschlange

Queue

Bedieneinheit

Server

Senke

Sink

Zeit

zufällige Startzeiten

zufällige Bedienzeiten




Bediensysteme – Typische Fragen

1. Welcher Anteil der Anforderungen wurde erfolgreich bedient?

2. Wie groß war die gesamte Bedienzeit? (im Mittel, in 90% der Fälle)

3. Wie groß war die Wartezeit? (im Mittel, in 90% der Fälle)

4. Wie viele Server werden benötigt?

5. Wie viele Warteplätze werden benötigt?

6. Was ist die geeignete Abfertigungsmethode? (FIFO, Priorität)

7. Lassen sich die einzelnen Prozesse besser anordnen? (parallel)

8. Wie lässt sich die „Kundenzufriedenheit“ steigern?

9. Wie stark sind die einzelnen Server ausgelastet?

10. …




Bediensysteme - Beispiele

1. Postschalter

2. Funkkanalzuteilung

3. Zugriffe auf ein BUS-System

4. Dimensionierung & Abfertigungsstrategien bei Routern

5. Dimensionierung von Servern

6. Lagerhaltung, Bestückung von Automaten

7. Geschäftsprozesse

8. Produktionsprozesse

9. …




R, tX t

Menge „gleichartiger“ Zufallsvariablen, die i.Allg. abhängig voneinander sind

Index wird vielfach als Zeit t (diskretisierte Zeit n) interpretiert

Wertebereich der Zufallsvariablen: Zustandsraum

0N, nX ndiskret: kontinuierlich

Markov-Kette Markov-Prozess

Spezialfall:

• Zustandsraum diskret

• Gedächtnislosigkeit

• Memoryless (M)

Stochastischer Prozess




Beispiele für stochastische Prozesse

1. Anzahl der Sonnenstunden in Meschede am Tag n im Jahr

2. Globale Durchschnittstemparatur im Jahr n

3. Tabellenplatz von Schalke 04 am Spieltag n

4. Position des Fußballs im Spiel zur Zeit t

5. Buchstabe in einem Text an Position

6. Windgeschwindigkeit an einem Standort zur Zeit t

7. Anzahl laufender Telefonate in einer Funkzelle zur Zeit t

8. Anzahl von insgesamt abgewickelter Telefonate am Tag n

9. Anzahl der Zugriffe auf einen Server

10. Aufgerufene Internetseite in einer Session nach Klick n

11. Bestand eines bestimmten Bauteils in einem Lager am Tag n

12. Ausfall eines Bauteils zur Zeit t




4.2 Markov-Ketten




)(

)|(BP

BAPBAP

Wahrscheinlichkeit, dass A eintritt unter der Bedingung, dass B vorliegt.

)(

)(

)(

)()|( AP

BP

BPAP

BP

BAPBAP

Falls A und B unabhängig sind:

Bedingte Wahrscheinlichkeit:

Wiederholung: Bedingte Wahrscheinlichkeiten

1 2 3 4 5 6

1 2 3 4 5 6 7

2 3 4 5 6 7 8

3 4 5 6 7 8 9

4 5 6 7 8 9 10

5 6 7 8 9 10 11

6 7 8 9 10 11 12

Bsp: Würfeln mit 2 Würfeln

A: Augenzahlsumme 7

12/736/21)( AP

B: erster Würfel zeigt 6

6/1

6/1

)()|(

BP

BAPBAP

B: erster Würfel gerade

3

2

2/1

36/12

)()|(

BP

BAPBAP




A

)(

)|(BP

BAPBAP


)(

)(

)(

)()|( AP

BP

BPAP

BP

BAPBAP

Falls A und B unabhängig sind:


)()|( BPBAPBAP

i

ii

i

i BPBAPBAPAP )()|(

B1 B2 B3

B4 B5 B6

B7 B8 B9

Wiederholung: Bedingte Wahrscheinlichkeiten

Zerlegung in disjunkte Mengen Bi




Diskreter stochastischer Prozess mit diskretem Zustandsraum und folgender Eigenschaft

für die bedingten Wahrscheinlichkeiten:

)|(),,...,,( 1100112211 nnnnnnnnnn jXjXPjXjXjXjXjXP

Nur der unmittelbar vorhergehende Zustand spielt eine Rolle, nicht die Vorgeschichte.

Übergangswahrscheinlichkeiten von Zustand i j: pi j

)|()( 1 iXjXPnp nnji

Homogene Kette: Übergangswahrscheinlichkeiten unabhängig von n: pi j

1 2 k 0

10p 21p 32p 1, kkp

1,1kp32p21p10p

00p

Definition einer Markov-Kette

Übergangsgraph

für die Zustände 0, 1, 2, … k, k+1, …




Die Übergangswahrscheinlichkeiten erlauben es, die Wahrscheinlichkeit für den

Zustand j zur Zeit n iterativ zu berechnen:

Zustandswahrscheinlichkeiten zur Zeit n als Zeilenvektor:

Iterative Berechnung der Zustandswahrscheinlichkeiten:

)()()()( 11 iPpiXPpjXPjP n

i

jin

i

jinn

............

...

...

...

222120

121110

020100

ppp

ppp

ppp

TP

...),2(),1(),0( nnnn PPPP

Tnn PP P 1

n

Tn PP P 0

Übergangswahrscheinlichkeiten als Elemente der

Übergangs-Matrix (Transission T):

Stationärer Zustand: Tstst PP P Eigenvektor der Übergangsmatrix z. Eigenwert 1

Wann existiert ein stationärer Zustand? Wann gilt: ? n

Tn

st PP P

0lim

Zeilensumme = ? Warum?

Markov-Kette: Übergangsmatrix u. stationärer Zustand

)|( 1 iXjXPp nnji




Beispiel 1: System mit Ausfall und Reparatur

1 0

r

a 1-a 1-r

rr

aapppp

1

1),(),( 1010

rpapp 100 )1( aprpp 011 )1(

rpap 100 aprp 010

Gerät in

Ordnung

Gerät

defekt Stationärer Zustand

in Matrix-Form

Wahrscheinlichkeit

• für Ausfall an Tag n: a

• für Reparatur an Tag n: r Gleichgewicht

bei 0

GGW

bei 1

p0 a = p1 r Gleichgewicht auf Link

B.1: MATLAB-Übung

[V, D] = eig(A)

D: Diagonalmatrix mit Eigenwerten von A

V: Matrix mit Eigenvektoren von A

als Spalten

Normierung des 1. Eigenvektors (Warum?)

VN = V(:, 1) / sum( V(:, 1) )

• Berechne: Eigenvektor

• Berechne: PT^n, n = 1, 2, 4, 8, …

01 pr

ap

00101 pr

appp

ra

rp

0

ra

ap

1

Normierung:

Ergebnis:

Vergleich für 2 Parametersätze




Beispiel 2: Erweiterte Systeme mit Reparatur

1 0

r1

a1

1- a1 – a2

3 2

1 0

r

a 1-a

Gerät in

Ordnung

Gerät

defekt

2

1 - r

l

1 - l

Nicht

reparabel

a2 r2

Wahrscheinlichkeit

• für Ausfall an Tag n: a

• für Reparatur an Tag n: r

• für Neubeschaffung am Tag n: l

r1

a1

a2 r2

1- r1 – a2

1- r1 – r2 1- a1 – r2

Ausfallwahrscheinlichkeit

• falls beide Teile erforderlich: 1 – p0

• falls Teile redundant sind: p3

Ausfall eines Geräts aus zwei Teilen Nicht-reparabler Defekt




Bediensystem mit einem Server und Warteschlange

• Diskreter Zeittakt

• unbegrenzte Warteplätze

• Zustand:

• Anzahl Kunden k im System

• Neuankunft Kunde: pa = l

• Bedienzeit-Ende: pe = µ

• Last: a = l / m < 1

1 Server, 1 Warteplatz

1 2 0

m m m

3

m

l l l l

1lm 1lm 1lm

1l

l pk = µ pk+1 a pk = pk+1

pk = ak p0 p0 = 1 – a

1 2 0

m m

l l

1lm

1l 1m

Mittlere Kundenzahl im System: E{K} = a/(1 – a)

B.2 MATLAB-Übung

• Übergangsmatrix PT aufstellen

• Eigenwerte / -vektoren berechnen

• PT^n berechnen

• Wartewahrscheinlichkeit?

• Serverauslastung?

• Mittlere Kundenzahl




Explizite Simulation einer Markov-Kette (Bsp. M/M/1) nmax = 100000; % Maximale Laufzeit

lambda = 0.06; % Erzeugungswahrscheinlichkeit

mu = 0.1; % Ende-Wahrscheinlichkeit

a = lambda/mu; % Last

k = 0; % aktuelle Zahl: Kunden im System

K = []; % speichert Kunden in jedem Zeitschritt

N = 0; % aktuelle Anzahl erzeugter Kundenwünsche

for n = 1:nmax, % Schleife über Zeitschritte

if(rand(1,1) < lambda + mu), % Änderung der Kundenzahl

% mit p = lambda + mu

if( rand(1, 1) > a /(1 + a) ) % Verringerung Kundenzahl (Ende)

k = max(0, k - 1);

else % Erhöhung der Kundenzahl

k = k + 1;

end

end

K = [K; k]; % Füge aktuelle Kundenzahl zu

end % Kundenzahlvektor hinzu

Statistische Auswertung, Diagramme




Zeitlicher Verlauf der Kundenzahl und Histogramm

der Zustandswahrscheinlichkeiten für M/M/1-System

0 500 1000 1500 20000

1

2

3

4

5

6Kundenzahl im System als Funktion der Zeit

Zeitschritt

Kundenzahl im

Syste

m

0 2 4 6 8 100

0.5

1

1.5

2

2.5

3

3.5

4x 10

4 Histogramm der Kundenzahl im System

Kundenzahl

Häufigkeit

Parameter wie auf vorheriger Folie

zeitlicher

Mittelwert




Aufg. B.3: Simulationsaufgaben zur Markov-Kette (M/M/1)

• Wie groß muss man nmax wählen, um ca. 20000 Ankünfte zu simulieren?

• Ergänzen Sie das Simulationsprogramm diskretMM1 um

– die Ausgabe der tatsächlichen Ankünfte

– ein Diagramm zum zeitlichen Verlauf der Kundenzahl sowie mit dem Histogramm

– die Ausgabe des Ensemble-Mittelwerts E{K} und des zeitlichen Mittelwerts der Kundenzahl

• Simulieren Sie das System für zwei verschiedene Werte der Last und vergleichen

Sie die Ergebnisse mit den Werten für die Formeln zum stationären Zustand.

• Simulieren Sie das System zweimal mit der Last a = 1! Was stellen Sie fest?

• Ergänzen Sie das Programm um die Möglichkeit des Abweisens (Blockierung B)

eines Kundenwunsches bei einer endlichen Warteschlangenlänge (z.B. nur 1 Platz)!

• Bauen Sie einen Zähler für die Anzahl von Blockierungen ein und geben Sie die die

Blockierwahrscheinlichkeit aus.

• Simulieren Sie das System für zwei verschiedene Werte der Last und vergleichen

Sie die Ergebnisse mit den Werten für die Formeln zum stationären Zustand.




Weitere Beispiele für Markov-Ketten

0 1

q

p

3 2

p

q

q q p p

01 11

00 10

21

20

02 12 22

0

3 4

0,5

0,5

1 2 0,5

0,5

0,75

0,25

0,5 0,5

0,5 0,5

B.4 MATLAB-Übung

• Erzeuge PT

• Berechne: PT^n,

• n = 50, 51, 52, 53

• Was fällt auf?

• Woran liegt das?

• ….

Zweidimensionaler Prozess

„Random Walk“

Überlegung

• Was geschieht langfristig?

• Pn im stationärem Zustand?

• ….




Internet-Suche als Markov-Kette

Einfaches Beispiel mit 5 verlinkten Seiten

5

2 4

3

1

p1j = (0, ¼, ¼ , ¼, ¼)

p2j = (½ , 0, 0 , ½, 0)

p3j = (0, 0, 0 , 1, 0)

p4j = (0, 1, 0 , 0, 0)

p5j = (1/3 ,1/3, 0 ,1/3, 0)

Übergangswahrscheinlichkeiten

aus dem Zustand i heraus: 1/L

L: Anzahl der Links

Wie sieht die stationäre Verteilung aus?

Ranking der Seiten nach

Wahrscheinlichkeit (stationär)

B.5a: MATLAB-Übung

• Übergangsmatrix PT einlesen

• liegt als Excel-Tabelle vor

• PT = xlsread(`www5.xlsx´)

• PT^n berechnen




Internet-Suche als Markov-Kette

Einfaches Beispiel mit 6 verlinkten Seiten

5

2 4

3

1

p1j = (1- a) p1j + a / N

Problem?

6 Absorbierende

Zustände

Prozess läuft irgendwann in die

absorbierenden Zustände 3 und 6

und kehrt nicht mehr zurück

Lösung?

„Steuersatz“ a bei allen Zuständen

Link p34 gestrichen

Jetzt: immer ein Link zurück

B.5b: MATLAB-Übung

• Übergangsmatrix PT einlesen

• PT = xlsread(`www6.xlsx´)

• PT^n berechnen

• Mit Steuersatz modifizieren




Für eine (positiv) rekurrente, aperiodische und irreduzible Markov-Kette gilt:

- Der folgende Grenzwert für die Übergangsmatrix existiert:

)(/1 kPstk

- Es existiert ein eindeutiger stationärer Zustand (Eigenvektor): Tstst PP P

Wichtige Ergebnisse für Markov-Ketten

n

Tn

T PP

lim

- Alle Zeilen von sind identisch und gleich

TP stP

- Mittlere Zeit für die Rückkehr in den Zustand k:

- Das Verhalten der Markov-Kette ist ergodisch:

Erwartungswert/Ensemble-Mittel = zeitliches Mittel

k

st kPkKE )(

N

n

nN

t KN

K1

1lim

z.B. Mittlere

Kundenzahl

Gilt auch für andere Größen.




4.3 Markov-Prozesse

(und ähnliche Prozesse)

- zeit-kontinuierlich




Markov-Kette (diskrete Zeitschritte)

)|(),,...,,( 1100112211 nnnnnnnnnn jXjXPjXjXjXjXjXP

Nur der unmittelbar vorhergehende Zustand spielt eine Rolle, nicht die Vorgeschichte.

Übergangswahrscheinlichkeiten von Zustand i j: )|( 1 iXjXPp nnji

Markov-Prozess und Übergangsraten (1)

)()|()( iXPiXjXPjXP t

i

ttttt

)(

)|()()(iXP

t

iXjXP

t

jXPjXPt

i

jitttttt

)()( tPtdt

dPi

i

ji

j

Markov-Prozess (Zeit kontinuierlich)

0t

)()()()( 11 iPpiXPpjXPjP n

i

jin

i

jinn 1 Zeilensumme j

jip

jjji

jii

ji PP

j

jzz

zji

jii

ji PP

0z

zi

Übergangsraten: ij




in Zustand j

hinein

aus Zustand j

heraus

Markov-Prozess und Übergangsraten (2)

)()|()( iXPiXjXPjXP t

i

ttttt

)(

)|()()(iXP

t

iXjXP

t

jXPjXPt

i

jitttttt

dt

dPj

Übergangsraten i j : ij

0t

j

jzz

zji

jii

ji PP

0eZeilensumm0 j

ji

„Wahrscheinlichkeitsfluss“

0: dt

dPjZustandrstationäre




)(

)|(BP

BAPBAP


)(1

)()|],[(

tP

ttpüberlebttZeitzurbistttinAusfallPP

T

Tausfall

)(1

)(lim)(

0 tP

tp

t

Pt

T

Tausfall

t

a

z.B. Wahrscheinlichkeit, dass Bauteil bis zur Zeit t ausfällt:

Verteilung: F(t) = PT( T ≤ t), Dichte: pT(t) = f(t) = F‘(t)

Ausfallrate:


Durch Reparatur bzw. Ausfall eines Bauteils / Ende oder Beginn eines Prozesses

ändert sich der Zustand eines Systems von i auf j = i ± 1

Bestimmung der Übergangsrate:

Bedingte Wahrscheinlichkeiten u. Ausfallrate / Ende-Rate




)(1

)()(

tP

tpt

T

T

aAusfallrate:

)exp()()( tµµtptf T

Exponential-Verteilung mit Parameter µ:

Dichte:

)exp(1)()( tµtPtF T

µµt

µtµ

tP

tpt

T

T

))exp(1(1

)exp(

)(1

)()(a

Ausfallrate

Sterberate

Ende-Rate

Bei einer Exponential-Verteilung ist die Ausfallrate konstant (µ) und hängt nicht

von der Vorgeschichte ab – Gedächtnislosigkeit (Memoryless)

Mittlere Lebensdauer (Mean Time Before Failure): MTBF = E{T} = = 1/µ

Analog bei Zwischenankunftszeiten/ Reparaturzeiten gemäß Exponential-Verteilung:

Konstante Ankunftsrate l / Reparaturrate r : r = 1/MTTR (Mean Time To Restore)

Ausfallrate: Gedächtnislosigkeit Exponential-Verteilung




)(1

)(lim)(

0 tP

tp

t

Ptr

T

Tausfall

t

Ausfallrate:

Weibull-Verteilung mit Parameter µ, k:

Dichte:

k

T tµtPtF exp1)()(

1

)(1

)()(

k

T

T tµµktP

tptr

Ausfallrate

Sterberate

Enderate

Einfluss von k ?

kk

T tµtµµktptf

exp)()(1

Ausfallrate (Weibull-Verteilung)

t

r Badewannenkurve

Früh-

ausfälle Normalbereich Alterung

k = 1

k > 1 k < 1

Ausfallrate bei Pareto?




Kendall-Notation für Bediensysteme A / B / s / q / K / SD

• A: Verteilung der Zwischen-Ankunftszeiten

• M: Memoryless (Exponential) D: Deterministisch

• G: General (beliebige Verteilung) P: Pareto ….

• B: Verteilung der Bedienzeiten

• Siehe A

• s: Anzahl der Bedieneinheiten 1, 2, 3, …

• q: Anzahl der Bedieneinheiten + Warteplätze s, s+1, s+2, …,

• K: maximale Anzahl der Kunden

• K = 1, 2, 3, …,

• SD: Sevice-Disziplin

• FIFO: First In First Out LIFO: Last In Last Out

• SJF: Shortest Job First PRIO: Prioriry …

Rot: Default-Werte

können weggelassen

werden, z.B.

• M/M/1 =

• M/M/1///FIFO

• M/M/8/8 =

• M/M/8/8//FIFO




Mehrere Bedienplätze s, unbegrenzte Warteplätze: M/M/s

1 2 s 0

l l l l l l

1m 2m 3m sm sm sm

l

sm

k

l l

km k1m

Kunden im System: k

Anzahl Server: s

Ankunftsrate: l

Mittlere Bearbeitungszeit:

Ende-Rate: µ = 1/

Gesamtlast A = l/µ = l

Last pro Server a = A/s

veränderlich

8

60 h-1

6 min = 0,1 h

10 h-1

60 h-1 0,1 h = 6 Erl

0,75 Erl pro Server

Parameter Zahlenbeispiel (kleines Callcenter)

0,95 Erl pro Server

Pseudo-Einheit

Erlang (Erl)

Agner Krarup Erlang

Dän. Mathematiker u. Ingenieur




Mehrere Bedienplätze s, unbegrenzte Warteplätze: M/M/s

1 2 s 0

l l l l l l

1m 2m 3m sm sm sm

l

sm

l pk + kµ pk = l pk-1 + (k+1)µ pk+1

l pk + sµ pk = l pk-1 + s µ pk+1

kµ pk = l pk-1

l p0 = 1µ p1

sµ pk = l pk-1

1 kk pk

Ap

0!

pk

Ap

k

k

1 kk pap

0!

pas

Ap sk

s

k

k = 1, 2, …., s – 1,

k = s, s+1, ….

1

11

0

0 )1(!!

a

s

A

k

Ap

ss

k

k

k

l l

km k1m Kunden im System: k

Anzahl Server: s



Gleichgewicht bei Gleichgewicht bei

Normierung

Mittlere Länge

Warteschlange

Mittlere Anzahl

Kunden im System

02)1(!p

a

a

s

AK

s

W

WKAK




M/M/s/s-System, keine Warteplätze Verlustsystem

1 2 s 0

l l l

1m 2m 3m

l

sm

0!

pk

Ap

k

k

1

0

0!

s

k

k

k

Ap

k

l l

km k1m

Kunden im System: k

Anzahl Server: s


0!

)|( ps

ApAnkunftskPB

s

s

Gesamtlast A = l / µ

Ankunftsrate l

Serverauslastung: A (1 – B)

Eintrittsrate: le l (1 – B)

System

Zustandswahrscheinlichkeiten für s ???

Blockierwahrscheinlichkeit (Erlang-B-Formel)

10-2

10-1

100

10-1

100

101

102

Last pro Server

Blo

ckin

g B

/ %

Neuankunft

Blockierung

s = 2 4 8 16 32

s =1




Poisson-Prozess

Ak

k ek

Ap

!Wahrscheinlichkeit für k Kunden:

System ohne Begrenzungen:

Zwischenankunftszeiten exponential-verteilt, Bedienzeiten exponential-verteilt

Anzahl Kunden Kt im System zur Zeit t ist ein Poisson-Prozess

tKEµ

A l

0 5 10 150

0.1

0.2

0.3

0.4

0.5

k

pk

0 5 10 150

0.1

0.2

0.3

0.4

0.5

k

pk

0 5 10 150

0.1

0.2

0.3

0.4

0.5

k

pk

0 5 10 150

0.1

0.2

0.3

0.4

0.5

k

pk

A = 0,75 A = 1,5

A = 3 A = 6




Theorem von Little Verweilzeiten im System

T1

T2

T3

T4

T5

T6

T7

T8

Ka

t1 t2 t3 t4

11

aK

j

jT

Ke

T

<<

GFTeK

j

j

1

1

I

i

ii tkGFFlächeGraue1

I

i

iiT

t tkT

KKE1

1lim

aa K

j

j

a

a

K

j

j TKT

KT

T 11

11

TE l

Ergodischer Prozess (ohne Verluste)

TEKKE t l

Mittlere Kunden-

Zahl im System

Mittlere

Verweilzeit

Ankunfts-

rate = x

t

k

QQ TEKE l

TEBKE l)1(

in der Queue

mit Blockierung:

benötigt Markov-Eigenschaft nicht!!!




Warte- & Verweilzeiten bei M/M/s über Theorem von Little

1 2 s 0

l l l l l l

1m 2m 3m sm sm sm

l

sm

1

11

0

0 )1(!!

a

s

A

k

Ap

ss

k

k

02)1(!p

a

a

s

AK

s

W

WKAK

k

l l

km k1m Kunden im System: k

Anzahl Server: s

Gesamtlast A = l / µ


Normierung

Mittlere Länge

Warteschlange

Mittlere Anzahl

Kunden im System

l/WW KT

WW T

µKAKT

1// ll

Mittlere Wartezeit Mittlere Verweildauer

Verteilung der Wartezeit Verteilung der Verweildauer

TttTP /exp1 WWW TtptTP /exp1




Jackson-Netz: Netz aus M/M/1-Systemen

µ1 µ1

µ1 µ2 µ1

µ2

µ3

µ1

µ2

µ3 l2

l1

l1

l2 l1

l3 l1 l2

l3 ql1

l2 1q)l1

l1

Last an Bedieneinheit i: ai = li / mi

1el

2el

3el

4el

43211 eeee lllll

1el

11 )1( eq ll

q

n

i

i

k

in aakkkp i

1

11 )1(),...,(Zustandswahrsch.

nges KKKK ...21Mittl. Kundenzahl:

emee

gesKT

lll ...21 Mittl. Verweilzeit:

Eintrittsraten der für externe (e) Quellen: ejl

1el1el

1el

2el




S-Bahn-Station

An einer belebten S-Bahn-Station treffen im Mittel

10 Bahnen pro Stunde ein.

Ein Fahrgast trifft zu zufälligen Zeit ein. Wie lange

muss er im Mittel auf eine Bahn warten?

0 1 2 3 t / h

0 1 2 3 t / h

18 min

Mittlere Wartezeit: 0,3 1 min + 0,7 21 min = 15 min




Pollaczek–Khinchine-Formel für M/G/1-System

2

)1( 2

,,

CTT MWGW

Mittlere

Wartedauer

Zwischenankunftszeiten: exponentialverteilt mit Rate l

Bedienzeiten T mit beliebiger Verteilung G = G(T)

• Mittelwert

• Standardabweichung s

• Variationskoeffizient C = s /

• Last a = l

für ein beliebiges G

kein Markov-Prozess

• Index M: M/M/1-System

• Index G: M/G/1-System

2

)1( 2

,,

CKK MWGW

Mittlere Länge

Warteschlange

WG TTMittlere

Verweilzeit

aKK GWG ,

Mittlere

Kundenzahl

Exponentialverteilung

C = ?




Bedien-Disziplinen

ohne Pre-emption mit Pre-emption

Ohne Prioritäten

First In First Out (FIFO) =

First Come First Serve (FCFS)

Last In First Out (LIFO)

Random (RAND)

Priority

Shortest Remaining

Processing Time (SRPT)

Ohne Prioritäten

Processor Sharing (PS)

= Time Div. Multiplex (TDM)

= Round Robin (RR)

Forground Background (FB)

Priority (PRIO)

• Shortest Job First (SJF)

• Job Type

• Transmission Quality

• Proportional Fair

• Nearest Deadline

Pre-emption: Vorkaufsrecht, Verdrängung

i.A. keine analytischen Formeln

Simulation




4.4 Ereignisorientierte Simulation

(von Bediensystemen)

Komplexere Systeme mit beliebigen Bedienzeit-

Verteilungen und komplexen Zuteilungsstrategien

nicht analytisch behandelbar




Ereignisorientierte Simulation

Nachteile von festen, diskreten Zeitschritten (Perioden-Orientierung) bei Zufallsprozessen

• Zeitschritt groß: mehrere „gleichzeitige“ Ereignisse in einem Zeitschlitz

• Zeitschlitz klein: in den meisten Zeitschlitzen geschieht nichts (ineffizient)

Lösung: Ereignisorientierte Simulation (Discrete Event (System) Simulation)

• Nur die Zeitpunkte werden betrachtet, an denen sich der Zustand des Systems ändert

• Ereignis: Zustandsänderung

• Abzuarbeitende Ereignisse werden in eine dynamische EventList eingetragen

• Kennzeichen eines Events: Eintrittszeitpunkt, Typ, zugeordneter Prozess (Kunde, Paket), …

• Abhängig von Typ wird zum Eintrittszeitpunkt bestimmte Routine (Funktion) abgearbeitet

• Je nach Ergebnis können neue Events erzeugt werden (zu zufälligen Zeiten)

• Abzuarbeitende Prozesse / erzeugte Pakete speichern ihren Lebenslauf

• Wann ist was passiert?




Illustration: Ereignisorientierte Simulation

E1 E2 E3 E4 E5 E7 E8 E9 E10 E11

Packet 3

• pnr

• state

• type

• length

• origin

• Statistics

(Lebenslauf)

P 1 P 2 P 4

Event 6

• time

• type

• pnr

Generator

• RandomGenerator1

- verschiedene Param.

Ankunftszeit

• RandomGenerator2

- verschiedene Param.

Paketlänge

• weitere Parameter

Prozessor belegt:

Einreihen in Warteschlange

Prozessor frei:

Start der Bearbeitung

Ende festlegen

Bearbeitung abgeschlossen

Neues Paket in Warteschlange?

Zeit t

Zwischenankunftszeit Bearbeitungszeit Wartezeit

P 3




Startzeiten von Prozessen / Paketen switch Generator.type

case 'deterministic'

NextTime = time + Generator.par1;

case 'uniformly_distributet'

dt = Generator.par2 - Generator.par1;

NextTime = time + Generator.par1 + dt*rand(1, 1);

case 'exponential'

lambda = Generator.par1; % hier nur par1 benötigt

NextTime = time - log(rand(1,1)) / lambda;

case 'weibull'

????????????;

case 'pareto'

????????????;

otherwise

NextTime = time

end

SetEvent(NextTime, 'generate', PacketNr)

Struktur Generator

Generator.type

Generator.par1

Generator.par2

analog für Paketlängen

Es lassen sich auch Ankunftszeiten

aus Dateien (Messungen) einbruingen.




Hauptprogramm und Abarbeitung von Events function t = ProcNextEvent()

…

…

% suche nächstes Event in der Liste –

% das mit dem niedrgsten Zeitwert t

% n ist die Position in der Liste

[t, n] = min([EventList.time]);

% führe die Verarbeitung für das Event

% gemäß Typ aus

switch EventList(n).type

case 'generate'

GeneratePacket(t);

case 'process'

ProcessPacket(t, EventList(n).Pnr);

case 'delete'

DeletePacket(t, EventList(n).Pnr);

otherwise

warning('Fehler Event.type')

end

% Lösche abgearbeitetes Event aus Liste

EventList(n) = [];

Initialisierung • Generator

• Bedieneinheit

• Event- und PacketList leer

MaxSimTime =1000; % maximale Zeit

time = 0; % aktuelle Zeit

pnr = 1; % 1. Paket

SetEvent(time, 'generate', pnr);

while time < MaxSimeTime;

time = ProcNextEvent();

end

Statistische Auswertung • Warte- und Verweilzeiten (Verteilung, MW, …)

• Erzeugte Anforderungen, Blockierungen

• Kundenzahl, Warteschlangenlänge (Vert., MW)

• Zeitlicher Verlauf, ….




Bearbeitung eines ankommenden Paketes

function ProcessPacket(time, PacketNr)

if load < channels, % falls Kanal/Server frei

load = load + 1; % ein Kanal mehr belegt

duration = PacketList(PacketNr).length / rate; % Dauer der Bearbeitung

SetEvent(time + duration, 'delete', PacketNr); % Beende nach Verarbeitung

Else % falls alle Kanäle/Server belegt

q = min(find(queue == 0)); % suche ersten freien Warteplatz

if length(q) == 0, % falls kein Warteplatz vorhanden

SetEvent(time, 'delete', PacketNr);

warning('blocking‘) % Blockierung

else % falls Warteplatz q frei

queue(q) = PacketNr; % trage Paketnummer dort ein

load = load + 1; % erhöhe Last im System

end

end

Processor.load = load; % Aktualisiere Zustand Bediensystem

Processor.queue = queue;




Beendigung eines Prozesses (Paket verarbeitet) function DeletePacket(time, PacketNr)

PacketList(PacketNr).end = time; % Verarbeitung von PacketNr. beendet

PacketList(PacketNr).state = 'deleted'; % Setze Zustand auf 'deleted'

p = queue(1); % Hole nächstes Paket aus Queue (FIFO)

if p == 0, % Falls keines in der Queue

noaction = 1; % keine weiteren Maßnahmen

else % Falls Paket p in der Queue

queue(1) = 0; % Warteplatz freigeben

queue = circshift(queue', -1)'; % die anderen Pakete nachrücken lassen

Packet = PacketList(p); % neu zu verarbeitendes Packet

duration = Packet.length / rate; % Dauer der Bearbeitung

SetEvent(time + duration, 'delete', Packet.nr); % Beende nach Verarbeitung

end

load = max(0, load - 1); % ein Kanal weniger belegt

Processor.load = load; % Aktualisiere Zustand Bediensystem

Processor.queue = queue;




Aufg. B.6: Simulation einfacher Bedienprozesse

a) Simulieren Sie mit Hilfe der vorliegenden SW ein M/M/1-System für zwei verschiedene

sinnvolle Werte für die Last. Stellen Sie MaxSimTime so ein, dass ca. 10000 Pakete erzeugt

werden.

– Wie groß ist die mittlere Zahl der Kunden im Systeme, wie groß die Auslastung pro Kanal?

– Speichern Sie die (komplementären) Verteilungsfunktionen für die Warte- und Verweilzeit in einem

Word-Dokument und notieren Sie dabei die eingestellte Last.

– Was fällt Ihnen beim Vergleich der beiden Diagramme auf? Wie kann man die

Wartewahrscheinlichkeit ablesen?

– Überprüfen Sie, inwieweit das Gesetz von Little erfüllt ist!

b) Simulieren Sie mit Hilfe der vorliegenden SW ein M/M/s-System für eine selbst gewählte

Kanalzahl s > 1. Erhöhen Sie dabei die Last auf das s-fache gegenüber (a). Vergleichen Sie

die Ergebnisse mit denen für (a)!

c) Simulieren Sie ein M/G/1-System für eine Verteilung G Ihrer Wahl (z.B. Pareto, nicht

exponentiell). Dazu müssen beim Paketlängen-Generator den entsprechenden Zufallszahlen-

generator implementieren. Vergleichen Sie die Ergebnisse mit denen für (a)!

d) Simulieren Sie ein M/M/1-System mit einer anderen Bediendisziplin als FIFO (z.B. LIFO).

Dazu sind Änderungen in „DeletePacket“ vorzunehmen. Vergleichen Sie die Ergebnisse mit

denen für (a)!

Modellierung und Simulation 1 - fh-swf.de · • Ergänzen Sie das Programm um die Möglichkeit des...

Documents

Transcript of Modellierung und Simulation 1 - fh-swf.de · • Ergänzen Sie das Programm um die Möglichkeit des...