FLÄCHENINHALT

%Plotte zunächst die Figur - sehr aufwendig.
hold on;
plot(1:0.05:3,1,'b.','MarkerSize',5);
plot(8:0.05:10,1,'b.','MarkerSize',5);
plot(3:0.05:6,2,'b.','MarkerSize',5);
plot(6:0.05:8,3,'b.','MarkerSize',5);
plot(5:0.05:6,4,'b.','MarkerSize',5);
plot(7:0.05:10,4,'b.','MarkerSize',5);
plot(1:0.05:2,5,'b.','MarkerSize',5);
plot(3:0.05:5,6,'b.','MarkerSize',5);
plot(7:0.05:9,6,'b.','MarkerSize',5);
plot(2:0.05:3,8,'b.','MarkerSize',5);
plot(6:0.05:9,9,'b.','MarkerSize',5);
plot(1,1:0.05:5,'b.','MarkerSize',5);
plot(2,5:0.05:8,'b.','MarkerSize',5);
plot(3,1:0.05:2,'b.','MarkerSize',5);
plot(3,6:0.05:8,'b.','MarkerSize',5);
plot(5,4:0.05:6,'b.','MarkerSize',5);
plot(6,2:0.05:3,'b.','MarkerSize',5);
plot(6,4:0.05:9,'b.','MarkerSize',5);
plot(7,4:0.05:6,'b.','MarkerSize',5);
plot(8,1:0.05:3,'b.','MarkerSize',5);
plot(9,6:0.05:9,'b.','MarkerSize',5);
plot(10,1:0.05:4,'b.','MarkerSize',5);
axis([0 10 0 10]);

n = 50;

for l=1:length(n)
    x = 10*rand(2,n);
    plot(x(1,:), x(2,:),'r.','MarkerSize',5);
    treffer = ((1<= x(1,:) & x(1,:) <= 3) & (1 <= x(2,:) & x(2,:) <= 5)) | ((2<= x(1,:) & x(1,:) <= 5) & (2 <= x(2,:) & x(2,:) <= 6)) | ((2 <= x(1,:) & x(1,:) <= 3) & (6 <= x(2,:) & x(2,:) <= 8)) | ((5 <= x(1,:) & x(1,:) <= 6) & (2 <= x(2,:) & x(2,:) <= 4)) | ((6 <= x(1,:) & x(1,:) <= 7) & (3 <= x(2,:) & x(2,:) <= 9)) | ((7 <= x(1,:) & x(1,:) <= 10) & (3 <= x(2,:) & x(2,:) <= 4)) | ((8 <= x(1,:) & x(1,:) <= 10) & (1 <= x(2,:) & x(2,:) <= 3));
    f_anteil = sum(treffer)/n;
    f_anteil
end;

---------------------

GEBURTS- UND TODESPROZESSE
%Hier wird der in Imke Rothers Ausarbeitung zu Geburts- und Todesprozessen vorgestellter Algorithmus programmiert.

%Wir simulieren: Zum Zeitpunkt 0 sei die Population bei i.
zeit_komplett = 0;
i = 5;

%Maximale Dauer der Simulation - Sonst kann es zu Endlosschleifen kommen.
dds = 100;

fprintf('\n Neue Zeit ti     Derzeitige Zeit t   Populationsgröße X(t) \n');
fprintf('-----------------------------------------------------\n');
fprintf('  0            %1.2f        %10d \n', zeit_komplett, i);

%Für den Verlauf der Population
hold on;
plot(0,i,'b.','MarkerSize',5);
max_i = i;


