Midterm Project

Green House
MATLAB Codes
Question: 03
% Wavelength range (µm) and solar intensity (W/m²/µm)

wavelength = linspace(0.3, 2.5, 100); % Example wavelengths from 0.3 to 2.5 µm

I_solar = exp(-0.5 * (wavelength - 1).^2) * 1000; % Example solar spectrum (simplified)

% Spectral properties for low iron glass

R_low_iron = 0.04 + 0.005 * (wavelength - 0.3); % Reflectance for low iron glass (example)
T_low_iron = 0.92 - 0.03 * (wavelength - 0.3); % Transmittance for low iron glass (example)

A_low_iron = 1 - R_low_iron - T_low_iron; % Absorptance for low iron glass

% Spectral properties for high iron glass

R_high_iron = 0.08 + 0.01 * (wavelength - 0.3); % Reflectance for high iron glass (example)

T_high_iron = 0.85 - 0.05 * (wavelength - 0.3); % Transmittance for high iron glass (example)

A_high_iron = 1 - R_high_iron - T_high_iron; % Absorptance for high iron glass

% Calculate energy contributions for low iron glass

E_reflected_low_iron = R_low_iron .* I_solar;

E_transmitted_low_iron = T_low_iron .* I_solar;

E_absorbed_low_iron = A_low_iron .* I_solar;

% Calculate energy contributions for high iron glass

E_reflected_high_iron = R_high_iron .* I_solar;

E_transmitted_high_iron = T_high_iron .* I_solar;

E_absorbed_high_iron = A_high_iron .* I_solar;

% Integrate total energies using the trapezoidal rule

E_total = trapz(wavelength, I_solar); % Total incident energy

% Total energies for low iron glass

E_reflected_total_low_iron = trapz(wavelength, E_reflected_low_iron);

E_transmitted_total_low_iron = trapz(wavelength, E_transmitted_low_iron);

E_absorbed_total_low_iron = trapz(wavelength, E_absorbed_low_iron);

% Total energies for high iron glass

E_reflected_total_high_iron = trapz(wavelength, E_reflected_high_iron);

E_transmitted_total_high_iron = trapz(wavelength, E_transmitted_high_iron);

E_absorbed_total_high_iron = trapz(wavelength, E_absorbed_high_iron);

% Calculate ratios for low iron glass

R_ratio_low_iron = E_reflected_total_low_iron / E_total;

T_ratio_low_iron = E_transmitted_total_low_iron / E_total;

A_ratio_low_iron = E_absorbed_total_low_iron / E_total;

% Calculate ratios for high iron glass

R_ratio_high_iron = E_reflected_total_high_iron / E_total;

T_ratio_high_iron = E_transmitted_total_high_iron / E_total;

A_ratio_high_iron = E_absorbed_total_high_iron / E_total;

% Display results for low iron glass

fprintf('Low Iron Glass:\n');

fprintf('Reflected Ratio: %.2f\n', R_ratio_low_iron);

fprintf('Transmitted Ratio: %.2f\n', T_ratio_low_iron);

fprintf('Absorbed Ratio: %.2f\n\n', A_ratio_low_iron);

% Display results for high iron glass

fprintf('High Iron Glass:\n');

fprintf('Reflected Ratio: %.2f\n', R_ratio_high_iron);

fprintf('Transmitted Ratio: %.2f\n', T_ratio_high_iron);

fprintf('Absorbed Ratio: %.2f\n', A_ratio_high_iron);

% Plot results for low iron glass

figure;
set(gcf, 'Position', [100, 100, 900, 400]); % Set the figure size (width x height)

subplot(1, 2, 1);

plot(wavelength, R_low_iron, 'r-', 'LineWidth', 1.5); hold on;

plot(wavelength, T_low_iron, 'g-', 'LineWidth', 1.5);

plot(wavelength, A_low_iron, 'b-', 'LineWidth', 1.5);

xlabel('Wavelength (µm)');

ylabel('Ratio');

legend('Reflectance', 'Transmittance', 'Absorptance');

title('Low Iron Glass Spectral Properties');

grid on;

yticks(0:0.2:1); % Set the y-axis ticks from 0 to 1 with a step of 0.2

hold off;
% Plot results for high iron glass

subplot(1, 2, 2);

plot(wavelength, R_high_iron, 'r-', 'LineWidth', 1.5); hold on;

plot(wavelength, T_high_iron, 'g-', 'LineWidth', 1.5);

plot(wavelength, A_high_iron, 'b-', 'LineWidth', 1.5);

xlabel('Wavelength (µm)');

ylabel('Ratio');
legend('Reflectance', 'Transmittance', 'Absorptance');

title('High Iron Glass Spectral Properties');

grid on;

yticks(0:0.2:1); % Set the y-axis ticks from 0 to 1 with a step of 0.2

hold off;
Results

Graphs

Question: 04
Case: 01 (Summer Daytime: Windspeed=0 m/s, Tsky=235 K)
Low-Iron glass

T_inf = 25 + 273.15;

T_s = 80 + 273.15;

