稳态解决方案

First, compute the steady-state solution. If the final time in the transient analysis is sufficiently large, the transient solution at the final time must be close to the steady state solution. By comparing these two results, you can check the accuracy of the transient analysis.

创建一个用于解决轴对称问题的稳态热模型。

thermalmodel = createpde('thermal','steadystate-axisymmetric');

The 2-D model is a rectangular strip whosex-dimension extends from the axis of symmetry to the outer surface andy-dimension extends over the actual length of the rod (from-1.5 m to 1.5 m). Create the geometry by specifying the coordinates of its four corners.

g = decsg（[[3 4 0 0 .2 .2 -1.5 1.5 1.5 -1.5]'）;

Include the geometry in the model.

geometryFromEdges(thermalmodel,g);

用边缘标签绘制几何形状。

figure pdegplot(thermalmodel,'edgelabels','on') axisequal

The rod is composed of a material with these thermal properties.

k = 40;%导热,W /(m*C)Rho = 7800;% Density, kg/m^3cp = 500;百分比特定热量，w*s/（kg*c）q = 20000;％热源，w/m^3

For a steady-state analysis, specify the thermal conductivity of the material.

Thermalproperties（热模型，'ThermalConductivity',k);

Specify the internal heat source.

InternalHootSource（ThermalModel，Q）;

Define the boundary conditions. There is no heat transferred in the direction normal to the axis of symmetry (edge 1). You do not need to change the default boundary condition for this edge. Edge 2 is kept at a constant temperatureT= 100°C.

thermalBC(thermalmodel,'Edge'，2，'Temperature'，100）;

指定外边界上的对流边界条件（边3）。外边界的周围温度为100°C, and the heat transfer coefficient is $50\mathrm{W}/\left(\mathrm{m}{\cdot }^{\circ }\mathrm{C}\right)$.

thermalBC(thermalmodel,'Edge'，3，...'ConvectionCoefficient',50,...'AmbientTemperature'，100）;

杆底部的热通量（边缘4）为 $5000\mathrm{W}/{\mathrm{m}}^2$.

thermalBC(thermalmodel,'Edge'，4，'Heatflux'，5000）;

Generate the mesh.

msh = generateMesh(thermalmodel); figure pdeplot(thermalmodel) axisequal

Solve the model and plot the result.

result = solve(thermalmodel); T = result.Temperature; figure pdeplot(thermalmodel,'XYData'，t，'轮廓','on') axisequaltitle'Steady-State Temperature'

Transient Solution

Switch the analysis type of the model to transient-axisymmetric.

hythmodel.analysistype ='transient-axisymmetric';

Specify the thermal conductivity, mass density, and specific heat of the material.

Thermalproperties（热模型，'ThermalConductivity',k,...'MassDensity',rho,...'femificheat',cp);

指定杆中的初始温度为0°C.

thermalIC(thermalmodel,0);

Compute the transient solution for solution times from t = 0 to t = 50000 seconds.

tfinal = 50000;tlist = 0：100：tfinal;结果= solve（ThermalModel，tlist）;

将温度分布绘制为t = 50000秒。

t =结果。图pdeplot（热模型，'XYData',T(:,end),'轮廓','on') axisequaltitle(sprintf([瞬态温度的...' at Final Time (%g seconds)'],tfinal))

Find the temperature at the bottom surface of the rod: first, at the center axis and then on the outer surface.

Tcenter = interpolateTemperature(result,[0.0;-1.5],1:numel(tlist)); Touter = interpolateTemperature(result,[0.2;-1.5],1:numel(tlist));

将杆左端的温度绘制为时间的函数。杆的外表面以100的恒定温度暴露于环境°C. When the surface temperature of the rod is less than 100°C，环境加热杆。外表面比内轴温暖。当表面温度大于100时°C，环境冷却杆。外表面比杆的内部凉爽。

figure plot(tlist,Tcenter) hold上情节（tlist，touter，'--'） 标题“底部的温度随时间的函数”XLABEL'Time, s'ylabel'Temperature, C'网格上legend(“中心轴”,'Outer Surface','地点','东南')