[Matlab] Ziektekiemen

Pagina: 1
Acties:
  • 713 views

Onderwerpen


Acties:
  • 0 Henk 'm!

  • kestemodon
  • Registratie: November 2008
  • Laatst online: 10:21
Ik moet een functie schrijven voor op de universiteit.
Mijn scripts en functies zijn klaar, leveren geen foutmeldingen meer op, maar toch klopt het van geen kanten.
Al mijn waardes worden herleid naar nul :(

We moeten aan de hand van ons nummer van de unif de fictieve beweging volgen van bacteriën.

De bacterie verspreiden zich als volgt:
Als het aantal bacteriën boven de targetgrens liggen gebeurt er niets, als het kleiner of gelijk is dan gaat de volgende stappen op.
Indien het een hoek is van de matrix, dan wordt het aantal bacteriën door 2 gedeeld en verspreid over de 2 aanliggende velden, als het aanwezige getal niet te delen is door 2 dan blijft de rest van de deling achter.
De methode voor de randen en het midden van de matrix is hetzelfde buiten dat er dan 3 of 4 buren zijn.

Ik vermoed dat mijn fout in de functie ziektekiemen zit.
Maar ik post toch maar even alles voor de zekerheid :P

Alvast bedankt

Het inleidend script
code:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
clear all

Mycode = [ 0 2 5 5 7 3 1];
Myvec=Mycode(4:7);
Mylimit=max(Mycode)*sum(Mycode);

N=zeros(21);
N(1,11)=Myvec(1);
N(21,11)=Myvec(2);
N(11,1)=Myvec(3);
N(11,21)=Myvec(4);

target=Mylimit;

% Bepaal gemiddeld aantal ziektekiemen en maxima
 
   [n,m] = size(N);
   Maxima = zeros(1,31);
   Gemiddelde = zeros(1,31);
   
  t = 0:30;
   % Bereken situatie na 30 dagen
   
   for p = 1:30
     M = ziektekiemen ( N, target);
     Maxima(1,p) = max(max(N(1:n)));
     Gemiddelde(1,p) = (sum(sum(N(1:n,:)))) /(n*m);
     M = N + M;
      
   end
   Gemiddelde(1,31) = (sum(sum(N(1:n,:)))) /(n*m);
   Maxima(1,31)     = max(max(N(1:n)));
   save matrix.dat M -ascii
   
   % Plot 'Gemiddelde' en 'Maxima' in functie van 'Tijd'
   clf
   plot(t, Gemiddelde, 'or-');
   hold on;
   plot(t, Maxima, '*b-');
   grid on;
   xlabel('Tijd (in dagen)');
   ylabel('Aantal ziektekiemen');
   legend('Gemiddelde','Maximum');
   title('Donaat  & 1');


Functie om het aantal buren te bepalen
code:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
function a=aantal_buren(n,m,r,k)

if r>n ||k>m
    disp 'Kan aantal buren niet berekenen, dimensies kloppen niet'
    a=0;

else
    if r==1 %voor de eerste rij geldt :
        if k>1 && k<m % als het het eerste element niet is, en het is ook niet het laatste dan zijn er 3 buren
         a=3;
        else 
            a=2;  %het is ofwel het eerste element of het laatste van de eerste rij
        end
    
    elseif r==n %voor de laatste rij geldt:
        if k>1&&k<m %als het het eerste element niet is, en het is ook niet het laatste dan zijn er 3 buren
            a=3;
        else %het is het eerste of het laatste element
            a=2;
        end
    
    
    else % voor alle rijen verschillend van de eerste of de laatste
      if m==k || n==r % het zit in de eerste of de laatste kolom
          a=3;
      else %het zit in het midden van de matrix
          a=4;
      end
    end
end


Functie die bepaald wat er met de bacteriën gebeurd
code:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
function [weg,blijf]=reproductie(z,a,target)

if z<=target         % Aantal kiemen kleiner dan de max dus
    z=2*z;             % Verdubbelen kiemen
    t=mod(z,a);     % Onderzoeken of z deelbaar is door a
    if t==0             % Het is deelbaar dus:
        weg=z/a;     % Het aantal kiemen dat weggaat is gelijk aan het # kiemen / # buurvelden
        blijf=0;         % Er blijven geen kiemen achter
    else                 % Z is niet deelbaar door a
        z=z-t;           % Er wordt een z berekend die wel deelbaar is door a
        blijf=t;          % Dat getal t is gelijk aan het aantal blijvers
        weg=z/a;     % Het deelbare getal / #buurvelden geet het aantal vertrekkers weer
    end
else                % Er zijn al teveel kiemen. Alle kiemen blijven in dit veld
    blijf=z;
    weg=0;
end


Functie die matrix vernieuwt (Hier zit de fout waarschijnlijk)
code:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
function M=ziektekiemen(N,target)

[n,m]=size(N);
M=zeros(n,m);
B=zeros(n,m);

for r=1:n
    for k=1:m
        if M(r,k) >0
            z=M(r,k);
            a=aantal_buren(n,m,r,k);
            [weg,blijf]=reproductie(z,a,target);
            B(r,k)=blijf;
        %end
        if a==2                             % Het is een hoek
            if r==1  && k==1                %linkerbovenhoek
                M(r+1,k)=weg+M(r+1,k);
                M(r,k+1)=weg+M(r,k+1);
                
            elseif r==1 && k==m             %rechterbovenhoek
                M(r+1,k)=weg+M(r+1,k);
                M(r,k-1)=weg+M(r,k-1);
                
            elseif r==n && k==1             %linkerbenedenhoek
                M(r-1,k)=weg+M(r-1,k);
                M(r,k+1)=weg+M(r,k+1);
                
            else                            %rechterbenedenhoek
                M(r-1,k)=weg+M(r-1,k);
                M(r,k-1)=weg+M(r,k-1);
            end
            
        elseif a==3                         % Het is een rand, maar geen hoek
            if r==1                         % linkerrand
                M(1,k-1)=weg+M(1,k-1);
                M(1,k+11)=weg+M(1,k+1);
                M(2,k)=weg+M(2,k);
                
            elseif r==n                     % rechterrand
                M(n,k-1)=weg+M(n,k-1);
                M(n,k+1)=weg+M(n,k+1);
                M(n-1,k)=weg+M(n-1,k);
                
            elseif k==1                     %bovenrand
                M(r-1,1)=weg+M(r-1,1);
                M(r+1,1)=weg+M(r+1,1);
                M(r,2)=weg+M(r,2);
                
            else                            % onderrand
                M(r-1,m)=weg+M(r-1,m);
                M(r+1,m)=weg+M(r+1,m);
                M(r,m-1)=weg+M(r,m-1);
            end
            
        else                           % binnenkant matrix
            M(r,k-1)=weg+M(r,k-1);
            M(r,k+1)=weg+M(r,k11);
            M(r-1,k)=weg+M(r-1,k);
            M(r+1,k)=weg+M(r+1,k);
        end
        end
    end
end

M=M+B;                                  %achterblijvers en verplaatsingen samen tellen

Acties:
  • 0 Henk 'm!

  • goeievraag
  • Registratie: Mei 2005
  • Niet online
Op regel 9-10 in de ziektekiemenfunctie moet volgens mij staan
code:
1
2
 if N(r,k) >0
            z=N(r,k);

want M heb je zelf op nul gezet. Verder moet op regels 17-59 allen M-en rechts van is tekens vervangen worden door N-en.

Verder in het aanroepscript moet regel 25 worden vervangen door:
code:
1
 M = ziektekiemen ( M, target);

en moet je ergens helemaal aan begin nog zeggen M=N. (En je hebt nog geen tijdlus ;))

