Нагрейте проводимость в многодоменной геометрии с неоднородным потоком тепла

В этом примере показано, как выполнить 3-D переходный анализ проводимости тепла полой сферы, сделанной из трех различных слоев материала.

Сфера подвергается неоднородному внешнему потоку тепла.

Физические свойства и геометрия этой проблемы описаны в Сингхе, джайне, и Rizwan-uddin (см. Ссылку), который также имеет аналитическое решение для этой проблемы. Внутренняя поверхность сферы имеет температуру нуля в любом случае. Внешнее полушарие с положительным y значение имеет неоднородный поток тепла, заданный

qouter=θ2(π-θ)2ϕ2(π-ϕ)2

0θπ,0ϕπ.

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

Создайте тепловую модель для переходного теплового анализа.

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

Создайте многослойную сферу с помощью multisphere функция. Присвойте получившуюся геометрию тепловой модели. Сфера имеет три слоя материала с полым внутренним ядром.

gm = multisphere([1,2,4,6],'Void',[true,false,false,false]);
thermalmodel.Geometry = gm;

Постройте геометрию и покажите метки ячейки и столкнитесь с метками. Используйте FaceAlpha из 0,25 так, чтобы метки внутренних слоев отобразились.

figure('Position',[10,10,800,400]);
subplot(1,2,1)
pdegplot(thermalmodel,'FaceAlpha',0.25,'CellLabel','on')
title('Geometry with Cell Labels')
subplot(1,2,2)
pdegplot(thermalmodel,'FaceAlpha',0.25,'FaceLabel','on')
title('Geometry with Face Labels')

Figure contains 2 axes. Axes 1 with title Geometry with Cell Labels contains 2 objects of type quiver, patch. Axes 2 with title Geometry with Face Labels contains 2 objects of type quiver, patch.

Сгенерируйте mesh для геометрии. Выберите размер mesh, который достаточно крупен, чтобы ускорить решение, но достаточно прекрасный, чтобы представлять геометрию обоснованно точно.

generateMesh(thermalmodel,'Hmax',1);

Задайте теплопроводность, массовую плотность и удельную теплоемкость для каждого слоя сферы. Свойства материала являются безразмерными значениями, не данными реалистическими свойствами материала.

thermalProperties(thermalmodel,'Cell',1,'ThermalConductivity',1, ...
                                        'MassDensity',1, ...
                                        'SpecificHeat',1);
thermalProperties(thermalmodel,'Cell',2,'ThermalConductivity',2, ...
                                        'MassDensity',1, ...
                                        'SpecificHeat',0.5);
thermalProperties(thermalmodel,'Cell',3,'ThermalConductivity',4, ...
                                        'MassDensity',1, ...
                                        'SpecificHeat',4/9);

Задайте граничные условия. Самая внутренняя поверхность имеет температуру нуля в любом случае.

thermalBC(thermalmodel,'Face',1,'Temperature',0);

Наружная поверхность сферы имеет внешний поток тепла. Используйте функциональную форму тепловых граничных условий, чтобы задать поток тепла.

function Qflux = externalHeatFlux(region,~)

[phi,theta,~] = cart2sph(region.x,region.y,region.z);

theta = pi/2 - theta; % transform to 0 <= theta <= pi

ids = phi > 0;

Qflux = zeros(size(region.x));

Qflux(ids) = theta(ids).^2.*(pi - theta(ids)).^2.*phi(ids).^2.*(pi - phi(ids)).^2;

end

Постройте поток на поверхности.

[phi,theta,r] = meshgrid(linspace(0,2*pi),linspace(-pi/2,pi/2),6);
[x,y,z] = sph2cart(phi,theta,r);
region.x = x;
region.y = y;
region.z = z;
flux = externalHeatFlux(region,[]);
figure
surf(x,y,z,flux,'LineStyle','none')
axis equal
view(130,10)
colorbar
xlabel 'x'
ylabel 'y'
zlabel 'z'
title('External Flux')

Figure contains an axes. The axes with title External Flux contains an object of type surface.

Включайте это граничное условие в модель.