L = 3;

T_avg = (T_inf + T_s) / 2;

k = 0.026;
nu = 15.89e-6;
alpha = 22.14e-6;

beta = 1 / T_avg;

g = 9.81;

Pr = nu / alpha;
Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

delta_x = 0.3;

k_material = 1;

A = 36;

R1 = (delta_x / (k_material * A)) + (1 / (h * A));

R2 = (0.005 / (2 * 18)) + (1 / (h * 18));

R = (1 / R1 + 1 / R2) ^ -1;

I0 = 1000 * 0.85;

t = linspace(0, 15, 100);

Delta_T_t = zeros(size(t));

T_s_time = zeros(size(t));

for i = 1:length(t)
T_avg = (T_inf + T_s) / 2;

beta = 1 / T_avg;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

q_t = I0 * sin(t(i) * pi / 15) * 18;

Delta_T_t(i) = q_t * R;
T_s = T_inf + Delta_T_t(i);

T_s_time(i) = T_s - 273.15;

end

figure;
plot(t, Delta_T_t + 25, 'LineWidth', 2);

xlabel('Time (s)');

ylabel('Room Temperature Variation (K)');

title('Room Temperature Variation vs Time');

grid on;

High-Iron glass

T_inf = 25 + 273.15;

T_s = 80 + 273.15;

L = 3;

T_avg = (T_inf + T_s) / 2;

k = 0.026;

nu = 15.89e-6;

alpha = 22.14e-6;

beta = 1 / T_avg;

g = 9.81;
Pr = nu / alpha;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

delta_x = 0.3;

k_material = 1;

A = 36;

R1 = (delta_x / (k_material * A)) + (1 / (h * A));

R2 = (0.005 / (2 * 18)) + (1 / (h * 18));

R = (1 / R1 + 1 / R2) ^ -1;

I0 = 1000 * 0.7;
t = linspace(0, 15, 100);

Delta_T_t = zeros(size(t));

T_s_time = zeros(size(t));

for i = 1:length(t)

T_avg = (T_inf + T_s) / 2;

beta = 1 / T_avg;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

q_t = I0 * sin(t(i) * pi / 15) * 18;

Delta_T_t(i) = q_t * R;

T_s = T_inf + Delta_T_t(i);

T_s_time(i) = T_s - 273.15;

end
figure;

plot(t, Delta_T_t + 25, 'LineWidth', 2);

xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('High-Iron Glass (day: v=0 m/s and Tsky=235 K)');

grid on;

Case: 02 (Summer Daytime: Windspeed=0 m/s, Tsky=285 K)

Low-Iron glass
T_inf = 25 + 273.15;

T_s = 80 + 273.15;

L = 3;

T_avg = (T_inf + T_s) / 2;

k = 0.026;
nu = 15.89e-6;

alpha = 22.14e-6;

beta = 1 / T_avg;

g = 9.81;
Pr = nu / alpha;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

delta_x = 0.3;

k_material = 1;

A = 36;

R1 = (delta_x / (k_material * A)) + (1 / (h * A));

R2 = (0.005 / (2 * 18)) + (1 / (h * 18));

R = (1 / R1 + 1 / R2) ^ -1;

I0 = 1000 * 0.85;
t = linspace(0, 15, 100);

Delta_T_t = zeros(size(t));

T_s_time = zeros(size(t));

for i = 1:length(t)

T_avg = (T_inf + T_s) / 2;

beta = 1 / T_avg;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

q_t = I0 * sin(t(i) * pi / 15) * 18;

Delta_T_t(i) = q_t * R;

T_s = T_inf + Delta_T_t(i);

T_s_time(i) = T_s - 273.15;

end
figure;

plot(t, Delta_T_t + 25, 'LineWidth', 2);

xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('Low-Iron glass (day: v=0 m/s, Tsky=285 K)');

grid on;

High-Iron glass
T_inf = 25 + 273.15;

T_s = 80 + 273.15;

L = 3;

T_avg = (T_inf + T_s) / 2;

k = 0.026;

nu = 15.89e-6;

alpha = 22.14e-6;
beta = 1 / T_avg;

g = 9.81;

Pr = nu / alpha;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

delta_x = 0.3;
k_material = 1;

A = 36;

R1 = (delta_x / (k_material * A)) + (1 / (h * A));

R2 = (0.005 / (2 * 18)) + (1 / (h * 18));

R = (1 / R1 + 1 / R2) ^ -1;

I0 = 1000 * 0.7;

t = linspace(0, 15, 100);

Delta_T_t = zeros(size(t));
T_s_time = zeros(size(t));

for i = 1:length(t)

T_avg = (T_inf + T_s) / 2;

beta = 1 / T_avg;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

q_t = I0 * sin(t(i) * pi / 15) * 18;

Delta_T_t(i) = q_t * R;

T_s = T_inf + Delta_T_t(i);

T_s_time(i) = T_s - 273.15;

end

figure;
plot(t, Delta_T_t + 25, 'LineWidth', 2);
xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('High-Iron glass (day: v=0 m/s, Tsky=285 K)');