I fell with my reet in the prikkeldreed. How do I frommel my hand in the cookiestrommel? The sun is shining over het randje van de heining. I'm shocking klem between de deuren van de tram.


Acties:
  • 0 Henk 'm!

  • ILUsion
  • Registratie: Augustus 2003
  • Laatst online: 23-01 08:12
Ik zou in de eerste plaats je code eens opkuisen, het leest langs geen kant: in je hoofdscript zijn de belangrijkste variabelen voor zover ik kan zien M en N terwijl hun naam niet duidelijk zegt wat ze zijn. MATLAB kan namen geven tot 63 tekens, gebruik die mogelijkheden dan ook. En ja, het typt wat trager, daarom is er ook een autocompletion (via de TAB-toets) ingebouwd; vroeger programmeerde ik ook met zo kort mogelijke variabelen, maar dat helpt je niet om je code te begrijpen.

Dat je itereert over p of t of het maakt eigenlijk niet veel uit, is niet zo'n probleem zolang het uit de context duidelijk is wat de betekenis van die index is. Als ik helemaal ga mierenneuken, doe ook iets aan die indentatie (gewoon selecteren en TAB of Shift-TAB). Meer aanbevelingen kan je trouwens vinden op http://www.mathworks.com/matlabcentral/fileexchange/2529.

