Attempt at making a Runge Kutta 4 script

Trying to create a simulation a of population using RK4 method, but I'm struggling to get it to run properly. The idea is to get y1 and y2 with 11 numbers, which will represent population amounts.

This is the main script:

close all
clear

alpha = 10^-3;
gamma = 6e-3;
c = 3;
a = 2;

Dt = 0.5;
lambda = sqrt(-6);


y1 = 600;
y2 = 1000;

t = 1:Dt:6;

for i = 1:10
    k11 = Dt*dy1(i,y1,y2);
    k12 = Dt*dy2(i,y1,y2);
    
    k21 = Dt*(dy1(i,y1,y2) + 0.5*k11);
    k22 = Dt*(dy2(i,y1,y2) + 0.5*k12);
    
    k31 = Dt*(dy1(i,y1,y2) + 0.5*k21);
    k32 = Dt*(dy2(i,y1,y2) + 0.5*k22);
    
    k41 = Dt*(dy1(i,y1,y2) + k31);
    k42 = Dt*(dy2(i,y1,y2) + k32);    
    
    y1(i+1) = y1(i) + (1/6)*(k11 + k21 + k31 + k41);
    y2(i+1) = y2(i) + (1/6)*(k12 + k22 + k32 + k42);
end

dy1, dy2 and W are these funtions:

function g = dy1(t,y1,y2)

alpha = 10^-3;
gamma = 6e-3;
c = 3;
a = 2;


lambda = sqrt(-6);

%for t = 1:11
g = (a - alpha.*y2).*y1 - W(t).*y1;
%end
end


function h = dy2(t,y1,y2)

alpha = 10^-3;
gamma = 6e-3;
c = 3;
a = 2;
%y2 = a/alpha;
%y1 = c/gamma;
Dt = 0.5;
lambda = sqrt(-6);
y1 = 600;
y2 = 1000;

h = (-c + gamma*y1)*y2;
end


function f = W(t)
    f = [];
    if 0 <= t && t < (3/12)
        f = 0;
    end
    if (3/12) <= t && t <= (8/12)
        f = 1;
    end
    if (8/12) < t && t < 1
        f = 0;
    end
    if t >= 1
        t = t - 1;
    end
end

No watter how i try to change the script or the functions, I keep getting differnt kinds of errors. I feel like the general idea is there, but I just don't understand what exactly I'm doing that matlab doesn't like. Any help or insight would be greatly appreciated.

[Jan]

The readers cannot guess, what you want to achieve. Why do you restrict the values of t at the end of the W(t) function? t is not used afterwards anymore, so you can omit the change to t-1.

You define t = 1:Dt:6, but do not use these values later on. I guess you want to replace the loop for i = 1:10 by

tVec = 1:Dt:6,
for k = 1:numel(tVec)
    t = tVec(k);
    k11 = Dt*dy1(t,y1,y2);
    %            ^ t instead of i
end

"I keep getting differnt kinds of errors" - then post the code and a copy of the complete error message. It is easier to fix an error than to guess, what the error is.

[Tom Oud]

The error message I'm getting right now is the following:

In an assignment A(:) = B, the number of elements in A and B must be the same.

Error in Sub2 (line 33)

y1(t+1) = y1(t) + (1/6)*(k11 + k21 + k31 + k41);

(I've gotten a couple different ones, which is why I thought posting it might just be confusing)

Thanks for the answer though!