Solving system of equations

2 次查看（过去 30 天）

EldaEbrithil 2020-5-27

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/535078-solving-system-of-equations

评论： darova 2020-5-30

Hi all

i have a question about solving this system of equations. Tt, Pt and M are related to space and time due to i and j; i want to solve the system maintaining that dependence, so the result will be a matrix respectively for Tt, Pt and M. When i try to solve, i obtain "Out of range subscript." error. gamma, deltax and deltat are constant

Thanks to all

Tt=zeros(length(x),length(t));
Pt=zeros(length(x),length(t));
M=zeros(length(x),length(t));
Tt(1,1)=3.000555630247608e+02;
Pt(1,1)=2.201018491400215e+05;
M(1,1)=0.023565919700319;
for j=1:length(t)-1
for i=2:length(x)-1
Alla = cell(length(x),length(t));
Allb = cell(length(x),length(t));
Allc = cell(length(x),length(t));
syms Tt Pt M 
[sola,solb,solc]=vpasolve(Tt(i,j+1)==0.5*(Tt(i+1,j)-Tt(i-1,j))+((1+((gamma-1)/2)*M(i,j)^2)^(gamma/(gamma-1)))*((Tt(i+1,j)-Tt(i-1,j))*deltat/(2*deltax))+((1+((gamma-1)/2)*M(i,j)^2))*((Pt(i+1,j)-Pt(i-1,j))*deltat/(2*deltax)),...
Pt(i,j+1)==0.5*(Pt(i+1,j)-Pt(i-1,j))+2*((1+((gamma-1)/2)*M(i,j)^2)^(gamma/(gamma-1)))*((Tt(i+1,j)-Tt(i-1,j))*deltat/(2*deltax))+3*((1+((gamma-1)/2)*M(i,j)^2))*((Pt(i+1,j)-Pt(i-1,j))*deltat/(2*deltax)),...
M(i,j+1)==0.5*(M(i+1,j)-M(i-1,j))+2*((1+((gamma-1)/2)*M(i,j)^2)^(gamma/(gamma-1)))*((Tt(i+1,j)-Tt(i-1,j))*deltat/(2*deltax))+3*((1+((gamma-1)/2)*M(i,j)^2))*((Pt(i+1,j)-Pt(i-1,j))*deltat/(2*deltax)));
Alla{i,j} = sola;
Allb{i,j} = solb;
Allc{i,j} = solc;
end
end

17 个评论
显示 15更早的评论隐藏 15更早的评论

darova 2020-5-27

I can't explain it here

can be re-written as (P(i,j+1)-P(i,j))/dt

can be re-written as (P(i+1,j)-P(i,j))/dx

you what i mean?

Read about this method. Read about "Method of lines"

EldaEbrithil 2020-5-28

Laxwendroff.PNG

Yes i understand, but i think it is what similar to what i have done in my code, the only difference is related to the typology of discretization: you have used a forward discretiation in space and time, i have used a Forward Time Centered Space, FTCS discretization. Thi is the only difference, but the problem i have is easier than you think: i do not understand how to write the code for solving the system of equations practically.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

darova 2020-5-28

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/535078-solving-system-of-equations#answer_441251

Here is a simple example. I hope it's clear enough. TR, TL, TD - boundary conditions (right, left and down boundaries)

2 个评论
显示无隐藏无

EldaEbrithil 2020-5-30

I have tried to implement the method for the equation tht you give me in the.m file but i am not very confident about the results

clc,clear
% problem definition and discretization
dx = 0.01;
dt = 0.008;
xdomain = [0 1];
tdomain = [0 1];
nx = round((xdomain(2)-xdomain(1))/dx);
nt = round((tdomain(2)-tdomain(1))/dt);
x = linspace(xdomain(1),xdomain(2),nx);
t = linspace(tdomain(1),tdomain(2),nt);
u = zeros(nt,nx);
% du/dt - 2*t*du/dx = 0
u(1,:) = sin(2*pi*x);
for k = 1:nt-1
    for i = 1:nx-1
        % Predictor step
        u(k+1,i) = 2*t(k)*dt/dx*(u(k,i+1)-u(k,i)) + u(k,i);
    end    
   
end
figure(1);set(gcf,'Visible', 'off')
plot(x,u(85,:))
figure(4);set(gcf,'Visible', 'off')
surf(x,t,u)
%%%%%LAX WENDROFF%%%%%
dx2 = 0.01;
dt2 = 8e-4;
xdomain2 = [0 1];
tdomain2 = [0 1];
nx2 = round((xdomain2(2)-xdomain2(1))/dx2);
nt2 = round((tdomain2(2)-tdomain2(1))/dt2);
x2= linspace(xdomain2(1),xdomain2(2),nx2);
t2 = linspace(tdomain2(1),tdomain2(2),nt2);
u2 = zeros(nx2,nt2);
u2(:,1) = sin(2*pi*x2);%initial condition
for i=2:nx2-1
    for j=1:nt2-1
  u2(i,j+1)=u2(i,j)+(2*t2(j)*dt2/(2*dx2))*(u2(i+1,j)-u2(i-1,j))+((dt2^2)/(2*dx2))*(u2(i+1,j)-u2(i-1,j))+2*((dt2^2)/(dx2^2))*(t2(j)^2)*(u2(i+1,j)-2*u2(i,j)+u2(i-1,j));     
  stab(j)=2*t2(j)*dt2/(2*dx2);%always less then one
     end
end
figure(2);set(gcf,'Visible', 'off')
plot(x2,u2(:,85))
figure(3);set(gcf,'Visible', 'off')
surf(t2,x2,u2)