# How to solve systems of ode in matlab?

1 view (last 30 days)
RITIKA Jaiswal on 10 Apr 2023
Edited: Walter Roberson on 11 Apr 2023
I am trying to solve following sets of odes for different cases of control volume .I have written the code here and after running the code I am getting error .Please help me to fix the error.I have considered the value of of all K as constant in the code for now.
Lx = 0.1;
Ly = 0.1;
dx = 0.01;
dy = 0.01;
nx = Lx/dx;
ny = Ly/dy;
Tin = 500;
x = dx/2:dx:Lx+(dx/2);
y = dy/2:dy:Ly+(dy/2);
Nx = nx+1;
Ny = ny+1;
M=zeros(Nx,Ny);
Tini = 600;
W=Tini*(1+M);
Told = W ;
Counter = 0;
for i = 1:Nx
for j = 1:Ny
if (i==1) && (j==Ny)
dvdt=@(t,W)[3*Tin-2*W(i,j)+4*W(i,j-1)];
end
if (i==1) && (j==1)
dvdt=@(t,W)[3*Tin + 4*W(i,j)-3*W(i,j+1)];
end
if (i==Nx) && (j==1)
dvdt=@(t,W)[3*T(i-1,j) + 4*W(i,j)-3*W(i,j+1)];
end
if (i==Nx) && (j==Ny)
dvdt=@(t,W)[3*T(i-1,j-1) + 4*W(i,j)-3*W(i,j+1)];
end
if (j==Ny)
dvdt=@(t,W)[2*W(i-1,j)+4*W(i,j-1)-3*W(i,j)];
end
if (j==1)
dvdt=@(t,W)[4*W(i+1,j+1)-3*W(i,j)];
else
dvdt=@(t,W)[3*T(i-1,j) + 4*W(i+1,j-1)-3*W(i-1,j+1)];
end
[t,W]=ode45(dvdt,[0 0.4],Told)
end
end
plot(t,W)
below are the Following errors which I am getting after running my code:
Attempted to access W(2,2); index out of bounds because size(W)=[121,1].
Error in @(t,W)[4*W(i+1,j+1)-3*W(i,j)]
Error in odearguments (line 87)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 113)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in FVM (line 41)
[t,W]=ode45(dvdt,[0 0.4],Told)

Torsten on 10 Apr 2023
Your code should somehow look like this. Fill in the details.
nx = 11;
ny = 11;
Lx = 0.1;
Ly = 0.1;
dx = Lx/(nx-1);
dy = Ly/(ny-1);
for ix = 1:nx
for iy = 1:ny
T(ix,iy) = 600;
end
end
for ix = 1:nx
for iy = 1:ny
Y0((ix-1)*ny + iy ) = T(ix,iy);
end
end
[time,Y] = ode15s(@(time,Y)fun(time,Y,nx,ny,dx,dy),tspan,Y0)
function dYdt = fun(time,Y,nx,ny,dx,dy)
% Write solution vector Y in matrix
T = zeros(nx,ny);
for ix = 1:nx
for iy = 1:ny
T(ix,iy) = Y((ix-1)*ny + iy);
end
end
% Allocate workspace for the time derivatives in the grid points
dTdt = zeros(nx,ny);
% Set the dTdt expressions of your attached paper (Don't use function
% handles, but just the variables here !)
...
% Write back the dTdt matrix into a column vector to return it to ODE15S
for ix = 1:nx
for iy = 1:ny
dYdt((ix-1)*ny + iy) = dTdt(ix,iy);
end
end
end

### Categories

Find more on Programming in Help Center and File Exchange

R2021a

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!