grid on;

Case: 03 (Summer Daytime: Windspeed=10 m/s, Tsky=235 K)

Low-Iron glass
T_inf=25+273.15;

T_s=80+273.15;

L=3;
T_avg=(T_inf+T_s)/2;

k=0.026;

nu=15.89e-6;

alpha=22.14e-6;

beta=1/T_avg;

g=9.81;

Pr=nu/alpha;
U_inf=10;
Re_L=U_inf*L/nu;

if Re_L<5e5

Nu_avg=0.664*Re_L^0.5*Pr^(1/3);

else
Nu_avg=0.037*Re_L^(4/5)*Pr^(1/3);

end

h=Nu_avg*k/L;

delta_x=0.3;

k_material=1;

A=36;

R1=(delta_x/(k_material*A))+(1/(h*A));

R2=(0.005/(2*18))+(1/(h*18));
R=(1/R1+1/R2)^-1;

I0=1000*0.85;

t=linspace(0,15,100);

Delta_T_t=zeros(size(t));

T_s_time=zeros(size(t));

for i=1:length(t)

T_avg=(T_inf+T_s)/2;
beta=1/T_avg;

Gr=(g*beta*(T_s-T_inf)*L^3)/(nu^2);

Ra=Gr*Pr;

Nu=0.68+(0.67*Ra^(1/4))/((1+(0.492/Pr)^(9/16))^(4/9));

h=(Nu*k)/L;

q_t=I0*sin(t(i)*pi/15)*18;

Delta_T_t(i)=q_t*R;
T_s=T_inf+Delta_T_t(i);
T_s_time(i)=T_s-273.15;

end

figure;

plot(t,Delta_T_t +25,'LineWidth',2);
xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('Low-class iron (day: 10 m/s, Tsky=235 K)');

grid on;

High-Iron glass
T_inf=25+273.15;

T_s=80+273.15;
L=3;

T_avg=(T_inf+T_s)/2;

k=0.026;

nu=15.89e-6;

alpha=22.14e-6;

beta=1/T_avg;

g=9.81;
Pr=nu/alpha;

U_inf=10;

Re_L=U_inf*L/nu;

if Re_L<5e5
Nu_avg=0.664*Re_L^0.5*Pr^(1/3);

else

Nu_avg=0.037*Re_L^(4/5)*Pr^(1/3);

end

h=Nu_avg*k/L;

delta_x=0.3;

k_material=1;

A=36;
R1=(delta_x/(k_material*A))+(1/(h*A));

R2=(0.005/(2*18))+(1/(h*18));

R=(1/R1+1/R2)^-1;

I0=1000*0.7;

t=linspace(0,15,100);

Delta_T_t=zeros(size(t));

T_s_time=zeros(size(t));
for i=1:length(t)

T_avg=(T_inf+T_s)/2;

beta=1/T_avg;

Gr=(g*beta*(T_s-T_inf)*L^3)/(nu^2);

Ra=Gr*Pr;

Nu=0.68+(0.67*Ra^(1/4))/((1+(0.492/Pr)^(9/16))^(4/9));

h=(Nu*k)/L;
q_t=I0*sin(t(i)*pi/15)*18;
Delta_T_t(i)=q_t*R;

T_s=T_inf+Delta_T_t(i);

T_s_time(i)=T_s-273.15;

end
figure;

plot(t,Delta_T_t +25,'LineWidth',2);

xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('High-class iron (day: 10 m/s, Tsky=235 K)');

grid on;

Case: 04 (Summer Daytime: Windspeed=10 m/s, Tsky=285 K)

Low-Iron glass
T_inf=25+273.15; T_s=80+273.15; L=3; T_avg=(T_inf+T_s)/2; k=0.026; nu=15.89e-6;
alpha=22.14e-6; beta=1/T_avg; g=9.81; Pr=nu/alpha; U_inf=10; Re_L=U_inf*L/nu;

if Re_L<5e5

Nu_avg=0.664*Re_L^0.5*Pr^(1/3);

else

Nu_avg=0.037*Re_L^(4/5)*Pr^(1/3);
end

h=Nu_avg*k/L;

delta_x=0.3; k_material=1; A=36;

R1=(delta_x/(k_material*A))+(1/(h*A)); R2=(0.005/(2*18))+(1/(h*18)); R=(1/R1+1/R2)^-1;

I0=1000*0.85; t=linspace(0,15,100); Delta_T_t=zeros(size(t)); T_s_time=zeros(size(t));

for i=1:length(t)

T_avg=(T_inf+T_s)/2; beta=1/T_avg; Gr=(g*beta*(T_s-T_inf)*L^3)/(nu^2); Ra=Gr*Pr;

Nu=0.68+(0.67*Ra^(1/4))/((1+(0.492/Pr)^(9/16))^(4/9)); h=(Nu*k)/L;

q_t=I0*sin(t(i)*pi/15)*18; Delta_T_t(i)=q_t*R;

T_s=T_inf+Delta_T_t(i); T_s_time(i)=T_s-273.15;

