The bottom two equations in result are storing the derivitive of pressure and flowrate as state variables so that I can reference the derivitives in the first two equations. I'm not sure if this was the correct approach to this problem
Trouble using ODE45 for coupled non-linear differential equation
7 次查看(过去 30 天)
显示 更早的评论
Hi,
I am attempting to model the back filling of a large vacuum vessel with Nitrogen. The flow of the Nitrogen is controlled by the valve with a certain valve flow coefficient. P2 is constant and approximately one atmosphere. I have derivived the following equations to model this.
So I wrote the following code with ODE45 to monitor the pressure over time.
% IC Vector
IC_BF = [P01, Q0, dP0, dQ0];
t = [0 10];
options = odeset('RelTol',1e-12 );
[t_ode45,Result] = ode45(@BF_dyn, t, IC_BF, options)
function result = BF_dyn(t, x)
% Constants
R = 8.314; % J / mol·K
C_v = 0.28*(6.309e-5/1); % (gallons/min) ( 6.309e-5 (m^3/s)/ 1 (gallons/min)) = m^3/s
Tamb = 293; % K
G = 0.967;
dewar_vol = 4.4; % m^3
rho_N2 = 1.25; %kg/m^3
M_N2 = 28.02*(1/1000); % g/mol * (1 kg/1000g) = kg/mol
P02 = 6894.76; % Pa
b = (sqrt(G)/(C_v))^2;
a = rho_N2*R*Tamb/(M_N2*dewar_vol);
result = [
-2*x(2)*x(4)*b;
(1/a)*x(4);
x(3);
x(4);
];
end
However, when I plot the results I get the following, which goes significantly higher than P2 before attempting to go negative when in reality it should asymptote out to P2.
Am I implement this system into ODE45 wrong? Is there any way to incorporate the realtionship between P1 and P2 into ODE45 by adding additional arguement to my results vector?
Thank you for any assistance you can offer!
回答(1 个)
Torsten
2024-10-2
编辑:Walter Roberson
2024-10-2
I can't recover the four differential equations you solve from the two equations you included in mathematical form.
Your equations can be solved analytically:
syms a b P01 Q0 dP0 dQ0
syms t x1(t) x2(t) x3(t) x4(t)
odes = [diff(x1,t) == -2*b*x2*x4,diff(x2,t) ==x4/a,diff(x3,t)==x3,diff(x4,t)==x4];
conds = [x1(0)==P01,x2(0)==Q0,x3(0)==dP0,x4(0)==dQ0];
sol = dsolve(odes,conds);
sol.x1
sol.x2
sol.x3
sol.x4
1 个评论
Walter Roberson
2024-10-2
(Unfortunately due to a glitch, the output is not rendering here on MATLAB Answers. If you run the code on your own display then it will work.)
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!