Iteration of multiple nonlinear functions
4 views (last 30 days)
Show older comments
I have an assignment where I need to clalculate the temperature distribution in a window with elctric heaters installed at nine different points.
I have made the correct temperature profile equation at each point. Now I need using iteration to find the temperature at a given time, and that is where my problem starts.
My code looks like this so far:
%Given values
k = 0.84;
Qh = 25;
Ti = -3;
To = -3;
hi = 6;
ho = 20;
dx = 0.002;
dy = 0.01;
n_steps = 1000;
tau = 0.078;
T1 = zeros(n_steps,1);
T1(1,1) = Ti;
T2 = zeros(n_steps,1);
T2(1,1) = Ti;
T3 = zeros(n_steps,1);
T3(1,1) = Ti;
T4 = zeros(n_steps,1);
T4(1,1) = Ti;
T5= zeros(n_steps,1);
T5(1,1) = Ti;
T6 = zeros(n_steps,1);
T6(1,1) = Ti;
T7 = zeros(n_steps,1);
T7(1,1) = Ti;
T8 = zeros(n_steps,1);
T8(1,1) = Ti;
T9 = zeros(n_steps,1);
T9(1,1) = Ti;
for i = 1:(n_steps-1)
syms T1 T2 T_ T4 T5 T6 T7 T8 T9
eqn1 = T1 == (1-2*hi*dy*tau/k-10.4*tau)*T1(i) + 2*hi*dy*tau*Ti/k + 10*tau*T2(i) + 0.4*tau*T4(i);
eqn2 = T2 == (1-10.4*tau)*T2(i) + 5*(T1(i)+T3(i))*tau + 0.4*tau*T5(i);
eqn3 = T3 == (1-2*ho*dy*tau/k-10.4*tau)*T3(i) + 2*ho*dy*tau*To/k + 10*tau*T2(i)+0.4*tau*T6(i);
eqn4 = T4 == (1-2*hi*dy*tau/k-10.4*tau)*T4(i) + 2*hi*dy*tau*Ti/k + 0.2*tau*(T1(i)+T7(i)) + 10*tau*T5(i);
eqn5 = T5 == (1-10.4*tau)*T5(i) + 0.2*(T2(i)+T8(i))*tau + 5*(T4(i)+T6(i))*tau;
eqn6 = T6 == (1-2*ho*dy*tau/k-10.4*tau)*T6(i) + 2*ho*dy*tau*To/k+0.2*(T3(i)+T9(i))*tau + 10*tau*T5(i);
eqn7 = T7 == (1-2*hi*dy*tau/k-10.4*tau)*T7(i) + 2*Qh*tau/k+2*hi*dy*tau*Ti/k+0.4*tau*T4(i) + 10*tau*T8(i);
eqn8 = T8 == (1-10.4*tau)*T8(i) + 0.4*tau*T5(i) + 5*(T7(i)+T9(i))*tau;
eqn9 = T9 == (1-2*ho*dy*tau/k-10.4*tau)*T9(i) + 2*ho*dy*tau*To/k + 0.4*tau*T6(i) + 10*T8(i)*tau;
s = solve([eqn1 eqn2 eqn3 eqn4 eqn5 eqn6 eqn7 eqn8 eqn9],[T1 T2 T3 T4 T5 T6 T7 T8 T9])
s = struct2array(s);
s = double(s);
T1(i+1,1) = s(1);
T2(i+1,1) = s(2);
T3(i+1,1) = s(3);
T4(i+1,1) = s(4);
T5(i+1,1) = s(5);
T6(i+1,1) = s(6);
T7(i+1,1) = s(7);
T8(i+1,1) = s(8);
T9(i+1,1) = s(9);
end
%Temperature after 15 min
res_15min = [T1(225,1);T2(225,1);T3(225,1);T4(225,1);T5(225,1);T6(225,1);T7(225,1);T8(225,1);T9(225,1)]
%Temperature at steady-state
res_ss = [T1(525,1);T2(525,1);T3(525,1);T4(525,1);T5(525,1);T6(525,1);T7(525,1);T8(525,1);T9(525,1)]
I am unsure about the command "solve".
Does anyone know a better way of iterating these equation and then later plot them?
0 Comments
Answers (2)
Torsten
on 20 Jan 2023
%Given values
k = 0.84;
Qh = 25;
Ti = -3;
To = -3;
hi = 6;
ho = 20;
dx = 0.002;
dy = 0.01;
n_steps = 1000;
tau = 0.078;
T1 = zeros(n_steps,1);
T1(1,1) = Ti;
T2 = zeros(n_steps,1);
T2(1,1) = Ti;
T3 = zeros(n_steps,1);
T3(1,1) = Ti;
T4 = zeros(n_steps,1);
T4(1,1) = Ti;
T5 = zeros(n_steps,1);
T5(1,1) = Ti;
T6 = zeros(n_steps,1);
T6(1,1) = Ti;
T7 = zeros(n_steps,1);
T7(1,1) = Ti;
T8 = zeros(n_steps,1);
T8(1,1) = Ti;
T9 = zeros(n_steps,1);
T9(1,1) = Ti;
for i = 1:(n_steps-1)
T1(i+1,1) = (1-2*hi*dy*tau/k-10.4*tau)*T1(i) + 2*hi*dy*tau*Ti/k + 10*tau*T2(i) + 0.4*tau*T4(i);
T2(i+1,1) = (1-10.4*tau)*T2(i) + 5*(T1(i)+T3(i))*tau + 0.4*tau*T5(i);
T3(i+1,1)= (1-2*ho*dy*tau/k-10.4*tau)*T3(i) + 2*ho*dy*tau*To/k + 10*tau*T2(i)+0.4*tau*T6(i);
T4(i+1,1)= (1-2*hi*dy*tau/k-10.4*tau)*T4(i) + 2*hi*dy*tau*Ti/k + 0.2*tau*(T1(i)+T7(i)) + 10*tau*T5(i);
T5(i+1,1)= (1-10.4*tau)*T5(i) + 0.2*(T2(i)+T8(i))*tau + 5*(T4(i)+T6(i))*tau;
T6(i+1,1)= (1-2*ho*dy*tau/k-10.4*tau)*T6(i) + 2*ho*dy*tau*To/k+0.2*(T3(i)+T9(i))*tau + 10*tau*T5(i);
T7(i+1,1)= (1-2*hi*dy*tau/k-10.4*tau)*T7(i) + 2*Qh*tau/k+2*hi*dy*tau*Ti/k+0.4*tau*T4(i) + 10*tau*T8(i);
T8(i+1,1)= (1-10.4*tau)*T8(i) + 0.4*tau*T5(i) + 5*(T7(i)+T9(i))*tau;
T9(i+1,1)= (1-2*ho*dy*tau/k-10.4*tau)*T9(i) + 2*ho*dy*tau*To/k + 0.4*tau*T6(i) + 10*T8(i)*tau;
end
plot((1:n_steps).',[T1,T2,T3,T4,T5,T6,T7,T8,T9])
2 Comments
Sam Chak
on 30 May 2024
If the beautiful solution you needed has been provided, please consider clicking 'Accept' ✔ on the answer to close the issue. Additionally, voting for @Alan Stevens' helpful answer serves as tokens of appreciation.
Alan Stevens
on 20 Jan 2023
Here's a quick and dirty way. I'll leave you to modify the following to record all the temperatures at every step.
%Given values
k = 0.84;
Qh = 25;
Ti = -3;
To = -3;
hi = 6;
ho = 20;
dx = 0.002;
dy = 0.01;
n_steps = 1000;
tau = 0.078;
V = [2*hi*dy*tau*Ti/k; 0; 2*ho*dy*tau*To/k; 2*hi*dy*tau*Ti/k; 0; 2*ho*dy*tau*To/k;
2*Qh*tau/k+2*hi*dy*tau*Ti/k; 0; 2*ho*dy*tau*To/k];
M = [(1-2*hi*dy*tau/k-10.4*tau), 10*tau, 0, 0.4*tau, 0, 0, 0, 0, 0;
5*tau, (1-10.4*tau), 5*tau, 0, 0.4*tau, 0,0,0,0;
0, 10*tau, (1-2*ho*dy*tau/k-10.4*tau), 0, 0, 0.4*tau, 0,0,0;
0.2*tau, 0, 0, (1-2*hi*dy*tau/k-10.4*tau), 10*tau, 0,0.2*tau,0,0;
0, 0.2*tau, 0, 5*tau, (1-10.4*tau), 5*tau, 0,0,0;
0, 0, 0.2*tau, 0, 10*tau, (1-2*ho*dy*tau/k-10.4*tau), 0, 0, 0.2*tau;
0, 0, 0, 0.4*tau, 0, 0, (1-2*hi*dy*tau/k-10.4*tau), 10*tau, 0;
0, 0, 0, 0, 0.4*tau, 0, 5*tau, (1-10.4*tau), 5*tau;
0, 0, 0, 0, 0, 0.4*tau, 0, 10*tau, (1-2*ho*dy*tau/k-10.4*tau)];
T = Ti*ones(9,1);
for i = 1:(n_steps-1)
T = M*T + V;
if i == 225
disp('Temps after 15mins')
disp(T')
figure
plot(1:9,T,'o-')
end
if i == 525
disp('Steady-state temps')
disp(T')
figure
plot(1:9,T,'o-')
end
end
0 Comments
See Also
Categories
Find more on Calculus in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!