How to solve interdependent differtional equations

Question

Ahmed Aburakhia 2023-12-9

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2058299-how-to-solve-interdependent-differtional-equations

评论： Sulaymon Eshkabilov 2023-12-9

在 MATLAB Online 中打开

Hi

I have a system of three differetional equations

I defined my function as follow

function dydx = StratNoObs(x, y)    
D=1; %depth [m]
rho_bs=2000;%[kg/m3]
rho_ts=1000;%[kg/m3]
drho_dy=rho_bs-rho_ts/D;
rho_t = rho_ts/(rho_bs*drho_dy); % 
rho_b = rho_bs/(rho_bs*drho_dy); % [kg/m3]
U=1;
Ri = (9.81*D^2*drho_dy)/(U^2*rho_bs); 
Re = 100; 
diffusifty = 0.000025; %[m2/s]
mu =  0.001002; % [kg/(m.s)]
Sc= (mu/(rho_bs*diffusifty)); %the ratio of the coefficient of diffusion to the kinematic viscosity
d_rho=rho_t/rho_b;
b = y(1);
f = y(2);
e = y(3);
dbdx = (280.*e+(1400/3).*f)./(Ri*Re);
fA = -(2.*f)./(Re);
fB = -Ri*(-4/735).*dbdx;
fC = (1/Re).*(-(1/3).*e-(2/9).*f);
fD = -((13/11340).*e+(31/17010).*f+(1/90));
fF = -dbdx*((-21*.e-41*.f-1404)./(882.*(rho_b-rho_t+(1/42).*b)));
fG = -(1080.*b)./(Re.*Sc.*(42.*rho_b-42.*rho_t+(1/42).*b));
fI = -(1470.*rho_b-1470.*rho_t+41.*b)/(882.*rho_b-882.*rho_t+21.*b);
fJ = ((31/17010).*e+(76/25515).*f+(1/54));
dfdx = (fA+fB+fC+fD.*(fF+fG)./(-fD.*fI+fJ);
dedx = fI.*dfdx+fF+fG;
dydx = [dbdx; dfdx; dedx];
end

I am not getting the expected answer so is I am questioning if it is ok to call the dbdx and dfdx in the defintion dedx?

Thank you

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Answer 1

Sulaymon Eshkabilov 2023-12-9

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2058299-how-to-solve-interdependent-differtional-equations#answer_1368429

在 MATLAB Online 中打开

Here is the corrected code. There were a few syntx errs and missing paranthesis.

x= linspace(0, pi);
y = [0 1 -1];
SOL=StratNoObs(x, y)   
SOL = 3×1
1.0e+05 *

    0.0000
   -4.9307
    8.2301
function dydx = StratNoObs(x, y)    
D=1; %depth [m]
rho_bs=2000;%[kg/m3]
rho_ts=1000;%[kg/m3]
drho_dy=rho_bs-rho_ts/D;
rho_t = rho_ts/(rho_bs*drho_dy); % 
rho_b = rho_bs/(rho_bs*drho_dy); % [kg/m3]
U=1;
Ri = (9.81*D^2*drho_dy)/(U^2*rho_bs); 
Re = 100; 
diffusifty = 0.000025; %[m2/s]
mu =  0.001002; % [kg/(m.s)]
Sc= (mu/(rho_bs*diffusifty)); %the ratio of the coefficient of diffusion to the kinematic viscosity
d_rho=rho_t/rho_b;
b = y(1);
f = y(2);
e = y(3);
dbdx = (280.*e+(1400/3).*f)./(Ri*Re);
fA = -(2.*f)./(Re);
fB = -Ri*(-4/735).*dbdx;
fC = (1/Re).*(-(1/3).*e-(2/9).*f);
fD = -((13/11340).*e+(31/17010).*f+(1/90));
fF = -dbdx*((-21.*e-41.*f-1404)./(882.*(rho_b-rho_t+(1/42).*b)));
fG = -(1080.*b)./(Re.*Sc.*(42.*rho_b-42.*rho_t+(1/42).*b));
fI = -(1470.*rho_b-1470.*rho_t+41.*b)/(882.*rho_b-882.*rho_t+21.*b);
fJ = ((31/17010).*e+(76/25515).*f+(1/54));
dfdx = (fA+fB+fC+fD.*(fF+fG))./(-fD.*fI+fJ);
dedx = fI.*dfdx+fF+fG;
dydx = [dbdx; dfdx; dedx];
end

2 个评论
显示无隐藏无

Ahmed Aburakhia 2023-12-9

thanks for correcting my code. I am more concern with the methodolagy as I am not famialir with a system of ODEs. I took the matlab course but their ODEs are not interdependant. In my case I am using dbdx in solving for dfdx and using both dbdx and dfdx in solving for dedx. Is this the correct way?

Sulaymon Eshkabilov 2023-12-9