FEM to solve internal heat source generated

21 views (last 30 days)
Eddy
Eddy on 8 Sep 2023
Commented: Bradley on 6 Oct 2023
Hi,
I am learning to use the PDE toolbox to solve a 2D heat transfer problem. The model is adapted from an example in Matlab. It is a rectangular surface with a another circular surface on it. The circular surface generates heat at the rate of 500 W/m3 and has a thermal conductivity of 1/6. The rectangular surface has a thermal conductivity of 1, three boundaries at 25 degree C and one boundary with a convection coefficient of 7.371 and ambient temperature of 30 degree C. The model works as expected and I am attaching the code below.
Now, I would like to find, using FEM, the heat generated by the circular plate for it to reach 100 deg.C. How can this be acheived?
Any help would be appreciated. Thanks
model = createpde;
R1 = [3,4,0,10,10,0,10,10,0,0]';
C1 = [1,5,8,0.2]';
% C1 = [1,0.5,-0.25,0.25]';
C1 = [C1;zeros(length(R1) - length(C1),1)];
gm = [R1,C1];
sf = 'R1+C1';
% Create the geometry.
ns = char('R1','C1');
ns = ns';
g = decsg(gm,sf,ns);
% Include the geometry in the model and plot it.
thermalModel = createpde("thermal","steadystate");
geometryFromEdges(thermalModel,g);
figure(1);
pdegplot(thermalModel,'FaceLabels','on','EdgeLabels','on')
axis equal
xlim([0,10])
mesh = generateMesh(thermalModel,"Hmax",0.25,"Hedge",{[5 6 7 8],0.01});
figure(2)
pdeplot(thermalModel)
thermalProperties(thermalModel,'ThermalConductivity',1,'Face',1);
thermalProperties(thermalModel,'ThermalConductivity',1/6,'Face',2);
internalHeatSource(thermalModel,500,'Face',2);
thermalBC(thermalModel,'Edge',[1,2,3],'Temperature',25);
thermalBC(thermalModel,'Edge',4,'ConvectionCoefficient',7.371,'AmbientTemperature',30);
thermalresults = solve(thermalModel);
T = thermalresults.Temperature;
set(gcf,'color','white')
figure (3)
pdeplot(thermalModel,'XYData',T,'Contour','on')
axis equal
title(['\fontsize{16}Steady-State Temperature. Max temp. ',num2str(max(T)),' \fontsize{16}^0C'])

Answers (1)

Shubham
Shubham on 20 Sep 2023
Hey @Eddy,
I understand that you want to find the heat generated by the circular plate for it to reach 100°C using FEM (Finite Elements Method).
You can define the target temperature and modify the heat generation rate until the current temperature reaches the target temperature. You can use a tolerance value for getting the desired temperature and adjust accordingly to control the accuracy of the result.
I have modified the code as below:
model = createpde;
R1 = [3,4,0,10,10,0,10,10,0,0]';
C1 = [1,5,8,0.2]';
C1 = [C1;zeros(length(R1) - length(C1),1)];
gm = [R1,C1];
sf = 'R1+C1';
% Create the geometry.
ns = char('R1','C1');
ns = ns';
g = decsg(gm,sf,ns);
% Include the geometry in the model and plot it.
thermalModel = createpde("thermal","steadystate");
geometryFromEdges(thermalModel,g);
figure(1);
pdegplot(thermalModel,'FaceLabels','on','EdgeLabels','on')
axis equal
xlim([0,10])
% Generate mesh
mesh = generateMesh(thermalModel,"Hmax",0.25,"Hedge",{[5 6 7 8],0.01});
figure(2)
pdeplot(thermalModel)
% Set thermal properties
thermalProperties(thermalModel,'ThermalConductivity',1,'Face',1);
thermalProperties(thermalModel,'ThermalConductivity',1/6,'Face',2);
% Define the internal heat generation
heatGeneration = 500;
internalHeatSource(thermalModel, heatGeneration, 'Face', 2);
% Set boundary conditions
thermalBC(thermalModel,'Edge',[1,2,3],'Temperature',25);
thermalBC(thermalModel,'Edge',4,'ConvectionCoefficient',7.371,'AmbientTemperature',30);
% Solve the steady-state heat transfer problem
thermalresults = solve(thermalModel);
T = thermalresults.Temperature;
% Plot the temperature distribution
set(gcf,'color','white')
figure(3)
pdeplot(thermalModel,'XYData',T,'Contour','on')
axis equal
title(['\fontsize{16}Steady-State Temperature. Max temp. ',num2str(max(T)),' \fontsize{16}^0C'])
% Find the heat generation required to reach the target temperature
targetTemperature = 100;
tolerance = 0.01;
while abs(max(T) - targetTemperature) > tolerance
% Adjust the heat generation based on the current temperature
heatGeneration = heatGeneration * (targetTemperature / max(T));
% Update the internal heat generation
internalHeatSource(thermalModel, heatGeneration, 'Face', 2);
% Solve the heat transfer problem
thermalresults = solve(thermalModel);
T = thermalresults.Temperature;
end
disp(['Heat generation required to reach ', num2str(targetTemperature), '°C: ', num2str(heatGeneration), ' W/m^3']);
The code has produced the following output:
Heat generation required to reach 100°C: 608.4878 W/m^3
You can try referring to:
Hope this helps!!

Products


Release

R2023a

Community Treasure Hunt

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

Start Hunting!