Solving system of ODEs with ode45 and multiple time dependent variables

19 vues (au cours des 30 derniers jours)
I was given a system of 2 ODE's to solve for 2 different temperatures in the span of 24 hours. I was able to solve for the temperature for a constant Temperature Ta and constant I but in reality these values change for every hour of the day. I was given a table with the values of Ta and I for every hour t. I do not know how to implement the values into my function for every hour t. I am fairly new to matlab so I'm sure somebody out there can help. Here is my code so far. I will also provide a picture of the table I was given.
syms I h1 Tw Tg Mw Cw dTwdt Ta h2
eqn1 = 0.043 * I == h1 * (Tw-Tg) +(Mw*Cw*dTwdt)
DTwdtEqn = solve(eqn1,dTwdt)
eqn2 = 0.048 * I + h1*(Tw-Tg) == h2*(Tg-Ta)
TgEqn = solve(eqn2,Tg)
dTgdt = diff(TgEqn)
tRange = [0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23]
T0 = [293; 298];
[tSol, TSol] = ode45(@ODEfun,tRange,T0);
waterTemp = TSol(:,1);
glassTemp = TSol(:,2);
plot(tSol,waterTemp,"LineWidth",2)
hold on
plot(tSol,glassTemp,"LineWidth",2)
hold off
function dTdt = ODEfun(t,T)
Ta = 200;
I = 20;
Mw = 5;
Cw = 4180;
h1 = 30;
h2 = 100;
Tw = T(1);
Tg = T(2);
TwPrime = (I.*(43/100) +h1 *(Tg - Tw))/(Mw*Cw);
TgPrime = Tw/(h1+h2) - ( I.*(6/125) + Ta.* h2 + Tw *h1)/(h1+h2)^2;
dTdt = [TwPrime;TgPrime];
end

Réponse acceptée

Bjorn Gustavsson
Bjorn Gustavsson le 25 Fév 2021
In your ODE-fun you can "simply" interpolate your time-dependent coefficients. Something like this:
function dTdt = ODEfun(t,T,t4Ta,Ta,t4I,I)
Ta = interp1(t4Ta,Ta,t,'pchip'); % Sensiblest spline-interpolation scheme, I think
I = interp1(t4I,I,t,'pchip');
Mw = 5;
Cw = 4180;
h1 = 30;
h2 = 100;
Tw = T(1);
Tg = T(2);
TwPrime = (I.*(43/100) +h1 *(Tg - Tw))/(Mw*Cw);
TgPrime = Tw/(h1+h2) - ( I.*(6/125) + Ta.* h2 + Tw *h1)/(h1+h2)^2;
dTdt = [TwPrime;TgPrime];
end
Then you send in the relevant time-parameter-pairings for t4Ta, Ta, t4I and I:
[tSol, TSol] = ode45(@(t,T) ODEfun(t,T,0:23,Ta,0:23,I),tRange,T0);
HTH

Plus de réponses (0)

Catégories

En savoir plus sur Programming dans Help Center et File Exchange

Produits


Version

R2020b

Community Treasure Hunt

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

Start Hunting!

Translated by