Iterative of Nonlinear Functions
26 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
Someone please help, whenever i run my code I keep getting the same errors. My code and a screenshoot of the error page is attached.
% Given Values
k = 120;
alpha = 3.91e-6;
Ti=1073.15;
Tw=288.15;
Ta=288.15;
hw=100;
ha=10;
dt=10;
n_steps=60;
l=0.1;
tau=0.00391;
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;
T10=zeros(n_steps,1);
T10(1,1)=Ti;
T11=zeros(n_steps,1);
T11(1,1)=Ti;
T12=zeros(n_steps,1);
T12(1,1)=Ti;
T13=zeros(n_steps,1);
T13(1,1)=Ti;
T14=zeros(n_steps,1);
T14(1,1)=Ti;
T15=zeros(n_steps,1);
T15(1,1)=Ti;
for i=1:(n_steps-1)
syms t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 t11 t12 t13 t14 t15
eqn1 = l*k*(ha+hw)*(Ta-T1(i))+(T6-T1(i))+(T2(i)-T1(i))==(t1-T1(i))/tau;
eqn2 = ha*l*(Ta-T2(i)+((T7(i)-T2(i))+(T3(i)-T2(i))+(T1(i)-T2(i))))==(t2-T2(i))/tau;
eqn3 = ha*l/k*(Ta-T3(i))+((T4(i)-T3(i))+(T8(i)-T3(i)+(T2(i)-T3(i))))==(t3-T1(i))/tau;
eqn4 = ha*l*(Ta-T4(i)+((T9(i)-T4(i))+(T3(i)-T4(i))+(T5(i)-T4(i))))==(t4-T4(i))/tau;
eqn5 = l*k*(ha+hw)*(Ta-T5(i))+(T10-T5(i))+(T4(i)-T5(i))==(t5-T5(i))/tau;
eqn6 = hw*l/k*(Ta-T6(i))+((T1(i)-T6(i))+(T7(i)-T6(i))+(T11(i)-T6(i)))==(t6-T6(i))/tau;
eqn7 = (T6(i)-T7(i))+(T8(i)-T7(i))+(T2(i)-T7(i))+(T12(i)-T7(i))==(t7-T7(i))/tau;
eqn8 = (T7(i)-T8(i))+(T9(i)-T8(i))+(T3(i)-T8(i))+(T13(i)-T8(i))==(t8-T8(i))/tau;
eqn9 = (T8(i)-T9(i))+(T10(i)-T9(i))+(T4(i)-T9(i))+(T14(i)-T9(i))==(t9-T9(i))/tau;
eqn10 = hw*l/k*(Ta-T10(i))+((T5(i)-T10(i))+(T9(i)-T10(i))+(T15(i)-T10(i)))==(t10-T10(i))/tau;
eqn11 = hw*l/k*(Ta-T11(i))+((T6(i)-T11(i))+(T12(i)-T11(i)))==(t11-T11(i))/tau;
eqn12 = hw*l/k*(Ta-T12(i))+((T13(i)-T12(i)+(T7(i)-T12(i))+(T11(i)-T12(i))))==(t12-T12(i))/tau;
eqn13 = hw*l/k*(Ta-T13(i))+((T12(i)-T13(i)+(T8(i)-T13(i))+(T14(i)-T13(i))))==(t13-T13(i))/tau;
eqn14 = hw*l/k*(Ta-T14(i))+((T13(i)-T14(i)+(T5(i)-T14(i))+(T9(i)-T14(i))))==(t14-T14(i))/tau;
eqn15 = hw*l/k*(Ta-T15(i))+((T14(i)-T15(i))+(T10(i)-T15(i)))==(t15-T15(i))/tau;
s = solve([eqn1 eqn2 eqn3 eqn4 eqn5 eqn6 eqn7 eqn8 eqn9 eqn10 eqn11 eqn12 eqn13 eqn14 eqn15],...
[t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 t11 t12 t13 t14 t15]);
nodal_temperatures = double(struct2array(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);
T10(i+1,1) = s(10);
T11(i+1,1) = s(11);
T12(i+1,1) = s(12);
T13(i+1,1) = s(13);
T14(i+1,1) = s(14);
T15(i+1,1) = s(15);
end
% Display the results for the final temperatures
disp('Final Temperatures at Step 60:');
disp([T1(60,1), T2(60,1), T3(60,1), T4(60,1), T5(60,1), T6(60,1), ...
T7(60,1), T8(60,1), T9(60,1), T10(60,1), T11(60,1), T12(60,1),...
T13(60,1), T14(60,1), T15(60,1)]);
0 comentarios
Respuestas (2)
Torsten
el 2 de Nov. de 2025
Editada: Torsten
el 2 de Nov. de 2025
Why do you use "solve" to get t1,...,t15 ? You can directly solve for the unknowns because you use the explicit Euler method to advance in time. That's why your time step size tau might have to be chosen smaller than 0.00391 for that your temperatures don't diverge. Or you simply made a mistake in setting up your equations.
% Given Values
k = 120;
alpha = 3.91e-6;
Ti=1073.15;
Tw=288.15;
Ta=288.15;
hw=100;
ha=10;
dt=10;
n_steps=60;
l=0.1;
tau=0.00391;
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;
T10=zeros(n_steps,1);
T10(1,1)=Ti;
T11=zeros(n_steps,1);
T11(1,1)=Ti;
T12=zeros(n_steps,1);
T12(1,1)=Ti;
T13=zeros(n_steps,1);
T13(1,1)=Ti;
T14=zeros(n_steps,1);
T14(1,1)=Ti;
T15=zeros(n_steps,1);
T15(1,1)=Ti;
for i=1:(n_steps-1)
syms t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 t11 t12 t13 t14 t15
eqn1 = l*k*(ha+hw)*(Ta-T1(i))+(T6(i)-T1(i))+(T2(i)-T1(i))==(t1-T1(i))/tau;
eqn2 = ha*l*(Ta-T2(i)+((T7(i)-T2(i))+(T3(i)-T2(i))+(T1(i)-T2(i))))==(t2-T2(i))/tau;
eqn3 = ha*l/k*(Ta-T3(i))+((T4(i)-T3(i))+(T8(i)-T3(i)+(T2(i)-T3(i))))==(t3-T1(i))/tau;
eqn4 = ha*l*(Ta-T4(i)+((T9(i)-T4(i))+(T3(i)-T4(i))+(T5(i)-T4(i))))==(t4-T4(i))/tau;
eqn5 = l*k*(ha+hw)*(Ta-T5(i))+(T10(i)-T5(i))+(T4(i)-T5(i))==(t5-T5(i))/tau;
eqn6 = hw*l/k*(Ta-T6(i))+((T1(i)-T6(i))+(T7(i)-T6(i))+(T11(i)-T6(i)))==(t6-T6(i))/tau;
eqn7 = (T6(i)-T7(i))+(T8(i)-T7(i))+(T2(i)-T7(i))+(T12(i)-T7(i))==(t7-T7(i))/tau;
eqn8 = (T7(i)-T8(i))+(T9(i)-T8(i))+(T3(i)-T8(i))+(T13(i)-T8(i))==(t8-T8(i))/tau;
eqn9 = (T8(i)-T9(i))+(T10(i)-T9(i))+(T4(i)-T9(i))+(T14(i)-T9(i))==(t9-T9(i))/tau;
eqn10 = hw*l/k*(Ta-T10(i))+((T5(i)-T10(i))+(T9(i)-T10(i))+(T15(i)-T10(i)))==(t10-T10(i))/tau;
eqn11 = hw*l/k*(Ta-T11(i))+((T6(i)-T11(i))+(T12(i)-T11(i)))==(t11-T11(i))/tau;
eqn12 = hw*l/k*(Ta-T12(i))+((T13(i)-T12(i)+(T7(i)-T12(i))+(T11(i)-T12(i))))==(t12-T12(i))/tau;
eqn13 = hw*l/k*(Ta-T13(i))+((T12(i)-T13(i)+(T8(i)-T13(i))+(T14(i)-T13(i))))==(t13-T13(i))/tau;
eqn14 = hw*l/k*(Ta-T14(i))+((T13(i)-T14(i)+(T5(i)-T14(i))+(T9(i)-T14(i))))==(t14-T14(i))/tau;
eqn15 = hw*l/k*(Ta-T15(i))+((T14(i)-T15(i))+(T10(i)-T15(i)))==(t15-T15(i))/tau;
s = solve([eqn1 eqn2 eqn3 eqn4 eqn5 eqn6 eqn7 eqn8 eqn9 eqn10 eqn11 eqn12 eqn13 eqn14 eqn15],...
[t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 t11 t12 t13 t14 t15]);
s = double(struct2array(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);
T10(i+1,1) = s(10);
T11(i+1,1) = s(11);
T12(i+1,1) = s(12);
T13(i+1,1) = s(13);
T14(i+1,1) = s(14);
T15(i+1,1) = s(15);
end
% Display the results for the final temperatures
disp('Final Temperatures at Step 60:');
disp([T1(60,1), T2(60,1), T3(60,1), T4(60,1), T5(60,1), T6(60,1), ...
T7(60,1), T8(60,1), T9(60,1), T10(60,1), T11(60,1), T12(60,1), ...
T13(60,1), T14(60,1), T15(60,1)]);
3 comentarios
Torsten
el 2 de Nov. de 2025
Editada: Torsten
el 2 de Nov. de 2025
You get undefined results from "solve" for i >= 14:
% Given Values
k = 120;
alpha = 3.91e-6;
Ti=1073.15;
Tw=288.15;
Ta=288.15;
hw=100;
ha=10;
dt=10;
n_steps=60;
l=0.1;
tau=0.0039;
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;
T10=zeros(n_steps,1);
T10(1,1)=Ti;
T11=zeros(n_steps,1);
T11(1,1)=Ti;
T12=zeros(n_steps,1);
T12(1,1)=Ti;
T13=zeros(n_steps,1);
T13(1,1)=Ti;
T14=zeros(n_steps,1);
T14(1,1)=Ti;
T15=zeros(n_steps,1);
T15(1,1)=Ti;
%for i=1:(n_steps-1)
for i = 1:14
syms t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 t11 t12 t13 t14 t15
eqn1 = 2*l/k*((ha+hw)*(Ta-T1(i)))+2*(T6(i)+T2(i)+2*T1(i))==(t1-T1(i))/tau;
eqn2 = 2*ha*l/k*(Ta-T2(i))+2*(T7(i)-T2(i))+(T3(i)-T2(i))+(T1(i)-T2(i))==(t2-T2(i))/tau;
eqn3 = 2*ha*l/k*(Ta-T3(i))+(T2(i)-T3(i))+(T4(i)-T3(i))+2*(T8(i)-T3(i))==(t3-T3(i))/tau;
eqn4 = 2*ha*l/k*(Ta-T4(i))+2*(T9(i)-T4(i))+(T3(i)-T5(i))+(T5(i)-T4(i))==(t4-T4(i))/tau;
eqn5 = 2*l/k*((ha+hw)*(Ta-T5(i)))+2*(T10(i)+T4(i)+2*T5(i))==(t5-T5(i))/tau;
eqn6 = 2*hw*l/k*(Ta-T6(i))+2*(T7(i)*T6(i))+(T11(i)-T6(i))==(t6-T6(i))/tau;
eqn7 = (T2(i)-T7(i))+(T8(i)-T7(i))+(T12(i)-T7(i))+(T6(i)-T7(i))==(t7-T7(i))/tau;
eqn8 = (T3(i)-T8(i))+(T9(i)-T8(i))+(T13(i)-T8(i))+(T7(i)-T8(i))==(t8-T8(i))/tau;
eqn9 = (T4(i)-T9(i))+(T10(i)-T9(i))+(T14(i)-T9(i))+(T8(i)-T9(i))==(t9-T9(i))/tau;
eqn10 = 2*hw*l/k*(Ta-T10(i))+2*(T9(i)*T10(i))+(T15(i)-T10(i))==(t10-T10(i))/tau;
eqn11 = 4*hw*l/k*(Ta-T11(i))+2*(T6(i)-T11(i))+2*(T12(i)-T11(i))==(t11-T11(i))/tau;
eqn12 = 2*hw*l/k*(Ta-T12(i))+2*(T7(i)-T12(i))+(T13(i)-T12(i))+(T11(i)-T12(i))==(t12-T12(i))/tau;
eqn13 = 2*hw*l/k*(Ta-T13(i))+2*(T8(i)-T13(i))+(T14(i)-T13(i))+(T12(i)-T13(i))==(t13-T13(i))/tau;
eqn14 = 2*hw*l/k*(Ta-T14(i))+2*(T9(i)-T14(i))+(T15(i)-T14(i))+(T13(i)-T14(i))==(t14-T14(i))/tau;
eqn15 = 4*hw*l/k*(Ta-T15(i))+2*(T10(i)-T15(i))+2*(T14(i)-T15(i))==(t15-T15(i))/tau;
s = solve([eqn1 eqn2 eqn3 eqn4 eqn5 eqn6 eqn7 eqn8 eqn9 eqn10 eqn11 eqn12 eqn13 eqn14 eqn15],...
[t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 t11 t12 t13 t14 t15]);
s = double(struct2array(s));
if i==14
s
end
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);
T10(i+1,1) = s(10);
T11(i+1,1) = s(11);
T12(i+1,1) = s(12);
T13(i+1,1) = s(13);
T14(i+1,1) = s(14);
T15(i+1,1) = s(15);
end
% Display the results for the final temperatures
%disp('Final Temperatures at Step 60:');
%disp([T1(60,1), T2(60,1), T3(60,1), T4(60,1), T5(60,1), T6(60,1), ...
% T7(60,1), T8(60,1), T9(60,1), T10(60,1), T11(60,1), T12(60,1), ...
% T13(60,1), T14(60,1), T15(60,1)]);
Torsten
el 2 de Nov. de 2025
Movida: Torsten
el 2 de Nov. de 2025
My suggestion: Implicit Euler method.
And maybe there is a mistake in equation (3) (see below).
% Given Values
k = 120;
alpha = 3.91e-6;
Ti=1073.15;
Tw=288.15;
Ta=288.15;
hw=100;
ha=10;
dt=10;
n_steps=60;
l=0.1;
tau=0.00391;
Tnew = sym('Tnew',[15 1]);
Told = sym('Told',[15,1]);
eqn(1) = l*k*(ha+hw)*(Ta-Tnew(1))+(Tnew(6)-Tnew(1))+(Tnew(2)-Tnew(1))==(Tnew(1)-Told(1))/tau;
eqn(2) = ha*l*(Ta-Tnew(2)+((Tnew(7)-Tnew(2))+(Tnew(3)-Tnew(2))+(Tnew(1)-Tnew(2))))==(Tnew(2)-Told(2))/tau;
%eqn(3) = ha*l/k*(Ta-Tnew(3))+((Tnew(4)-Tnew(3))+(Tnew(8)-Tnew(3)+(Tnew(2)-Tnew(3))))==(Tnew(3)-Told(3))/tau;
eqn(3) = ha*l/k*(Ta-Tnew(3))+((Tnew(4)-Tnew(3))+(Tnew(8)-Tnew(3)+(Tnew(2)-Tnew(3))))==(Tnew(3)-Told(1))/tau;
eqn(4) = ha*l*(Ta-Tnew(4)+((Tnew(9)-Tnew(4))+(Tnew(3)-Tnew(4))+(Tnew(5)-Tnew(4))))==(Tnew(4)-Told(4))/tau;
eqn(5) = l*k*(ha+hw)*(Ta-Tnew(5))+(Tnew(10)-Tnew(5))+(Tnew(4)-Tnew(5))==(Tnew(5)-Told(5))/tau;
eqn(6) = hw*l/k*(Ta-Tnew(6))+((Tnew(1)-Tnew(6))+(Tnew(7)-Tnew(6))+(Tnew(11)-Tnew(6)))==(Tnew(6)-Told(6))/tau;
eqn(7) = (Tnew(6)-Tnew(7))+(Tnew(8)-Tnew(7))+(Tnew(2)-Tnew(7))+(Tnew(12)-Tnew(7))==(Tnew(7)-Told(7))/tau;
eqn(8) = (Tnew(7)-Tnew(8))+(Tnew(9)-Tnew(8))+(Tnew(3)-Tnew(8))+(Tnew(13)-Tnew(8))==(Tnew(8)-Told(8))/tau;
eqn(9) = (Tnew(8)-Tnew(9))+(Tnew(10)-Tnew(9))+(Tnew(4)-Tnew(9))+(Tnew(14)-Tnew(9))==(Tnew(9)-Told(9))/tau;
eqn(10) = hw*l/k*(Ta-Tnew(10))+((Tnew(5)-Tnew(10))+(Tnew(9)-Tnew(10))+(Tnew(15)-Tnew(10)))==(Tnew(10)-Told(10))/tau;
eqn(11) = hw*l/k*(Ta-Tnew(11))+((Tnew(6)-Tnew(11))+(Tnew(12)-Tnew(11)))==(Tnew(11)-Told(11))/tau;
eqn(12) = hw*l/k*(Ta-Tnew(12))+((Tnew(13)-Tnew(12)+(Tnew(7)-Tnew(12))+(Tnew(11)-Tnew(12))))==(Tnew(12)-Told(12))/tau;
eqn(13) = hw*l/k*(Ta-Tnew(13))+((Tnew(12)-Tnew(13)+(Tnew(8)-Tnew(13))+(Tnew(14)-Tnew(13))))==(Tnew(13)-Told(13))/tau;
eqn(14) = hw*l/k*(Ta-Tnew(14))+((Tnew(13)-Tnew(14)+(Tnew(5)-Tnew(14))+(Tnew(9)-Tnew(14))))==(Tnew(14)-Told(14))/tau;
eqn(15) = hw*l/k*(Ta-Tnew(15))+((Tnew(14)-Tnew(15))+(Tnew(10)-Tnew(15)))==(Tnew(15)-Told(15))/tau;
vars = Tnew;
[A,b] = equationsToMatrix(eqn,vars);
Anum = double(A);
Tsol = zeros(15,n_steps);
Tsol(:,1) = Ti*ones(15,1);
for i = 2:n_steps
bnum = double(subs(b,Told,Tsol(:,i-1)));
Tsol(:,i) = Anum\bnum;
end
plot(Tsol(:,end))
0 comentarios
Ver también
Categorías
Más información sobre Calculus en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!

