Radiation heat transfer using PDE toolbox

Question

Ruben 2020-8-21

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/582698-radiation-heat-transfer-using-pde-toolbox

评论： Paul Schmitz 2023-9-8

Hello,

I'm trying to obtain the temperature at a given point of my layered 2D cross-section geometry considering thermal conductivity, convection and radiation. Inside the geometry there are horizontal and vertical edges.

The code is working fine when no convection or radiation are used. However when I include them, the results don't vary much (only on the 6th decimal) although they should do it (for radiation) since I'm considering high temperatures.

I'm defining what edges should consider the radiation and convection in the following way:

eps = @(region,state) obj.epsilon*region

outerCC = @(region,state) obj.conv_coef*region

thermalBC(obj.tm_SS,'Edge', 31,33],'ConvectionCoefficient',outerCC,'Emissivity',eps,'AmbientTemperature',Ta,'Vectorized','on');

Do I need to specify the pde coefficients for each edge or face? or the solver extracts them when using thermalBC?

In the code, I'm not using commands such as ApplyBoundaryCondition or specifyCoefficients . If I try to use specifyCoefficients, I define the coefficients in the following way:

Ta = 25;

c = length*K;

a = @(~,state) 2*obj.conv_coef + 2*obj.epsilon*obj.tm_SS.StefanBoltzmannConstant*state.u.^3;

f = 2*obj.conv_coef*Ta + 2*obj.epsilon*obj.tm_SS.StefanBoltzmannConstant*Ta^4;

d = length*Density*SpecificHeat;

coefs= specifyCoefficients(obj.tm_SS,'m',0,'d',d,'c',c,'a',a,'f',f,'Face',2);

but always obtain the following error:

Check for missing argument or incorrect argument data type in call to function 'specifyCoefficients'.

The model is defined as follows:

obj.tm_SS = createpde('thermal','steadystate');

gdm = [r0;r9;r1;r2;r3;r4;r5;r6;r7;r8]';

