Yesterday I got an answer from a guy (@Lutz Lehmann) that resolved my problem in using Heun method using matlab. Here is he's answer:
You write the right side of your equation in Matlab format
function dotx = ballbeam(t,x)
dx2 = x(1)*x(4)^2 - g*sin(x(3));
dx4 = (-2*m*x(1)*x(2)*x(4)-m*g*x(1)*cos(x(3))+tau)/(m*x(1)^2+J);
dotx = [ x(2), dx2, x(4), dx4 ];
end
and then feed this into a generic Heun loop
function x = odeheun(f,t,x0)
x = [x0];
for i=1:(length(t)-1)
h = t(i+1)-t(i);
k1 = h*f(t(i),x(i,:));
k2 = h*f(t(i+1),x(i,:)+k1);
x(i+1,:) = x(i,:)+0.5*(k1+k2);
end%for
end%function
and then call it like (tested in octave)
x0 = [x1_ini, 0, 0, 0 ];
t = t0:h:tf;
x = odeheun(@(t,x)ballbeam(t,x),t,x0); % possibly only @ballbeam in matlab
plot(t,x(:,1));
It seems fine and easy but I never used matlab even if professor says to me I have to, so the question is how can I put that code in matlab and see the result? I tried to make a file called genericHeunMethod in which I putted only the Heun loop and after I made another file where I putted my equation and then from this last file in the command window I inserted :
x0 = [x1_ini, 0, 0, 0 ];
t = t0:h:tf;
x = odeheun(@(t,x)ballbeam(t,x),t,x0); % possibly only @ballbeam in matlab
plot(t,x(:,1));
but there are a lot of error in recognising variables etc and I don't know how to do... can someone help me? ty in andvance.