# Coupled differential equation to solve fluids

3 views (last 30 days)
Dhruv Sirohi on 29 Oct 2019
Commented: Dhruv Sirohi on 30 Oct 2019
x0_t=@(t) (-3+0.2*sin(0.2*t));
x1_t=@(t) (3+0.1*sin(0.1*t));
x2_t=@t (0.0+0.13*sin(0.13*t));
y0_t= @(t) (1.5+0.12*sin(0.12*t));
y1_t= @(t) (1.5+0.15*sin(0.14*t));
y2_t= @(t) (-3.5+0.05*sin(0.14*t))
r0=@(x,y,t) (((x-x0_t)^2 + (y-y0_t)^2)+0.2);
r1= @(x,y,t) (((x-x1_t)^2+(y-y1_t)^2)+0.2);
r2= @(x,y,t) (((x-x2_t)^2+(y-y2_t)^2)+0.2);
ode1=diff(x)==-5*((y-y0_t(t))/(r0(x,y,t))^2)+5*((y-y1_t(t))/(r1(x,y,t))^2)+5*((x-x2_t(t))/(r2(x,y,t))^2)
ode2=diff(y)==5*((x-x0_t(t))/(r0(x,y,t))^2)-5*((x-x1_t(t))/(r1(x,y,t))^2)+5*((y-y2_t(t))/(r2(x,y,t))^2)
odes=[ode1;ode2]
>> cond1=x(0)==0;
>> cond2=y(0)==0;
>> conds=[cond1;cond2];
>> [xsol(t),ysol(t)]=dsolve(odes,conds)
Warning: Unable to find explicit solution.
> In dsolve (line 190)
Error using sym/subsindex (line 855)
Invalid indexing or function definition. Indexing must follow MATLAB indexing. Function arguments must be
symbolic variables, and function body must be sym expression.
I am a little out of touch with matlab, and need to solve a fluid mechanics problem with it, particularly finding the pathline. Can someone help me with this warning+error. Should I be using ode45 instead of dsolve?

Stephan on 29 Oct 2019
Edited: Stephan on 29 Oct 2019
syms x(t) y(t) x0_t x1_t x2_t y0_t y1_t y2_t r0 r1 r2
f1 = x0_t == (-3+0.2*sin(0.2*t));
f2 = x1_t == (3+0.1*sin(0.1*t));
f3 = x2_t == (0.0+0.13*sin(0.13*t));
f4 = y0_t == (1.5+0.12*sin(0.12*t));
f5 = y1_t == (1.5+0.15*sin(0.14*t));
f6 = y2_t == (-3.5+0.05*sin(0.14*t));
f7 = r0 == (((x-x0_t)^2 + (y-y0_t)^2)+0.2);
f8 = r1 == (((x-x1_t)^2+(y-y1_t)^2)+0.2);
f9 = r2 == (((x-x2_t)^2+(y-y2_t)^2)+0.2);
ode2=diff(y,t)==5*((x-x0_t)/(r0)^2)-5*((x-x1_t)/(r1)^2)+5*((y-y2_t)/(r2)^2);
ode(1) = subs(ode2,[r0, r1, r2],[rhs(f7), rhs(f8), rhs(f9)]);
ode(1) = subs(ode(1),[x0_t, x1_t, x2_t, y0_t, y1_t, y2_t],[rhs(f1),...
rhs(f2), rhs(f3), rhs(f4), rhs(f5), rhs(f6)]);
ode1=diff(x,t)==-5*((y-y0_t)/(r0)^2)+5*((y-y1_t)/(r1)^2)+5*((x-x2_t)/(r2)^2);
ode(2) = subs(ode1,[r0, r1, r2],[rhs(f7), rhs(f8), rhs(f9)]);
ode(2) = subs(ode(2),[x0_t, x1_t, x2_t, y0_t, y1_t, y2_t],[rhs(f1),...
rhs(f2), rhs(f3), rhs(f4), rhs(f5), rhs(f6)]);
[S,V] = odeToVectorField(ode);
fun = matlabFunction(S,'vars',{'t','Y'});
tspan = [0 50];
conds = [0 0];
[t,res] = ode45(fun, tspan, conds);
figure(1)
plot(t,res,'LineWidth',2)
figure(2)
plot3(t,res(:,1),res(:,2),'LineWidth',2)
grid on
xlabel('Time');
ylabel('x-values');
zlabel('y-values');

#### 1 Comment

Dhruv Sirohi on 30 Oct 2019
Thank! I only needed the 2D plot, and now I have it with the time axis!!