Решите 2-D нестационарную тепловую проблему.

Создайте нестационарную тепловую модель для этой проблемы.

thermalmodel = createpde('thermal','transient');

Создайте геометрию и включите ее в модель.

SQ1 = [3; 4; 0; 3; 3; 0; 0; 0; 3; 3];
D1 = [2; 4; 0.5; 1.5; 2.5; 1.5; 1.5; 0.5; 1.5; 2.5];
gd = [SQ1 D1];
sf = 'SQ1+D1';
ns = char('SQ1','D1');
ns = ns';
dl = decsg(gd,sf,ns);
geometryFromEdges(thermalmodel,dl);
pdegplot(thermalmodel,'EdgeLabels','on','FaceLabels','on')
xlim([-1.5 4.5])
ylim([-0.5 3.5])
axis equal

Для квадратной области присвойте следующие тепловые свойства:

thermalProperties(thermalmodel,'ThermalConductivity',10, ...
                               'MassDensity',2, ...
                               'SpecificHeat',0.1, ...
                               'Face',1);

Для алмазной области присвойте следующие тепловые свойства:

thermalProperties(thermalmodel,'ThermalConductivity',2, ...
                               'MassDensity',1, ...
                               'SpecificHeat',0.1, ...
                               'Face',2);

Предположим, что ромбовидная область является источником тепла с плотностью $4\mathrm{}{\mathrm{}}^{\mathrm{W/m2}}$.

internalHeatSource(thermalmodel,4,'Face',2);

Прикладывайте постоянную температуру $$\mathrm{}{}^{}\mathrm{0\circ C}$$ к сторонам квадратной пластины.

thermalBC(thermalmodel,'Temperature',0,'Edge',[1 2 7 8]);

Установите начальную температуру 0 ° C.

thermalIC(thermalmodel,0);

Создайте сетку.

generateMesh(thermalmodel);

Динамика для этой проблемы очень быстрая. Температура достигает стационарного состояния примерно за 0,1 секунды. Чтобы зафиксировать интересную часть динамики, установите время решения на logspace(-2,-1,10). Эта команда возвращает 10 логарифмически разнесенных значений времени решения от 0,01 до 0,1.

tlist = logspace(-2,-1,10);

Решите уравнение.

thermalresults = solve(thermalmodel,tlist)

thermalresults = 
  TransientThermalResults with properties:

      Temperature: [1481x10 double]
    SolutionTimes: [1x10 double]
       XGradients: [1481x10 double]
       YGradients: [1481x10 double]
       ZGradients: []
             Mesh: [1x1 FEMesh]

Постройте график решения с изотермическими линиями, используя контурный график.

T = thermalresults.Temperature;
pdeplot(thermalmodel,'XYData',T(:,10),'Contour','on','ColorMap','hot')

Проанализируйте теплопередачу в пруте с круглым поперечным сечением и внутренним выделением тепла, упростив 3D осесимметричную модель до 2-й модели.

Создайте нестационарную тепловую модель для решения осесимметричной задачи.

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

Модель 2-D представляет собой прямоугольную полосу, x-размер которой простирается от оси симметрии до внешней поверхности и y-размер которой простирается по фактической длине стержня (от -от 1,5 м до 1,5 м). Создайте геометрию, указав координаты ее четырех углов. Для осесимметричных моделей панель инструментов предполагает, что ось вращения является вертикальной осью, проходящей через r = 0.

g = decsg([3 4 0 0 .2 .2 -1.5 1.5 1.5 -1.5]');

Включите геометрию в модель.

geometryFromEdges(thermalmodel,g);

Постройте график геометрии с метками кромок.

figure
pdegplot(thermalmodel,'EdgeLabels','on')
axis equal

Стержень выполнен из материала с этими термическими свойствами.

k = 40; % thermal conductivity, W/(m*C)
rho = 7800; % density, kg/m^3
cp = 500; % specific heat, W*s/(kg*C)
q = 20000; % heat source, W/m^3

Укажите теплопроводность, массовую плотность и удельную теплоту материала.

thermalProperties(thermalmodel,'ThermalConductivity',k,...
                                'MassDensity',rho,...
                                'SpecificHeat',cp);

Укажите внутренний источник тепла и граничные условия.

internalHeatSource(thermalmodel,q);

Определите граничные условия. Теплопередача в направлении, перпендикулярном оси симметрии (кромке 1), отсутствует. Изменение граничного условия по умолчанию для этой кромки не требуется. Кромку 2 поддерживают при постоянной температуре Т = 100 ° С.

thermalBC(thermalmodel,'Edge',2,'Temperature',100);

Задайте условие границы конвекции на внешней границе (кромка 3). Окружающая температура на внешней границе составляет 100 ° C, а коэффициент теплопередачи - $\mathrm{}\text{50}\mathrm{}Вт/\mathrm{(}{}^{}\mathrm{}$m⋅∘C).

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

Тепловой поток в нижней части стержня (кромка 4) составляет $5000\mathrm{}{\mathrm{}}^{\mathrm{W/m2}}$.

thermalBC(thermalmodel,'Edge',4,'HeatFlux',5000);

Укажите, что начальная температура в стержне равна нулю.

thermalIC(thermalmodel,0);

Создайте сетку.

generateMesh(thermalmodel);

Вычислите переходное решение для времени решения из t = 0 кому t = 50000 секунд.

tfinal = 50000;
tlist = 0:100:tfinal;
result = solve(thermalmodel,tlist)

result = 
  TransientThermalResults with properties:

      Temperature: [259x501 double]
    SolutionTimes: [1x501 double]
       RGradients: [259x501 double]
       ZGradients: [259x501 double]
             Mesh: [1x1 FEMesh]

Постройте график распределения температуры в t = 50000 секунд.

T = result.Temperature;

figure 
pdeplot(thermalmodel,'XYData',T(:,end),'Contour','on')
axis equal
title(sprintf('Transient Temperature at Final Time (%g seconds)',tfinal))