Heat Flow Equation Matlab

0
votes

I'm currently working on a project to build a solar food dryer and I need to model on Matlab how temperature of the product will change with respect to change in solar radiation, Q .

Q is given by ;

Q =960*(sin(2*pi*Time2/24)).^2; %W/m2

where 

Time2 = (1:t:12); %hours

The heat flow equation is given by

Q(t)A = mcp*(T2-T1) + (mw*lw)

where :

mw = 0.706; % Mass of water loss from product in hr (Kg/h)

m = 15; % Mass of product to dry (Kg)

lw = 2260000; % Latent heat of vaporisation of water (J/Kg)

A = 1; % Surface Area of collector (m^2)

cp= 3746; % Specific heat capacity of product (J/Kg/C)

T1 = temperature at t

T2 = temperature at t + dt

Manipulating the heat flow give T2 as ;

T2= (((Q*A*3600) -(mw*lw))/(m*cp)) + T1; % 3600 is there to convert j/s to J/h

however implementing this on Matlab is proving a challenge for me- I'm fairly new to Matlab

This is what I have so far ;

close all

clear;

mw = 0.706; % Mass of water loss from product in hr (Kg/h)

m = 15; % Mass of product to dry (Kg)

lw = 2260000; % Latent heat of vaporisation of water (J/Kg)

A = 1; % Surface Area of Collector (m^2)

cp= 3746; % Specific heat capacity of product (J/Kg/C)

t = 1; % Time step

T = 24; % Initial Temperature (°C)



Time2=(1:t:12); hours

Q=960*(sin(2*pi*Time2/24)).^2; % Solar irradiation in tropical regions at specific time (W/m2)

for j = 1:12

 T(j+1)= ((((QQ2(j)*A*j*3600))-(mw*lw))/(m*cp))+ T(1);

end

figure(2)

plot(T)

title('Temperature')

xlabel('time (hours)')

ylabel('Temperature (°C)')

This seems wrong since the mass, m should decrease by mw after each hr and the temperature profile should follow the profile of the solar radiation. i.e peak at the same time

I have been spending days to get my head around this but i'm pretty bad at Matlab so I haven't made any meaningful progress. Any help will be appreciated

1
matlab

1 Answers

0
votes

So I am not sure if this is what you are looking for Ran, but there are a few points that looked like typos (?) and these would change behaviour. You have QQ2(j) suddenly appears in the middle of your script... I presumed that was just Q[j). Each cycle of your loop you add to t(1) and I think you meant t(j)? And surely the loop should also decrease the m? So I modified to this...

for j = 1:12
     T(j+1)= ((((Q(j)*A*j*3600))-(mw*lw))/(m*cp))+  T(j);
     m=m-mw*t;
end

Now T peaks at 12noon.

All that said, I would think a problem like this would be better solved with a differential equation solver like ‘ode45()’, but I would need to see that differential equation before advising how to do that in Matlab, but it shouldn’t be too tricky!

 clear;
global conv;
conv=60*60;
mw = 0.706/conv; % Mass of water loss from product in hr (Kg/h)
m = 15.0; % Mass of product to dry (Kg)
lw = 2260000.0; % Latent heat of vaporisation of water (J/Kg)
A= 1.0; % Surface Area of Collector (m^2)
cp= 3746.0; % Specific heat capacity of product (J/Kg/C)
T = 24.0; % Initial Temperature (°C)
%%
%% for j = 1:12
%% T(j+1)= ((((Q(j)*A*j*3600))-(mw*lw))/(m*cp))+ T(j);
%% m=m-mw*t;
%% end

tspan = [0, 12*conv];
yo=[T;m];
[t,y] = ode45(@(t,y) ode(t, y,[A,mw,lw,cp]), tspan, yo);

yfinal=y;
figure (1)
plot(t./conv,yfinal(:,1))
title('Temperature')
xlabel('time (hours)')
ylabel('Temperature (°C)')
figure (2)
tqs=linspace(0,12);
Qt=960.0*(sin(2.0.*pi.*tqs./(24.0))).^2;
plot(tqs,Qt);
function dydt=ode(t,y,x)
    global conv;
    A =x(1);
    mw =x(2);
    lw =x(3);
    cp =x(4);
    m=y(2);
    Q=960*(sin(2*pi*t/(24*conv)))^2; % Solar irradiation in tropical regions at specific time (W/m2)
    dydt=[((Q * A)-(mw*lw))/(m*cp);-mw];
end

This gives this output, I think the equations will need a fiddle, but hopefully the structure of the ODE helps? P.s., I am not convinced that peak temperature would match peak heat input since usually there is a lag. Hottest part of the day is often 3pm in the UK... but I haven’t looked at heat transfer since High School... THIS

Regards R