while ((i > 0) & (dds > 0))
    %Geburts- und Todesprozessparameter werden festgelegt in m bzw. n
    m = 1/(i+1);
    n = 1/(i+1);

    %Gleichverteilte Zufallszahl generieren, um an eine exponentialverteilte Zeit t zu kommen.
    u = rand(1,1);
    zeit_neu = -log(u)/(m+n);

    plot(zeit_komplett:0.05:zeit_komplett+zeit_neu,i,'b.','MarkerSize',5);

    %Der in der Ausarbeitung beschrieben 'Münzwurf' wird nun simuliert.
    u = rand(1,1);
    if (u <= m/(m+n))
        i = i+1; %Die Population vergrößert sich.
    else
        i = i-1; %Die Population verringert sich.
    end;

    if(max_i < i)
        max_i = i;
    end;

    %Zeichnen des Verlaufs

    zeit_komplett = zeit_komplett + zeit_neu;

    %Verringern der max. Schleifendurchläufe
    dds = dds - 1;

    fprintf('  %1.2f         %1.2f        %10d \n',zeit_neu, zeit_komplett, i);
end;


plot(zeit_komplett:0.05:zeit_komplett+5,i,'b.','MarkerSize',5);


axis([0 zeit_komplett+5 0 max_i+1]);
xlabel('Zeit t','Fontsize',14);
ylabel('Populationsgröße','Fontsize',14);
title('Verlauf eines Geburts- und Todesprozesses','Fontsize',18);

----------------

HEIMWEGPROBLEM MIT GRAPHISCHER AUSGABE
%erst Stadt plotten, dann Verlauf des Mathematikers

mathe_start_x = 0;
mathe_start_y = 0;
wx = 5;
wy = 0;

hold on;
for y = -5:5
    x = -5:0.01:5;
    plot(x,y,'k-');
end;
for x = -5:5
    y = -5:0.01:5;
    plot(x,y,'k-');
end;
plot(mathe_start_x,mathe_start_y,'r.','MarkerSize',5);
plot(wx,wy,'g.','MarkerSize',5);
axis([-6 6 -6 6]);

N  = 1;
M  = 0;
for l=1:N
    mx = mathe_start_x;
    my = mathe_start_y;
    while ((mx ~= 5) & (mx ~= -5) & (my ~= 5) & (my ~= -5))
        ziehung = 4*rand(1,1);
        if ziehung <= 1
            plot(mx,my:0.05:my + 1,'r.','MarkerSize',5);
            my = my + 1;
        elseif ziehung <= 2
            plot(mx:0.05:mx + 1,my,'r.','MarkerSize',5);
            mx = mx + 1;           
        elseif ziehung <= 3
            plot(mx,my:-0.05:my - 1,'r.','MarkerSize',5);
            my = my - 1;
        else
            plot(mx:-0.05:mx - 1,my,'r.','MarkerSize',5);
            mx = mx - 1;
        end;
    end;
    if (mx == wx) & (my == wy)
        M = M + 1;
    end;
plot(mx,my,'r.','MarkerSize',5);
end;

--------------------------------

HEIMWEGPROBLEM
% Generiere Stadt, in der sich der Mathematiker aufhält
%figure(1);
%plot();
%axis([-6 -6 -6 -6]);


%maximal 10^6 ! Sonst dauert es zu lange.
%auch mal verschiedene Wohnungsstandorte probieren und etwas Verblüffendes feststellen.

N  = 100000;
M  = 0;
mathe_start_x = 0;
mathe_start_y = 0;
wx = 5;
wy = 0;
for l=1:N
    mx = mathe_start_x;
    my = mathe_start_y;
    while ((mx ~= 5) & (mx ~= -5) & (my ~= 5) & (my ~= -5))
        ziehung = 4*rand(1,1);
        if ziehung <= 1
            my = my + 1;
        elseif ziehung <= 2
            mx = mx + 1;
        elseif ziehung <= 3
            my = my - 1;
        else
            mx = mx - 1;
        end;
    end;
    if (mx == wx) & (my == wy)
        M = M + 1;
    end;
end;
ausgabe = ['P(Betrunkener findet selbst nach Hause) = ',num2str(M/N)];
ausgabe

-----------------

PRODUKTIONSPROZESS-VERTEILUNG

k = [0 1 2 3 4 5 6];
l = 2;            %Lamdba
summe = 0;
n = 24;            %Arbeitsstunden
a2 = 0;
prod = 2.5;        %Anzahl fertiger Erzeugnise pro Stunde