end

figure; plot(t,Delta_T_t+25,'LineWidth',2);
xlabel('Time (Hr)');

ylabel(' Room Temperature Variation (K)');

title('Low-class iron (day: 10 m/s, Tsky=285 K)');

grid on;
High-Iron glass

T_inf=25+273.15; T_s=80+273.15; L=3; T_avg=(T_inf+T_s)/2; k=0.026; nu=15.89e-6;

alpha=22.14e-6; beta=1/T_avg; g=9.81; Pr=nu/alpha; U_inf=10; Re_L=U_inf*L/nu;

if Re_L<5e5

Nu_avg=0.664*Re_L^0.5*Pr^(1/3);

else

Nu_avg=0.037*Re_L^(4/5)*Pr^(1/3);
end

h=Nu_avg*k/L;
delta_x=0.3; k_material=1; A=36;

R1=(delta_x/(k_material*A))+(1/(h*A)); R2=(0.005/(2*18))+(1/(h*18)); R=(1/R1+1/R2)^-1;

I0=1000*0.7; t=linspace(0,15,100); Delta_T_t=zeros(size(t)); T_s_time=zeros(size(t));

for i=1:length(t)

T_avg=(T_inf+T_s)/2; beta=1/T_avg; Gr=(g*beta*(T_s-T_inf)*L^3)/(nu^2); Ra=Gr*Pr;

Nu=0.68+(0.67*Ra^(1/4))/((1+(0.492/Pr)^(9/16))^(4/9)); h=(Nu*k)/L;

q_t=I0*sin(t(i)*pi/15)*18; Delta_T_t(i)=q_t*R;

T_s=T_inf+Delta_T_t(i); T_s_time(i)=T_s-273.15;

end

figure; plot(t,Delta_T_t+25,'LineWidth',2);

xlabel('Time (Hr)');
ylabel(' Room Temperature Variation (K)');

title('High-class iron (day: 10 m/s, Tsky=285 K)');

grid on;
Case: 05 (Summer Night: Windspeed=0 m/s, Tsky=235 K)
% Parameters

k1=1; k2=2; A1=27; A2=18; A3=45; m=27000; Cp=800; Tsky=235; sigma=5.67e-8;

delta_x1=0.3; delta_x2=0.005;

Troom=298; total_time=32400; time_step=1; num_steps=total_time/time_step;

time=linspace(0,total_time,num_steps); Troom_history=zeros(1,num_steps);
Ts_history=zeros(1,num_steps);

for t=1:num_steps

Ts=(sigma*A3*(Troom^4-Tsky^4))/((k1*A1/delta_x1)+(k2*A2/delta_x2))+Troom;

Q=sigma*A3*(Troom^4-Tsky^4);
Troom=Troom-(Q*time_step)/(m*Cp);

Troom_history(t)=Troom; Ts_history(t)=Ts;

end

figure;

plot(time,Troom_history,'b','LineWidth',1.5);

hold on;

plot(time,Ts_history,'r','LineWidth',1.5);
hold off;
grid on;

xlabel('Time (s)');

ylabel('Temperature (K)');

title('Room and Surface Temperature vs Time');

legend('T_{room}','T_s');

Case: 06 (Summer Night: Windspeed=0 m/s, Tsky=285 K)

% Parameters

k1 = 1; k2 = 2; A1 = 27; A2 = 18; A3 = 45; m = 27000; Cp = 800; Tsky = 285;

sigma = 5.67e-8; delta_x1 = 0.3; delta_x2 = 0.005;

Troom = 298; total_time = 32400; time_step = 1; num_steps = total_time / time_step;

time = linspace(0, total_time, num_steps); Troom_history = zeros(1, num_steps); Ts_history =

zeros(1, num_steps);

% Time-stepping loop

for t = 1:num_steps
Ts = (sigma * A3 * (Troom^4 - Tsky^4)) / ((k1 * A1 / delta_x1) + (k2 * A2 / delta_x2)) +
Troom;

Q = sigma * A3 * (Troom^4 - Tsky^4);

delta_T = -(Q * time_step) / (m * Cp);

Troom = Troom + delta_T;

 Troom_history(t) = Troom; Ts_history(t) = Ts;

end

% Plot results

figure; plot(time, Troom_history, 'b', 'LineWidth', 1.5); hold on;

plot(time, Ts_history, 'r', 'LineWidth', 1.5); hold off;

grid on; xlabel('Time (s)'); ylabel('Temperature (K)');

title('Room and Surface Temperature vs Time'); legend('T_{room}', 'T_s');

Case: 07 (Summer Night: Windspeed=10 m/s, Tsky=235 K)

% Parameters

k1 = 1; k2 = 2; A1 = 27; A2 = 18; A3 = 45; m = 27000; Cp = 800; Tsky = 235;

sigma = 5.67e-8; delta_x1 = 0.3; delta_x2 = 0.005;

Troom = 298; total_time = 32400; time_step = 1; num_steps = total_time / time_step;

time = linspace(0, total_time, num_steps); Troom_history = zeros(1, num_steps); Ts_history =

