В этом примере показан процесс квазистатического осесимметричного анализа теплового напряжения путем воспроизведения результатов упрощенной модели дискового тормоза, описанной в [1]. Дисковые тормоза поглощают механическую энергию через трение и преобразуют её в тепловую, которая затем рассеивается. В примере используется упрощенная модель дискового тормоза в одном процессе торможения от постоянной начальной угловой скорости до остановки. Поток операций состоит из двух этапов:
Нестационарный тепловой анализ для расчета распределения температуры в диске с использованием теплового потока от тормозных колодок
Квазистатический структурный анализ для расчета тепловых напряжений в несколько раз с использованием ранее полученного распределения температур для определения тепловых нагрузок
Полученные графики показывают распределение температуры, радиальное напряжение, напряжение кольца и напряжение по Мизесу для соответствующего времени решения.
Основываясь на допущениях, использованных в [1], пример уменьшает область анализа до прямоугольной области, соответствующей осесимметричному сечению кольцевого диска. Из-за геометрической и нагрузочной симметрии диска пример моделирует только половину толщины диска и эффект одной подушки. На следующем рисунке левый край соответствует внутреннему радиусу диска rd. Правый край соответствует внешнему радиусу диска и также совпадает с внешним радиусом площадки . Диск испытывает давление от подушки, которая генерирует тепловой поток. Вместо того чтобы явно моделировать прокладку, включите ее влияние в тепловой анализ, указав этот тепловой поток в качестве граничного условия от внутреннего радиуса прокладки до внешнего радиуса прокладки .