fprintf('\n  k  P(K=k)  F(k) \n');
fprintf('-------------------------\n');

for i = 1:length(k)
    wkeit = (l^k(i)/gamma(k(i)+1)) * exp(-l);
    summe = summe + wkeit;
    verteilung(i) = summe;
    fprintf('%3d  %1.2f    %1.2f \n',k(i),wkeit,summe);
end;   


fprintf('\n I   U      A1  A2 \n');
fprintf('-------------------------\n');

for m = 1:n
    x = rand(1,1);
    if(x <= verteilung(1))
        a1 = 0;
    elseif x <= verteilung(2)
        a1 = 1;
    elseif x <= verteilung(3)
        a1 = 2;
    elseif x <= verteilung(4)
        a1 = 3;
    elseif x <= verteilung(5)
        a1 = 4;
    elseif x <= verteilung(6)
        a1 = 5;
    else
        a1 = 6;
    end;
    if(a2 + a1 - prod <= 0)
        a2 = 0;
    else
        a2 = a2 + a1 - prod;
    end;
    fprintf('%2d   %1.2f  %2d   %1.1f \n',m,x,a1,a2);
end;

------------------------

WIENER-PROZESS

%f(x) = sqrt(2/pi)*exp(-(x^2)/2)
%Hilfsfunktion hierfür ist h(x) = exp(-x), also Parameter lambda = 1!
%Generiere ein Y mit Verteilung h
%Generiere ein u als gleichverteilte Zufallsvariable aus ]0,1[
%Wenn u <= f(Y)/(c*g(Y)), setze Zufallszahl = Y
%Ein c hierfür wäre c = sqrt(2*exp(1)/pi)
%Entscheide dann mit 50% W'keit ob das Vorzeichen + oder - ist.

%Konstante C
c = sqrt(2*exp(1)/pi);

hold on;
weg = 0;
max_weg = 0;
min_weg = 0;

for i=0:0.1:100

        %Zunächst das Y
        z = rand(1,1);
        y = -log(z)/1;

        %jetzt prüfen, ob das y 'passend' ist
        u = rand(1,1);
        f_y = sqrt(2/pi)*exp(-(y^2)/2);
        h_y = exp(-y);

        %Wenn es der Überprüfung nicht standhält, wird ein neues generiert.
        while (u > f_y/(c*h_y))
                z = rand(1,1);
                y = -log(z)/1;

                u = rand(1,1);
                f_y = sqrt(2/pi)*exp(-(y^2)/2);
                h_y = exp(-y);
        end;


        %Vorzeichenwahl der endgültigen Zufallszahl
        v = rand(1,1);
        if (v <= 0.5)
                x = -y;
                plot(i,weg:-0.1:weg+x,'b.','MarkerSize',5);
        else
                x = y;
                plot(i,weg:0.1:weg+x,'b.','MarkerSize',5);
        end;

        weg = weg + x;
        if(max_weg < weg)
                max_weg = weg;
        end;
        if(min_weg > weg)
                min_weg = weg;
        end;
end;

axis([-1 101 min_weg-1 max_weg+1]);

----------------------------

PI APPROXIMIEREN

hold off;
n = [10 100 1000 10000];
for l = 1: length(n);
    figure(n(l));
    z = 0:0.01:1;
    y = sqrt(1-z.^2);
    x=rand(2,n(l));
    plot(x(1,:),x(2,:),'.',z,y, 'r');
    axis([0 1.2 0 1.2]);
    treffer = sqrt(x(1,:).^2 + x(2,:).^2) <= 1;
    pi_naeherung = 4*sum(treffer)/n(l);
    ausgabe = ['\pi = ',num2str(pi_naeherung)];
    text(0.4,-0.1,ausgabe,'Fontsize',14);
    ausgabe2 = ['Für n = ',num2str(n(l)),' ergibt sich Pi = ',num2str(pi_naeherung)]
end;