Динамика ослабленного консольного луча

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

Модель затухания является основным вязким затуханием, распределенным однородно через объем луча. Луч деформирован путем применения внешней загрузки в совете луча и затем выпущен во время t=0. Этот пример не использует дополнительной загрузки, таким образом, смещение уменьшений луча в зависимости от времени из-за затухания. Пример использует плоское напряжение модальные, статические модели, и анализа переходных процессов в его рабочем процессе с тремя шагами:

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

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

  3. Выполните анализ переходных процессов с и без затухания.

Затухание обычно описывается как процент критического затухания, ξ, для выбранной частоты вибрации. Этот пример использует ξ=0.03, который составляет три процента критического затухания.

Пример задает значения параметров с помощью имперской системы модулей. Можно заменить их на значения, заданные в метрической системе. Если вы делаете так, гарантируете, что задаете все значения в примере с помощью той же системы.

Модальный анализ

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

modelM = createpde('structural','modal-planestress');

Создайте геометрию и включайте ее в модель. Предположим, луч 5 дюймов длиной и 0,1 дюйма толщиной.

width = 5; 
height = 0.1;

gdm = [3;4;0;width;width;0;0;0;height;height];
g = decsg(gdm,'S1',('S1')');
geometryFromEdges(modelM,g);

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

figure; 
pdegplot(modelM,'EdgeLabels','on'); 
axis equal
title 'Geometry With Edge Labels Displayed'

Figure contains an axes. The axes with title Geometry With Edge Labels Displayed contains 5 objects of type line, text.

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

hmax = height/5;
msh = generateMesh(modelM,'Hmax',hmax);

Задайте модуль Молодежи, отношение Пуассона и массовую плотность стали.

E = 3.0e7; 
nu = 0.3; 
rho = 0.3/386;
structuralProperties(modelM,'YoungsModulus',E, ...
                            'PoissonsRatio',nu, ...
                            'MassDensity',rho);

Укажите, что левый край луча является фиксированным контуром.

structuralBC(modelM,'Edge',4,'Constraint','fixed');

Решите задачу для частотного диапазона от 0 к 1e5. Рекомендуемый подход должен использовать значение, которое немного меньше, чем ожидаемая самая низкая частота. Таким образом используйте -0.1 вместо 0.

res = solve(modelM,'FrequencyRange',[-0.1,1e5]')
res = 
  ModalStructuralResults with properties:

    NaturalFrequencies: [8x1 double]
            ModeShapes: [1x1 FEStruct]
                  Mesh: [1x1 FEMesh]

По умолчанию решатель возвращает круговые частоты.

modeID = 1:numel(res.NaturalFrequencies);

Опишите получившиеся частоты в Гц путем деления их на . Отобразите частоты в таблице.

tmodalResults = table(modeID.',res.NaturalFrequencies/(2*pi));
tmodalResults.Properties.VariableNames = {'Mode','Frequency'};
disp(tmodalResults)
    Mode    Frequency
    ____    _________

     1       126.94  
     2       794.05  
     3       2216.8  
     4       4325.3  
     5       7110.7  
     6       9825.9  
     7        10551  
     8        14623  

Вычислите аналитическую основную частоту (Гц) с помощью теории луча.

I = height^3/12;
freqAnalytical = 3.516*sqrt(E*I/(width^4*rho*height))/(2*pi)
freqAnalytical = 126.9498

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

freqNumerical = res.NaturalFrequencies(1)/(2*pi)
freqNumerical = 126.9416

Вычислите период, соответствующий самому низкому режиму вибрации.

longestPeriod = 1/freqNumerical
longestPeriod = 0.0079

Постройте y-компонент решения для самой низкой частоты луча.

figure;
pdeplot(modelM,'XYData',res.ModeShapes.uy(:,1))
title('Lowest Frequency Vibration Mode')
axis equal

Figure contains an axes. The axes with title Lowest Frequency Vibration Mode contains an object of type patch.

Начальное смещение из статического решения

Луч деформирован путем применения внешней загрузки в ее совете и затем выпущен во время t=0. Найдите начальное условие для анализа переходных процессов при помощи статического решения луча с вертикальной загрузкой в совете.

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

modelS = createpde('structural','static-planestress');

Используйте ту же геометрию и mesh, которую вы использовали для модального анализа.

geometryFromEdges(modelS,g);
modelS.Mesh = msh;

Задайте те же значения для модуля Молодежи, отношения Пуассона и массовой плотности материала.

structuralProperties(modelS,'YoungsModulus',E, ...
                           'PoissonsRatio',nu, ...
                           'MassDensity',rho);
                       

Задайте то же ограничение на левый конец луча.

structuralBC(modelS,'Edge',4,'Constraint','fixed');

Примените статическую вертикальную нагрузку на правую сторону луча.

structuralBoundaryLoad(modelS,'Edge',2,'SurfaceTraction',[0;1]);

Решите статическую модель. Получившееся статическое решение служит начальным условием для анализа переходных процессов.

Rstatic = solve(modelS);

Анализ переходных процессов

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

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

modelT = createpde('structural','transient-planestress');

Используйте ту же геометрию и mesh, которую вы использовали для модального анализа.

geometryFromEdges(modelT,g);
modelT.Mesh = msh;

Задайте те же значения для модуля Молодежи, отношения Пуассона и массовой плотности материала.

structuralProperties(modelT,'YoungsModulus',E, ...
                            'PoissonsRatio',nu, ...
                            'MassDensity',rho);
                       

Задайте то же ограничение на левый конец луча.

structuralBC(modelT,'Edge',4,'Constraint','fixed');

Задайте начальное условие при помощи статического решения.

structuralIC(modelT,Rstatic)
ans = 
  NodalStructuralICs with properties:

    InitialDisplacement: [6511x2 double]
        InitialVelocity: [6511x2 double]

Решите незатухающую переходную модель в течение трех полных периодов, соответствующих самому низкому режиму вибрации.

tlist = 0:longestPeriod/100:3*longestPeriod;
resT = solve(modelT,tlist,'ModalResults',res);

Интерполируйте смещение в совете луча.

intrpUt = interpolateDisplacement(resT,[5;0.05]);

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

plot(resT.SolutionTimes,intrpUt.uy)
grid on
title('Undamped Solution')
xlabel('Time')
ylabel('Tip of beam displacement')

Figure contains an axes. The axes with title Undamped Solution contains an object of type line.

Теперь решите модель с затуханием равного 3% критического затухания.

zeta = 0.03;
omega = 2*pi*freqNumerical;
structuralDamping(modelT,'Zeta',zeta);
resT = solve(modelT,tlist,'ModalResults',res);

Интерполируйте смещение в совете луча.

intrpUt = interpolateDisplacement(resT,[5;0.05]);

Y-смещение в совете является синусоидальной функцией времени с амплитудой, экспоненциально уменьшающейся со временем.

figure
hold on
plot(resT.SolutionTimes,intrpUt.uy)
plot(tlist,intrpUt.uy(1)*exp(-zeta*omega*tlist),'Color','r')
grid on
title('Damped Solution')
xlabel('Time')
ylabel('Tip of beam displacement')

Figure contains an axes. The axes with title Damped Solution contains 2 objects of type line.

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