Создание переходной осесимметричной тепловой модели.
modelT = createpde('thermal','transient-axisymmetric');
Создайте геометрию с двумя смежными прямоугольниками. Верхний край длинного прямоугольника (справа) представляет контактную область диска.
R1 = [3,4, [ 66, 76.5, 76.5, 66, -5.5, -5.5, 0, 0]/1000]'; R2 = [3,4, [76.5, 113.5, 113.5, 76.5, -5.5, -5.5, 0, 0]/1000]'; gdm = [R1 R2]; ns = char('R1','R2'); g = decsg(gdm,'R1 + R2',ns');
Назначьте геометрию тепловой модели.
geometryFromEdges(modelT,g);
Постройте график геометрии с метками кромок и граней.
figure pdegplot(modelT,'EdgeLabels','on','FaceLabels','on')

Создайте сетку. Чтобы сопоставить сетку, используемую в [1], используйте линейный геометрический порядок вместо квадратичного порядка по умолчанию.
generateMesh(modelT,'Hmax',0.5E-04,'GeometricOrder','linear');
Укажите свойства термического материала диска.
alphad = 1.44E-5; % Diffusivity of disc Kd = 51; rhod = 7100; cpd = Kd/rhod/alphad; thermalProperties(modelT,'ThermalConductivity',Kd, ... 'MassDensity',rhod, ... 'SpecificHeat',cpd);
Задайте граничное условие теплового потока для учета области площадки. Для определения qFcn см. раздел Функция теплового потока.
thermalBC(modelT,'Edge',6,'HeatFlux',@qFcn);
Установите начальную температуру.
thermalIC(modelT,20);
Определите модель для времени, используемого в [1].
tlist = [0 0.1 0.2 1.0 2.0 3.0 3.96]; Rt = solve(modelT,tlist);
Постройте график изменения температуры со временем в трех ключевых радиальных местоположениях. Полученный график сравним с графиком, полученным в [1].
iTRd = interpolateTemperature(Rt,[0.1135;0],1:numel(Rt.SolutionTimes)); iTrp = interpolateTemperature(Rt,[0.0765;0],1:numel(Rt.SolutionTimes)); iTrd = interpolateTemperature(Rt,[0.066;0],1:numel(Rt.SolutionTimes)); figure plot(tlist,iTRd) hold on plot(tlist,iTrp) plot(tlist,iTrd) title('Temperature Variation with Time at Key Radial Locations') legend('R_d','r_p','r_d') xlabel 't, s' ylabel 'T,^{\circ}C'

Создайте осесимметричную модель статического структурного анализа.
model = createpde('structural','static-axisymmetric');
Назначьте геометрию и сетку, используемые для тепловой модели.
model.Geometry = modelT.Geometry; model.Mesh = modelT.Mesh;
Укажите структурные свойства диска.
structuralProperties(model,'YoungsModulus',99.97E9, ... 'PoissonsRatio',0.29, ... 'CTE',1.08E-5);
Наложение зависимости на модель для предотвращения жесткого движения.
structuralBC(model,'Edge',[3,4],'ZDisplacement',0);
Укажите исходную температуру, соответствующую состоянию нулевого теплового напряжения модели.
model.ReferenceTemperature = 20;
Укажите тепловую нагрузку с помощью результатов нестационарной тепловой нагрузки Rt. Время решения совпадает с временем анализа тепловой модели. Для каждого решения решите соответствующую задачу статического структурного анализа и постройте график распределения температуры, радиального напряжения, напряжения кольца и напряжения по Мизесу. Для определения plotResults см. раздел Функция результатов графика. Результаты сопоставимы с рис. 5 из [1].
for n = 2:numel(Rt.SolutionTimes) structuralBodyLoad(model,'Temperature',Rt,'TimeStep',n); R = solve(model); plotResults(model,R,modelT,Rt,n); end






Эта вспомогательная функция вычисляет переходное значение теплового потока от площадки к диску. Он использует эмпирическую формулу из [1].
function q = qFcn(r,s) alphad = 1.44E-5; % Diffusivity of disc Kd = 51; % Conductivity of disc rhod = 7100; % Density of disc cpd = Kd/rhod/alphad; % Specific heat capacity of disc alphap = 1.46E-5; % Diffusivity of pad Kp = 34.3; % Conductivity of pad rhop = 4700; % Density of pad cpp = Kp/rhop/alphap; % Specific heat capacity of pad f = 0.5; % Coefficient of friction omega0 = 88.464; % Initial angular velocity ts = 3.96; % Stopping time p0 = 1.47E6*(64.5/360); % Pressure only spans 64.5 deg occupied by pad omegat = omega0*(1 - s.time/ts); % Angular speed over time eta = sqrt(Kd*rhod*cpd)/(sqrt(Kd*rhod*cpd) + sqrt(Kp*rhop*cpp)); q = (eta)*f*omegat*r.r*p0; end
Эта вспомогательная функция отображает распределение температуры, радиальное напряжение, напряжение обруча и напряжение по Мизесу.
function plotResults(model,R,modelT,Rt,tID) figure subplot(2,2,1) pdeplot(modelT,'XYData',Rt.Temperature(:,tID), ... 'ColorMap','jet','Contour','on') title({'Temperature'; ... ['max = ' num2str(max(Rt.Temperature(:,tID))) '^{\circ}C']}) xlabel 'r, m' ylabel 'z, m' subplot(2,2,2) pdeplot(model,'XYData',R.Stress.srr, ... 'ColorMap','jet','Contour','on') title({'Radial Stress'; ... ['min = ' num2str(min(R.Stress.srr)/1E6,'%3.2f') ' MPa']; ... ['max = ' num2str(max(R.Stress.srr)/1E6,'%3.2f') ' MPa']}) xlabel 'r, m' ylabel 'z, m' subplot(2,2,3) pdeplot(model,'XYData',R.Stress.sh, ... 'ColorMap','jet','Contour','on') title({'Hoop Stress'; ... ['min = ' num2str(min(R.Stress.sh)/1E6,'%3.2f') ' MPa']; ... ['max = ' num2str(max(R.Stress.sh)/1E6,'%3.2f') ' MPa']}) xlabel 'r, m' ylabel 'z, m' subplot(2,2,4) pdeplot(model,'XYData',R.VonMisesStress, ... 'ColorMap','jet','Contour','on') title({'Von Mises Stress'; ... ['max = ' num2str(max(R.VonMisesStress)/1E6,'%3.2f') ' MPa']}) xlabel 'r, m' ylabel 'z, m' sgtitle(['Time = ' num2str(Rt.SolutionTimes(tID)) ' s']) end
[1] Adamowicz, Адам. «Осесимметричная модель КЭ для анализа тепловых напряжений в тормозном диске». Журнал теоретической и прикладной механики 53, выпуск 2 (апрель 2015): 357-370. https://doi.org/10.15632/jtam-pl.53.2.357.