zeros(1, num_steps);

% Time-stepping loop

for t = 1:num_steps

Ts = (sigma * A3 * (Troom^4 - Tsky^4)) / ((k1 * A1 / delta_x1) + (k2 * A2 / delta_x2)) +

Troom;

Q = sigma * A3 * (Troom^4 - Tsky^4);

 delta_T = -(Q * time_step) / (m * Cp);

Troom = Troom + delta_T;

Troom_history(t) = Troom; Ts_history(t) = Ts;

end
% Plot results

figure; plot(time, Troom_history, 'b', 'LineWidth', 1.5); hold on;

plot(time, Ts_history, 'r', 'LineWidth', 1.5); hold off;

grid on; xlabel('Time (s)'); ylabel('Temperature (K)');

title('Room and Surface Temperature vs Time'); legend('T_{room}', 'T_s');

Case: 08 (Summer Night: Windspeed=10 m/s, Tsky=285 K)

% Parameters

k1 = 1; k2 = 2; A1 = 27; A2 = 18; A3 = 45; m = 27000; Cp = 800; Tsky = 285;

sigma = 5.67e-8; delta_x1 = 0.3; delta_x2 = 0.005;

Troom = 298; total_time = 32400; time_step = 1; num_steps = total_time / time_step;

time = linspace(0, total_time, num_steps); Troom_history = zeros(1, num_steps); Ts_history =

zeros(1, num_steps);

% Time-stepping loop

for t = 1:num_steps
 Ts = (sigma * A3 * (Troom^4 - Tsky^4)) / ((k1 * A1 / delta_x1) + (k2 * A2 / delta_x2)) +
Troom;

Q = sigma * A3 * (Troom^4 - Tsky^4);

delta_T = -(Q * time_step) / (m * Cp);

Troom = Troom + delta_T;

Troom_history(t) = Troom; Ts_history(t) = Ts;

end
% Plot results

figure; plot(time, Troom_history, 'b', 'LineWidth', 1.5); hold on;

plot(time, Ts_history, 'r', 'LineWidth', 1.5); hold off;

grid on; xlabel('Time (s)'); ylabel('Temperature (K)');

title('Room and Surface Temperature vs Time'); legend('T_{room}', 'T_s');

Case: 09 (Winter Daytime: Windspeed=0 m/s, Tsky=235 K)

Low-Iron glass
T_inf = 25 + 273.15;

T_s = 80 + 273.15;

L = 3;

T_avg = (T_inf + T_s) / 2;

k = 0.026;

nu = 15.89e-6;

alpha = 22.14e-6;

beta = 1 / T_avg;
g = 9.81;

Pr = nu / alpha;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

delta_x = 0.3;

k_material = 1;
A = 36;

R1 = (delta_x / (k_material * A)) + (1 / (h * A));

R2 = (0.005 / (2 * 18)) + (1 / (h * 18));

R = (1 / R1 + 1 / R2) ^ -1;

I0 = 800 * 0.85;

t = linspace(0, 9, 100);

Delta_T_t = zeros(size(t));
T_s_time = zeros(size(t));

for i = 1:length(t)

T_avg = (T_inf + T_s) / 2;

beta = 1 / T_avg;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;
q_t = I0 * sin(t(i) * pi / 9) * 18;

Delta_T_t(i) = q_t * R;

T_s = T_inf + Delta_T_t(i);

T_s_time(i) = T_s - 273.15;

end

figure;

plot(t, Delta_T_t + 25, 'LineWidth', 2);

xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('Room Temperature Variation vs Time');

grid on;

High-Iron glass

T_inf = 25 + 273.15;
T_s = 80 + 273.15;

L = 3;

T_avg = (T_inf + T_s) / 2;

k = 0.026;

nu = 15.89e-6;
alpha = 22.14e-6;

beta = 1 / T_avg;

g = 9.81;

Pr = nu / alpha;
Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

delta_x = 0.3;

k_material = 1;

A = 36;

R1 = (delta_x / (k_material * A)) + (1 / (h * A));

R2 = (0.005 / (2 * 18)) + (1 / (h * 18));

R = (1 / R1 + 1 / R2) ^ -1;

I0 = 800 * 0.7;

t = linspace(0, 9, 100);

Delta_T_t = zeros(size(t));

T_s_time = zeros(size(t));

for i = 1:length(t)
T_avg = (T_inf + T_s) / 2;

beta = 1 / T_avg;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

q_t = I0 * sin(t(i) * pi / 9) * 18;

Delta_T_t(i) = q_t * R;
T_s = T_inf + Delta_T_t(i);

T_s_time(i) = T_s - 273.15;

end

figure;
plot(t, Delta_T_t + 25, 'LineWidth', 2);

xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('High-Iron Glass (day: v=0 m/s and Tsky=235 K)');

grid on;

Case: 10 (Winter Daytime: Windspeed=0 m/s, Tsky=285 K)

Low-Iron glass
T_inf = 25 + 273.15;

T_s = 80 + 273.15;
L = 3;

