ode45 for M'=c1A'+ c2B' ?

vkjbjlj ckjjb,mn.

2015 7 3

0 个回答

6 次查看（30 天）

0 个投票

Hi MATLAB community I am solving an ODE using ode45. The equation is dM/dx= (c1*dA/dx)/A + (c2*dT/dx) where A and T are some function of x.I have written a code but this is not working.

function dMdx=scramjetwithoutshock1(x,M)
gamma_b=1.36;
x=0:0.001:0.35; 
xi=0; %x3=xi;
x4=0.35;
X=((x-xi)/(x4-xi));
a3=0.0038;
A=a3*(1+X);
a=diff(A)/A;
tau_b=1.4;
theta=2;
Ttx=1+(tau_b-1)*((theta*X)/(1+(theta-1)*X));
b=D(Ttx)/Ttx;
dMdx = M*((1+((gamma_b-1)/2)*M^2)/(1-M^2))*(-a+((1+(gamma_b*M^2))/2)*b); 
end

and recalled this function to plot M vs x as follows-

x=0:0.001:0.35;  %Combuster axial length variation
initial_M=2;     % M2=M3=2 (assumed)
[x,M]=ode45( @scramjetwithoutshock, x, initial_M);

Any help will be appreciated.

2 个评论
显示无隐藏无

RahulTandon 2015-7-6

编辑：James Tursa 2015-7-6

在 MATLAB Online 中打开

i tried to work out your problem.See the code below. But, if you could send me the pure M(x) and all the related equations, i will be able to solve the problem.

syms x a b a3 gamma_b theta tau_b M(x) c1 A(x) c2 T(x) Ttx(x);
str1 = 'DM == (c1*DA)/A + (c2*DT)' ;
str2 = 'A(x) == a3*(1 + x)';
str4 = 'a3 == DA/A';
str5 = 'Ttx(x) == 1+(tau_b-1)*((theta*x)/(1+(theta-1)*x))';
str6 = 'DTtx == b*Ttx';
str7 = 'DM == M*((1+((gamma_b-1)/2)*M^2)/(1-M^2))*(-a+((1+(gamma_b*M^2))/2)*b)';
str8 = 'b==DTtx/Ttx(x)';
clc;
M = dsolve(str1,str2,str3,str5,str6,str8,'x');   % str4, str7 are % omitted!