Trouble using ODE45 for coupled non-linear differential equation

2 次查看（过去 30 天）
Jacob
Jacob 2024-10-2
编辑： Walter Roberson 2024-10-2
Ran in:
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];
Unrecognized function or variable 'P01'.
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!
  3 个评论
显示 1更早的评论
Sam Chak
Sam Chak 2024-10-2
@Jacob, For what purpose do you put these unused constants inside the BF_dyn(t, x) function?
delta_P0= 6894.76; % Pa
P01 = 1.33322e-5; % pa
dP01 = 0; % Pa
P02 = 6894.76; % Pa
Torsten
Torsten 2024-10-2
编辑：Torsten 2024-10-2
If you equate the two equations you posted, dQ/dt cancels out and you get
Q(t) = -R*T*M*rho/V / (2*G/Cv^2)
Does that make sense ?

请先登录，再进行评论。

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

回答（1 个）

Torsten
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
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.)

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Programming 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by