I have to solve numerically this set of Ordinary Differential Equations
$$ \frac{dx_1}{ds} = \frac{1}{x_1} \left[x_2 \left(a + \frac{x_2}{s}\right)-\alpha x_1 z\right]$$
$$ \frac{dx_2}{ds} = -\frac{1}{x_1} \left[x_2 \left(b + x_2\right)z\right] - \left(c+\frac{x_2}{s} \right)$$
where $z(s) = \sqrt{x_1^2 +\left(b + x_2 \right)^2}$. The equations are time-independent, $s$ is a distance, $a(s)$, $b(s)$ and $c(s)$ are constant vectors, i.e. $f(s)$; and $\alpha$ is a constant. I've tried to solve using Matlab's ODE45 to no avail because perhaps the fact that the constants are vectors.
Here's the code and some vector examples:
global a b c d
y = ode45(@odeeqns, s, [1 2],[0 0]) % -> I'm not sure about these because I'm working with distance, not time.
function dyds = odeeqns(s,y)
global a b c d
dyds(:,1) = 1./y(1) .*(y(2) .* (a + y(2)./ s) - alpha .* y(1) .* sqrt(y(1).^2 + (b + y(2)).^2));
dyds(:,2) = - (y(2) .* (b + y(2)) .* sqrt(y(1).^2 + (b + y(2)).^2))./y(1) - y(2)./s - c;
plot(dyds);hold on
return
alpha = .0167;
a = [0.0000 0.2985 1.4973 2.4266 2.7838 2.7917 2.6397 2.4295 2.2094 2.0000];
b = [0.0000 0.0298 0.2246 0.4853 0.6960 0.8375 0.9239 0.9718 0.9942 1.0000];
c = [0.0299 0.2615 0.9764 1.4490 1.5680 1.5098 1.3870 1.2499 1.1188 1.0058];
s = [1 2 3 4 5 6 7 8 9 10]; % km
This is the error message:
??? Error using ==> odearguments at 113
YPRIME must return a column vector.
Error in ==> ode45 at 173
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in ==> windfield_powell at 33
[T,sols] = ode45(@yprime, [1 2],[1 1]);
Any help please? Thanks!
You can try something like that, first define a function
function dx=f(s,x)
a=1.2; b=2.3; c=3.4; alpha=4.5;
dx=zeros(2,1);
z=sqrt(x(1)^2+(b+x(2))^2);
dx(1)=1/x(1)(x(2)(a+x(2)/s)-alpha*x(1)*z);
dx(2)=-1/x(1)(x(2)(b+x(2))*z)-(c+x(2)/s);
and then use ode45 with vectors:
[T,Y] = ode45(@f,[3 12],[1 1])
plot(T,Y(:,1),'-',T,Y(:,2),'-.')
I hope it helps although perhaps you don't want to define the parameters $a$, $b$, $c$, $\alpha$ inside the function $f$