T_avg = (T_inf + T_s) / 2;

k = 0.026;

nu = 15.89e-6;

alpha = 22.14e-6;
beta = 1 / T_avg;

g = 9.81;

Pr = nu / alpha;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

delta_x = 0.3;

k_material = 1;

A = 36;

R1 = (delta_x / (k_material * A)) + (1 / (h * A));

R2 = (0.005 / (2 * 18)) + (1 / (h * 18));

R = (1 / R1 + 1 / R2) ^ -1;

I0 = 800 * 0.85;

t = linspace(0, 9, 100);

Delta_T_t = zeros(size(t));

T_s_time = zeros(size(t));

for i = 1:length(t)

T_avg = (T_inf + T_s) / 2;

beta = 1 / T_avg;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

q_t = I0 * sin(t(i) * pi / 9) * 18;

Delta_T_t(i) = q_t * R;
T_s = T_inf + Delta_T_t(i);
T_s_time(i) = T_s - 273.15;

end

figure;

plot(t, Delta_T_t + 25, 'LineWidth', 2);

xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('Low-Iron glass (day: v=0 m/s, Tsky=285 K)');

grid on;

High-Iron glass

T_inf = 25 + 273.15;

T_s = 80 + 273.15;

L = 3;

T_avg = (T_inf + T_s) / 2;

k = 0.026;
nu = 15.89e-6;

alpha = 22.14e-6;

beta = 1 / T_avg;
g = 9.81;

Pr = nu / alpha;

Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;
Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

delta_x = 0.3;

k_material = 1;

A = 36;

R1 = (delta_x / (k_material * A)) + (1 / (h * A));

R2 = (0.005 / (2 * 18)) + (1 / (h * 18));

R = (1 / R1 + 1 / R2) ^ -1;
I0 = 800 * 0.7;

t = linspace(0, 9, 100);

Delta_T_t = zeros(size(t));

T_s_time = zeros(size(t));

for i = 1:length(t)

T_avg = (T_inf + T_s) / 2;

beta = 1 / T_avg;
Gr = (g * beta * (T_s - T_inf) * L^3) / (nu^2);

Ra = Gr * Pr;

Nu = 0.68 + (0.67 * Ra^(1/4)) / ((1 + (0.492 / Pr)^(9/16))^(4/9));

h = (Nu * k) / L;

q_t = I0 * sin(t(i) * pi / 9) * 18;

Delta_T_t(i) = q_t * R;

T_s = T_inf + Delta_T_t(i);

T_s_time(i) = T_s - 273.15;
end

figure;

plot(t, Delta_T_t + 25, 'LineWidth', 2);

xlabel('Time (hr)');
ylabel('Room Temperature Variation (K)');

title('High-Iron glass (day: v=0 m/s, Tsky=285 K)');

grid on;

Case: 11 (Winter Daytime: Windspeed=10 m/s, Tsky=235 K)

Low-Iron glass
T_inf=25+273.15;

T_s=80+273.15;

L=3;

T_avg=(T_inf+T_s)/2;

k=0.026;

nu=15.89e-6;

alpha=22.14e-6;
beta=1/T_avg;
g=9.81;

Pr=nu/alpha;

U_inf=10;

Re_L=U_inf*L/nu;
if Re_L<5e5

Nu_avg=0.664*Re_L^0.5*Pr^(1/3);

else

Nu_avg=0.037*Re_L^(4/5)*Pr^(1/3);

end

h=Nu_avg*k/L;

delta_x=0.3;

k_material=1;
A=36;

R1=(delta_x/(k_material*A))+(1/(h*A));

R2=(0.005/(2*18))+(1/(h*18));

R=(1/R1+1/R2)^-1;

I0=800*0.85;

t=linspace(0,9,100);

Delta_T_t=zeros(size(t));
T_s_time=zeros(size(t));

for i=1:length(t)

T_avg=(T_inf+T_s)/2;

beta=1/T_avg;

Gr=(g*beta*(T_s-T_inf)*L^3)/(nu^2);

Ra=Gr*Pr;

Nu=0.68+(0.67*Ra^(1/4))/((1+(0.492/Pr)^(9/16))^(4/9));
h=(Nu*k)/L;
q_t=I0*sin(t(i)*pi/9)*18;

Delta_T_t(i)=q_t*R;

T_s=T_inf+Delta_T_t(i);

T_s_time(i)=T_s-273.15;
end

figure;

plot(t,Delta_T_t +25,'LineWidth',2);

xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('Low-class iron (day: 10 m/s, Tsky=235 K)');

grid on;

High-Iron glass

T_inf=25+273.15;

T_s=80+273.15;
L=3;

T_avg=(T_inf+T_s)/2;

k=0.026;

nu=15.89e-6;

alpha=22.14e-6;
beta=1/T_avg;

g=9.81;

Pr=nu/alpha;

U_inf=10;
Re_L=U_inf*L/nu;

if Re_L<5e5

Nu_avg=0.664*Re_L^0.5*Pr^(1/3);

else

Nu_avg=0.037*Re_L^(4/5)*Pr^(1/3);