Maar goed, om in te gaan op jouw code, ik denk dat het sowieso al in je hoofdprogramma misloopt, veel verder ben ik nog niet gaan kijken:
code:
1
2
3
4
5
6
for p = 1:30
     M = ziektekiemen ( N, target);
     Maxima(1,p) = max(max(N(1:n)));
     Gemiddelde(1,p) = (sum(sum(N(1:n,:)))) /(n*m);
     M = N + M;
end


Qua functionaliteit loopt daar iets mis, wat daar staat is in pseudocode
code:
1
2
3
4
for i = 1:30
    M = functie(N);
    [maxima, gemiddelde] = berekenStatistieken(N);
    M = andereFunctie(M,N);

Je gebruikt op geen enkele manier dus een terugkoppeling in die code: N is invariant in die lus, dus zullen je statistieken van N ook invariant zijn. Let er ook op dat je vanboven in je code niet alle bronnen invoegt in je matrix.

Verder een klein trucje die je wat kan afhelpen van dubbele maxima/sommen, zoals de code er nu staat, zou ik hem schrijven als:
code:
1
2
Maxima(p) = max(N(:));
Gemiddelde(p) = mean(N(:));


Misschien handig om te weten, maar ik zou tijdens het debuggen in je lus ook zorgen dat je visualiseert wat er in je matrix gebeurt. Dat kan ofwel met de 3D-plots (surf, contour, ...) ofwel met de functie spy (die geeft aan welke elementen verschillend zijn van 0). Om alles automatisch te volgen op je hoog niveau kan je ook gebruik maken van pause(0.5) in je lus, daarmee pauzeer je een halve seconde elke iteratie.

Acties:
  • 0 Henk 'm!

  • kestemodon
  • Registratie: November 2008
  • Laatst online: 10:21
Sorry dat het niet goed leest, maar de variabelen zijn door de unif gekozen en moeten zo gebruikt worden...

Allebei al vreselijk bedankt, ik ga het deze avond eens bekijken :P

Acties:
  • 0 Henk 'm!

  • NMe
  • Registratie: Februari 2004
  • Laatst online: 09-09 13:58

NMe

Quia Ego Sic Dico.

En daarmee kun je het doen helaas. Je dumpt je vraagstelling, je dumpt je code. En vervolgens verwacht je dat wij je code doorspitten en je probleem voor je oplossen? Sorry, maar daar beginnen we niet aan: debuggen mag je zelf doen. Als je eenmaal gedebugd hebt kun je je vraag terugbrengen naar een handjevol regels, en dán kun je eventueel een gerichte vraag stellen op het forum. Met "dit is mijn code, fix het maar" ga je hier buiten de antwoorden die je al gehad hebt geen antwoorden meer krijgen.

Succes in elk geval.

'E's fighting in there!' he stuttered, grabbing the captain's arm.
'All by himself?' said the captain.
'No, with everyone!' shouted Nobby, hopping from one foot to the other.


Dit topic is gesloten.