Iterative of Nonlinear Functions

26 visualizaciones (últimos 30 días)
Jamari
Jamari el 2 de Nov. de 2025
Editada: Torsten el 3 de Nov. de 2025
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
Error using catMany>checkDimensions (line 53)
Dimensions of objects being concatenated must match.

Error in catMany (line 12)
[resz, ranges] = checkDimensions(sz, dim);
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
Error in sym/horzcat (line 20)
y = catMany(2, args{:});
^^^^^^^^^^^^^^^^^^^^^^^^
% 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)]);

Respuestas (2)

Torsten
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:');
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)]);
1.0e+39 * -2.9940 0.0017 0.7202 0.0017 -2.9940 0.0023 -0.0000 -0.0005 -0.0000 0.0023 -0.0000 0.0000 -0.0000 0.0023 -0.0000
  3 comentarios
Torsten
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
s = 1×15
1.0e+197 * 3.6785 0.0000 0.0000 0.0000 3.6784 NaN 1.8393 0.0000 1.8392 NaN 3.6785 0.0000 0.0000 0.0000 3.6784
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
% 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
Torsten el 2 de Nov. de 2025
Editada: Torsten el 3 de Nov. de 2025
Systems like yours usually take the form
dT_i/dt = h_ext * (T_ext - T_i) + sum_{j~=i} h_ij * (T_j - T_i)
with a symmetric matrix h_ij.
I cannot find this structure in your equations.

Iniciar sesión para comentar.


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

Community Treasure Hunt

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

Start Hunting!

Translated by