end

h=Nu_avg*k/L;

delta_x=0.3;
k_material=1;

A=36;

R1=(delta_x/(k_material*A))+(1/(h*A));

R2=(0.005/(2*18))+(1/(h*18));

R=(1/R1+1/R2)^-1;

I0=800*0.7;

t=linspace(0,9,100);
Delta_T_t=zeros(size(t));

T_s_time=zeros(size(t));

for i=1:length(t)

T_avg=(T_inf+T_s)/2;

beta=1/T_avg;

Gr=(g*beta*(T_s-T_inf)*L^3)/(nu^2);

Ra=Gr*Pr;
Nu=0.68+(0.67*Ra^(1/4))/((1+(0.492/Pr)^(9/16))^(4/9));
h=(Nu*k)/L;

q_t=I0*sin(t(i)*pi/9)*18;

Delta_T_t(i)=q_t*R;

T_s=T_inf+Delta_T_t(i);
T_s_time(i)=T_s-273.15;

end

figure;

plot(t,Delta_T_t +25,'LineWidth',2);

xlabel('Time (hr)');

ylabel('Room Temperature Variation (K)');

title('High-class iron (day: 10 m/s, Tsky=235 K)');

grid on;

Case: 12 (Winter Daytime: Windspeed=10 m/s, Tsky=285 K)

Low-Iron glass
T_inf=25+273.15; T_s=80+273.15; L=3; T_avg=(T_inf+T_s)/2; k=0.026; nu=15.89e-6;
alpha=22.14e-6; beta=1/T_avg; g=9.81; Pr=nu/alpha; U_inf=10; Re_L=U_inf*L/nu;

if Re_L<5e5
 Nu_avg=0.664*Re_L^0.5*Pr^(1/3);

else

Nu_avg=0.037*Re_L^(4/5)*Pr^(1/3);

end
h=Nu_avg*k/L;

delta_x=0.3; k_material=1; A=36;

R1=(delta_x/(k_material*A))+(1/(h*A)); R2=(0.005/(2*18))+(1/(h*18)); R=(1/R1+1/R2)^-1;

I0=800*0.85; t=linspace(0,9,100); Delta_T_t=zeros(size(t)); T_s_time=zeros(size(t));

for i=1:length(t)

T_avg=(T_inf+T_s)/2; beta=1/T_avg; Gr=(g*beta*(T_s-T_inf)*L^3)/(nu^2); Ra=Gr*Pr;

Nu=0.68+(0.67*Ra^(1/4))/((1+(0.492/Pr)^(9/16))^(4/9)); h=(Nu*k)/L;

q_t=I0*sin(t(i)*pi/9)*18; Delta_T_t(i)=q_t*R;
T_s=T_inf+Delta_T_t(i); T_s_time(i)=T_s-273.15;

end

figure; plot(t,Delta_T_t+25,'LineWidth',2);

xlabel('Time (Hr)');

ylabel(' Room Temperature Variation (K)');

title('Low-class iron (day: 10 m/s, Tsky=285 K)');

grid on;
High-Iron glass

T_inf=25+273.15; T_s=80+273.15; L=3; T_avg=(T_inf+T_s)/2; k=0.026; nu=15.89e-6;

alpha=22.14e-6; beta=1/T_avg; g=9.81; Pr=nu/alpha; U_inf=10; Re_L=U_inf*L/nu;

if Re_L<5e5

Nu_avg=0.664*Re_L^0.5*Pr^(1/3);

else
Nu_avg=0.037*Re_L^(4/5)*Pr^(1/3);

end

h=Nu_avg*k/L;

delta_x=0.3; k_material=1; A=36;

R1=(delta_x/(k_material*A))+(1/(h*A)); R2=(0.005/(2*18))+(1/(h*18)); R=(1/R1+1/R2)^-1;

I0=800*0.7; t=linspace(0,9,100); Delta_T_t=zeros(size(t)); T_s_time=zeros(size(t));

for i=1:length(t)
T_avg=(T_inf+T_s)/2; beta=1/T_avg; Gr=(g*beta*(T_s-T_inf)*L^3)/(nu^2); Ra=Gr*Pr;

Nu=0.68+(0.67*Ra^(1/4))/((1+(0.492/Pr)^(9/16))^(4/9)); h=(Nu*k)/L;

q_t=I0*sin(t(i)*pi/9)*18; Delta_T_t(i)=q_t*R;

T_s=T_inf+Delta_T_t(i); T_s_time(i)=T_s-273.15;

end
figure; plot(t,Delta_T_t+25,'LineWidth',2);

xlabel('Time (Hr)');

ylabel(' Room Temperature Variation (K)');

title('High-class iron (day: 10 m/s, Tsky=285 K)');

grid on;

Case: 13 (Winter Night: Windspeed=0 m/s, Tsky=235 K)

% Parameters

k1=1; k2=2; A1=27; A2=18; A3=45; m=27000; Cp=800; Tsky=235; sigma=5.67e-8;

delta_x1=0.3; delta_x2=0.005;

