A quick and dirty solution. Define current
as a global variable (you need to do this in both the base workspace and in the function) and use a for
loop, e.g.
current_vector=1e-7*(6:0.1:8.5);
global current
for k=1:length(current_vector)
current = current_vector(k);
[t{k},y{k}]=ode45(f,<tspan>,<y0>)
end
Replace <tspan>
and <y0>
by the appropriate values/variables. I assume the other constants are defined in the body of the function, in which case your function should look something like this:
function f=RSJ(t,y)
global current
M1 = ... % etc...
f=(current/P)-(M1/P)*sin(y)+(M2/P)*sin(y+E);
end
BTW, I don't see any explicit dependence on time t
in your function...