obj.g = decsg(gdm,'R0+(R9-R1)+R1+(R2-R3-R6-R7)+R3+R6+R7+(R4-R6-R7)+R6+R7+(R5-R8)+R8',['R0';'R9';'R1';'R2';'R3';'R4';'R5';'R6';'R7';'R8']');

obj.geometry=geometryFromEdges(obj.tm_SS,obj.g);

internalHeatSource(obj.tm_SS,obj.hfV_h,'Face',7);

Do you know what I'm doing wrong? or, Do you have any advice on how to procede?

Thank you.

2 个评论
显示无隐藏无

Ravi Kumar 2020-8-21

Hi Ruben,

Can you attache the complete working script to try out. Also which version of MATLAB are you using?

Regards,

Ravi

Ruben 2020-8-24

Dear Ravi,

thank you. I'm using Matlab R2020a. Find below the script:

clear

model = createpde('thermal','steadystate');

hf_h = 400e6;

RadandConv=1;

% geometry dimensions

width=1000e-6;

w_h=7e-6;

h_h=0.1e-6;

length_h=1e-3;

h_a= 10e-6;

h_i= 4.8e-6;

htr= h_i;

tro= 2.5e-6;

hw= 500e-6;

wt= 2e-6;

wtSi= 13e-6;

hSi= 0.5e-6;

wo= 0;

H_max= 10e-6;

H_min= 0.01e-6;

H_grad= 1.3;

hfV = hf_h / h_h; %HEAT FLUX W/m

%Material properties

ka=0.031; Da=1.225; Qa=1.006;

kh=310; Dh=19320; Qh=130;

ki=1.5; Di=2650; Qi=710;

kw=130; Dw=2500; Qw=710;

%Cross section definition

r0 = [3 4 -width/2 width/2 width/2 -width/2 0 0 h_a h_a];

r1 = [3 4 -w_h/2 w_h/2 w_h/2 -w_h/2 0 0 h_h h_h];

r2 = [3 4 -width/2 width/2 width/2 -width/2 -h_i -h_i 0 0];

r5= [3 4 -width/2 width/2 width/2 -width/2 -h_i-hw -h_i-hw -h_i -h_i];

r6 = [3 4 -(w_h/2+tro+wo) -(w_h/2+tro+wt-wo) -(w_h/2+tro+wt) -(w_h/2+tro) -htr -htr 0 0];

r7 = [3 4 (w_h/2+tro+wo) (w_h/2+tro+wt-wo) (w_h/2+tro+wt) (w_h/2+tro) -htr -htr 0 0];

r8= [3 4 -wtSi/2 wtSi/2 wtSi/2 -wtSi/2 -h_i -h_i -h_i-hSi -h_i-hSi];

gdm = [r0;r1;r2;r5;r6;r7;r8]';

g = decsg(gdm,'(R0-R1)+R1+(R2-R6-R7)+R6+R7+(R5-R8)+R8',['R0';'R1';'R2';'R5';'R6';'R7';'R8']');

geometry=geometryFromEdges(model,g);

e_bot = 5; %bottom edge

f_h = 9; %face of the internal heat source

internalHeatSource(model,hfV,'Face',f_h);

epsilon_i=0.79;

epsilon_Si= 0.67;

epsilon_h = 0.8;

conv_coef= 1;

ta = 25;

c = length_h*ki;

a = @(~,state) 2*conv_coef + 2*epsilon_i*model.StefanBoltzmannConstant*state.u.^3;

f = 2*conv_coef*ta + 2*epsilon_i*model.StefanBoltzmannConstant*ta^4;

d = length_h*Di*Qi;

thermalBC(model,'Edge',e_bot,'Temperature',ta);

if RadandConv==1

model.StefanBoltzmannConstant = 5.670367e-8;

eps_i = @(region,state) epsilon_i*region;

outerCC = @(region,state) conv_coef*region;

thermalBC(model,'Edge',[6,7,8,9,21,23,29],'ConvectionCoefficient',outerCC,'Emissivity',eps_i,'AmbientTemperature',ta,'Vectorized','on');

% coefs= specifyCoefficients(model,'m',1,'d',d,'c',c,'a',a,'f',f,'Edge',[6,7,8,9,21,23,29]);

% coefs= specifyCoefficients(model,'m',1,'d',d,'c',c,'a',a,'f',f,'Face',4);

eps_h = @(region,state) epsilon_h*region;

thermalBC(model,'Edge',[2,3,4],'ConvectionCoefficient',outerCC,'Emissivity',eps_h,'AmbientTemperature',ta,'Vectorized','on');

eps_Si = @(region,state) epsilon_Si*region;

thermalBC(model,'Edge',[10,11,12,27,31],'ConvectionCoefficient',outerCC,'Emissivity',eps_Si,'AmbientTemperature',ta,'Vectorized','on');

end

thermalProperties(model,'Face',[1 2 5 8],'ThermalConductivity',ka,'MassDensity',Da,'SpecificHeat',Qa);

thermalProperties(model,'Face',9,'ThermalConductivity',kh,'MassDensity',Dh,'SpecificHeat',Qh);

thermalProperties(model,'Face',[3 4 7],'ThermalConductivity',ki,'MassDensity',Di,'SpecificHeat',Qi);

thermalProperties(model,'Face',6,'ThermalConductivity',kw,'MassDensity',Dw,'SpecificHeat',Qw);

X_LIM = [-25e-6 25e-6]; Y_LIM= [-10e-6 10e-6];FIG_AXIS=[X_LIM Y_LIM];

nFIG=1; figH = figure(nFIG);

pdegplot(model,'EdgeLabels','on');

axis(FIG_AXIS);hold off

nFIG = nFIG + 1; figH = figure(nFIG);

pdegplot(model,'FaceLabels','on');

axis(FIG_AXIS);hold off

msh=generateMesh(model,'Hmax',H_max,'Hmin',H_min,'Hgrad',H_grad);

nFIG = nFIG + 1;

figH = figure(nFIG); set(figH,'NumberTitle','off','Name','2 - Geometry with Mesh');

pdeplot(model); axis(FIG_AXIS);hold off

R = solve(model);

Temp=R.Temperature;

% Temperature measurement

x=[-(w_h/2+tro):0.01e-6:(w_h/2+tro)];

y=-h_i/2*ones(size(x));

T_interp = interpolateTemperature(R,x,y);

T = mean(T_interp);

disp(['measured Temp = ' num2str(T) ' C']);

nFIG = nFIG + 1;

figH = figure(nFIG); set(figH,'NumberTitle','off','Name','5 - Temperature, Steady State Solution');

pdeplot(model,'XYData',Temp,'Contour','on','ColorMap','hot'); hold on; xlabel('m'); ylabel('m'); zlabel('Celsius');

pdegplot(model,'FaceLabels','on'); axis(FIG_AXIS);

title(['Tmeas=' num2str(T) ' ºC']);

hold off

Regards,

Ruben

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Ravi Kumar 2020-8-24

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/582698-radiation-heat-transfer-using-pde-toolbox#answer_484412

Hi Ruben,

Thanks for the code, it helped. All the BCs inside RadandConv == 1 block are applied on the edges that are interior to the domain. These boundary conditions, more accurately interface conditions, are not supported in PDE Toolbox:

https://www.mathworks.com/help/pde/ug/do-not-specify-boundary-conditions-between-subdomains.html

Based on your setup it looks like you need surface-to-surface radiation modelling capability, unfortunately, it is not supported in PDE Toolbox.

Regards,

Ravi