Зажатая квадратная изотропная пластина с равномерной нагрузкой под давлением

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

Дифференциальное уравнение в частных производных для тонкой изотропной пластины с нагрузкой давления

где $$D$$ - жесткость пластины при изгибе,

и $$E$$ - модуль упругости, $$\mathrm{start}$$- отношение Пуассона, $$h$$ - толщина пластины, $$w$$ - поперечное отклонение пластины, p $$$$- нагрузка давления.

Граничными условиями для зажатых границ являются $$w\mathrm{=}$$ 0 $${и}^{}\mathrm{}$$w′=0, $${\mathrm{где}}^{}$$w ′ - производная $$$$w в направлении, перпендикулярном границе.

Дифференциальное уравнение в частных производных Toolbox™ не может непосредственно решить это уравнение четвертого порядка. Преобразуйте уравнение четвертого порядка в эти два дифференциальных уравнения второго порядка в частных производных, где $$v$$ - новая зависимая переменная.

Невозможно непосредственно задать граничные условия для $$$$$${}^{}$$′ w и w в этой системе второго порядка. Вместо этого укажите, что $${}^w$$ ′ равно 0, и $${\mathrm{определите}}^{}$$v ′ так, $$$$чтобы w также равнялось 0 на границе. Чтобы задать эти условия, используйте жесткие «пружины», распределенные по границе. Пружины прикладывают поперечную силу сдвига к краю пластины. Определите силу сдвига вдоль границы, обусловленную этими пружинами, $$как$$n⋅D∇v=-kw, $$$$где n - нормаль к границе, $$$$а k - жесткость пружин. Это выражение является обобщенным граничным условием Неймана, поддерживаемым панелью инструментов. Значение $$$$k должно быть достаточно большим, чтобы $$$$w было приблизительно равно 0 во всех точках на границе. Он также должен быть достаточно маленьким, чтобы избежать числовых ошибок из-за плохо кондиционированной матрицы жесткости.

Панель инструментов использует зависимые переменные ${\mathit{}}_{\mathrm{u1}}$ и ${\mathit{}}_{\mathrm{u2}}$ вместо $$w$$ и $$$$v. Перезаписать две дифференциальные уравнения второго порядка с помощью переменных ${\mathit{}}_{\mathrm{u1}}$ и ${\mathit{}}_{\mathrm{u2}}$:

Создайте модель PDE для системы из двух уравнений.

model = createpde(2);

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

len = 10;
gdm = [3 4 0 len len 0 0 0 len len]';
g = decsg(gdm,'S1',('S1')');
geometryFromEdges(model,g);

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

figure
pdegplot(model,'EdgeLabels','on')
ylim([-1,11])
axis equal
title 'Geometry With Edge Labels Displayed'

Коэффициенты PDE должны быть указаны в формате, требуемом панелью инструментов. Для получения более подробной информации см.

Коэффициент c в этом примере является тензором, который может быть представлен в виде матрицы 2 на 2 из блоков 2 на 2:

Эта матрица дополнительно распрямляется в вектор-столбец из шести элементов. Записи в полной матрице 2 на 2 (определяющей коэффициент a) и векторе 2 на 1 (определяющей коэффициент f) следуют непосредственно из определения системы двух уравнений.

E = 1.0e6; % Modulus of elasticity
nu = 0.3; % Poisson's ratio
thick = 0.1; % Plate thickness
pres = 2; % External pressure

D = E*thick^3/(12*(1 - nu^2));

c = [1 0 1 D 0 D]';
a = [0 0 1 0]';
f = [0 pres]';
specifyCoefficients(model,'m',0,'d',0,'c',c,'a',a,'f',f);

Для определения граничных условий сначала задайте жесткость пружины.

k = 1e7;

Определите распределенные пружины на всех четырех кромках.

bOuter = applyBoundaryCondition(model,'neumann','Edge',(1:4),...
                                     'g',[0 0],'q',[0 0; k 0]);

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

generateMesh(model);

Решите модель.

res = solvepde(model);

Доступ к решению в узловых местоположениях.

u = res.NodalSolution;

Постройте график поперечного отклонения.

numNodes = size(model.Mesh.Nodes,2);
figure
pdeplot(model,'XYData',u(:,1),'Contour','on')
title 'Transverse Deflection'

Найдите поперечное отклонение в центре пластины.

numNodes = size(model.Mesh.Nodes,2);
wMax = min(u(1:numNodes,1))

wMax = -0.2763

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

wMax = -.0138*pres*len^4/(E*thick^3)

wMax = -0.2760