Troom=298; total_time= 54000; time_step=1; num_steps=total_time/time_step;

time=linspace(0,total_time,num_steps); Troom_history=zeros(1,num_steps);
Ts_history=zeros(1,num_steps);

for t=1:num_steps

Ts=(sigma*A3*(Troom^4-Tsky^4))/((k1*A1/delta_x1)+(k2*A2/delta_x2))+Troom;

Q=sigma*A3*(Troom^4-Tsky^4);

Troom=Troom-(Q*time_step)/(m*Cp);

Troom_history(t)=Troom; Ts_history(t)=Ts;
end
figure;

plot(time,Troom_history,'b','LineWidth',1.5);

hold on;

plot(time,Ts_history,'r','LineWidth',1.5);
hold off;

grid on;

xlabel('Time (s)');

ylabel('Temperature (K)');

title('Room and Surface Temperature vs Time');

legend('T_{room}','T_s');

Case: 14 (Winter Night: Windspeed=0 m/s, Tsky=285 K)

% Parameters

k1 = 1; k2 = 2; A1 = 27; A2 = 18; A3 = 45; m = 27000; Cp = 800; Tsky = 285;

sigma = 5.67e-8; delta_x1 = 0.3; delta_x2 = 0.005;

Troom = 298; total_time = 54000; time_step = 1; num_steps = total_time / time_step;

time = linspace(0, total_time, num_steps); Troom_history = zeros(1, num_steps); Ts_history =

zeros(1, num_steps);

% Time-stepping loop
for t = 1:num_steps
 Ts = (sigma * A3 * (Troom^4 - Tsky^4)) / ((k1 * A1 / delta_x1) + (k2 * A2 / delta_x2)) +
Troom;

Q = sigma * A3 * (Troom^4 - Tsky^4);

delta_T = -(Q * time_step) / (m * Cp);

Troom = Troom + delta_T;

Troom_history(t) = Troom; Ts_history(t) = Ts;

end
% Plot results

figure; plot(time, Troom_history, 'b', 'LineWidth', 1.5); hold on;

plot(time, Ts_history, 'r', 'LineWidth', 1.5); hold off;

grid on; xlabel('Time (s)'); ylabel('Temperature (K)');

title('Room and Surface Temperature vs Time'); legend('T_{room}', 'T_s');

Case: 15 (Winter Night: Windspeed=10 m/s, Tsky=235 K)

% Parameters

k1 = 1; k2 = 2; A1 = 27; A2 = 18; A3 = 45; m = 27000; Cp = 800; Tsky = 235;

sigma = 5.67e-8; delta_x1 = 0.3; delta_x2 = 0.005;

Troom = 298; total_time = 54000; time_step = 1; num_steps = total_time / time_step;

time = linspace(0, total_time, num_steps); Troom_history = zeros(1, num_steps); Ts_history =
zeros(1, num_steps);

% Time-stepping loop

for t = 1:num_steps

Ts = (sigma * A3 * (Troom^4 - Tsky^4)) / ((k1 * A1 / delta_x1) + (k2 * A2 / delta_x2)) +

Troom;

Q = sigma * A3 * (Troom^4 - Tsky^4);

delta_T = -(Q * time_step) / (m * Cp);

Troom = Troom + delta_T;

Troom_history(t) = Troom; Ts_history(t) = Ts;

end

% Plot results

figure; plot(time, Troom_history, 'b', 'LineWidth', 1.5); hold on;

plot(time, Ts_history, 'r', 'LineWidth', 1.5); hold off;

grid on; xlabel('Time (s)'); ylabel('Temperature (K)');

title('Room and Surface Temperature vs Time'); legend('T_{room}', 'T_s');

Case: 16 (Winter Night: Windspeed=10 m/s, Tsky=285 K)

% Parameters
k1 = 1; k2 = 2; A1 = 27; A2 = 18; A3 = 45; m = 27000; Cp = 800; Tsky = 285;

sigma = 5.67e-8; delta_x1 = 0.3; delta_x2 = 0.005;

Troom = 298; total_time = 54000; time_step = 1; num_steps = total_time / time_step;

time = linspace(0, total_time, num_steps); Troom_history = zeros(1, num_steps); Ts_history =

zeros(1, num_steps);

% Time-stepping loop

for t = 1:num_steps
Ts = (sigma * A3 * (Troom^4 - Tsky^4)) / ((k1 * A1 / delta_x1) + (k2 * A2 / delta_x2)) +
Troom;
Q = sigma * A3 * (Troom^4 - Tsky^4);

delta_T = -(Q * time_step) / (m * Cp);

Troom = Troom + delta_T;

Troom_history(t) = Troom; Ts_history(t) = Ts;

end

% Plot results

figure; plot(time, Troom_history, 'b', 'LineWidth', 1.5); hold on;

plot(time, Ts_history, 'r', 'LineWidth', 1.5); hold off;

grid on; xlabel('Time (s)'); ylabel('Temperature (K)');

title('Room and Surface Temperature vs Time'); legend('T_{room}', 'T_s');