You should build a simple case study outside of this function. It seems that you want to build an array whose rows contain the evolution of each component of y0. You should try to see (e.g. in the command window) whether you are able to initialize such an array, and then whether you are able to fill it in correctly.
To illustrate the process:
>> a = 2 ; b = 12 ; n = 10 ; h = (b-a)/n ;
>> y0 = [3, 4, 5] ;
>> x = [a zeros(1,n)] ;
>> y = [y0 zeros(length(y0),n)];
Error using horzcat
Dimensions of matrices being concatenated are not consistent.
well, it seems that when the user defines y0 as a row vector, it doesn't work. Therefore, we can use linear indexing to guarantee that we are dealing with a column vector:
>> y = [y0(:), zeros(length(y0),n)];
>> size(y)
ans =
3 11
>> y
y =
3 0 0 0 0 0 0 0 0 0 0
4 0 0 0 0 0 0 0 0 0 0
5 0 0 0 0 0 0 0 0 0 0
which seems to be working. Next, we define a simple f for the test:
>> f = @(x,y) 2*x + 0*y ;
Then we try to perform the first iteration..
>> i = 1 ;
>> x(i+1) = x(i)+h;
>> x
x =
2 3 0 0 0 0 0 0 0 0 0
this seems to be working.. let's see how it goes for y:
>> y(i+1) = y(i)+h*f(x(i),y(i));
>> y
y =
3 0 0 0 0 0 0 0 0 0 0
7 0 0 0 0 0 0 0 0 0 0
5 0 0 0 0 0 0 0 0 0 0
this doesn't work; the 4 from the initial condition was replaced with 7. But if we look a little closer, it seems that you are indexing y, which is a 2D array, as if it were a 1D array. Looking at the array, it seems that step i should define y(:,i+1) based on x(i) and y(:,i) .. and here I let you bring the update. Note that the upper boundary of the range defining the loop statement is n in your code, so the last loop will define the n+1-th column of x and y, and you might want hence to review this upper boundary.
I can't write the solution for you as I suspect that it is a homework, but this shows you how to proceed: simple steps where you display/test all the steps of the computation.