thermalBC(thermalmodel,'Face',4,'HeatFlux',@externalHeatFlux,'Vectorized','on');

Задайте начальную температуру, чтобы быть нулем во всех точках.

thermalIC(thermalmodel,0);

Задайте вектор такта и решите переходную тепловую задачу.

tlist = [0,2,5:5:50];
R = solve(thermalmodel,tlist);

Чтобы строить контуры в несколько раз, с уровнями контура, являющимися тем же самым для всех графиков, определяют область значений температур в решении. Минимальная температура является нулем, потому что это - граничное условие на внутренней поверхности сферы.

Tmin = 0;

Найдите максимальную температуру из итогового решения такта.

Tmax = max(R.Temperature(:,end));

Постройте контуры в области значений Tmin к Tmax в те времена в tlist.

h = figure;
for i = 1:numel(tlist)
    pdeplot3D(thermalmodel,'ColorMapData',R.Temperature(:,i))
    caxis([Tmin,Tmax])
    view(130,10)
    title(['Temperature at Time ' num2str(tlist(i))]);
    M(i) = getframe;
   
end

Чтобы посмотреть фильм контуров при выполнении этого примера на компьютере, выполните следующую линию:

movie(M,2)

Визуализируйте температурные контуры на поперечном сечении. Во-первых, задайте прямоугольную сетку точек на y-z плоскость, где x=0.

[YG,ZG] = meshgrid(linspace(-6,6,100),linspace(-6,6,100));
XG = zeros(size(YG));

Интерполируйте температуру в узлах решетки. Выполните интерполяцию для нескольких временных шагов, чтобы наблюдать эволюцию температурных контуров.

tIndex = [2,3,5,7,9,11];
varNames = {'Time_index','Time_step'};
index_step = table(tIndex.',tlist(tIndex).','VariableNames',varNames);
disp(index_step);
    Time_index    Time_step
    __________    _________

         2            2    
         3            5    
         5           15    
         7           25    
         9           35    
        11           45    
TG = interpolateTemperature(R,XG,YG,ZG,tIndex);

Задайте геометрические сферические слои на поперечном сечении.

t = linspace(0,2*pi);
ylayer1 = cos(t); zlayer1 = sin(t);
ylayer2 = 2*cos(t); zlayer2 = 2*sin(t);
ylayer3 = 4*cos(t); zlayer3 = 4*sin(t);
ylayer4 = 6*cos(t); zlayer4 = 6*sin(t);

Постройте контуры в области значений Tmin к Tmax для временных шагов, соответствующих индексам времени tIndex.

figure('Position',[10,10,1000,550]);
for i = 1:numel(tIndex)
    subplot(2,3,i)
    contour(YG,ZG,reshape(TG(:,i),size(YG)),'ShowText','on')
    colorbar
    title(['Temperature at Time ' num2str(tlist(tIndex(i)))]);
    hold on
    caxis([Tmin,Tmax])
    axis equal
    % Plot boundaries of spherical layers for reference.
    plot(ylayer1,zlayer1,'k','LineWidth',1.5)
    plot(ylayer2,zlayer2,'k','LineWidth',1.5)
    plot(ylayer3,zlayer3,'k','LineWidth',1.5)
    plot(ylayer4,zlayer4,'k','LineWidth',1.5)
end

Figure contains 6 axes. Axes 1 with title Temperature at Time 2 contains 5 objects of type contour, line. Axes 2 with title Temperature at Time 5 contains 5 objects of type contour, line. Axes 3 with title Temperature at Time 15 contains 5 objects of type contour, line. Axes 4 with title Temperature at Time 25 contains 5 objects of type contour, line. Axes 5 with title Temperature at Time 35 contains 5 objects of type contour, line. Axes 6 with title Temperature at Time 45 contains 5 objects of type contour, line.

Ссылка

[1] Сингх, Suneet, P. K. Джайн и Rizwan-uddin. "Аналитическое Решение для 3D, Неустойчивой Проводимости Тепла в Многоуровневой Сфере". ASME. J. Теплопередача. 138 (10), 2016, стр 101301-101301-11.

Для просмотра документации необходимо авторизоваться на сайте