两边简支梁有限元程序_MATLAB_FEM_简支梁

以下是一个用MATLAB编写的简支梁有限元程序：

clear all; clc;

% 梁的参数
L = 5; % 梁的长度
E = 2.1e11; % 弹性模量
I = 8.3333e-6; % 截面惯性矩
rho = 7800; % 密度
A = 0.01; % 横截面积

% 离散化参数
n = 20; % 节点数目
m = n-1; % 单元数目

% 计算单元长度
dx = L/m;

% 初始化矩阵
K = zeros(n, n); % 刚度矩阵
M = zeros(n, n); % 质量矩阵

% 计算局部刚度矩阵和质量矩阵
k = E*I/dx^3*[12 6*dx -12 6*dx; 6*dx 4*dx^2 -6*dx 2*dx^2; -12 -6*dx 12 -6*dx; 6*dx 2*dx^2 -6*dx 4*dx^2];
m = rho*A*dx/420*[156 22*dx 54 -13*dx; 22*dx 4*dx^2 13*dx -3*dx^2; 54 13*dx 156 -22*dx; -13*dx -3*dx^2 -22*dx 4*dx^2];

% 组装全局刚度矩阵和质量矩阵
for i = 1:m
    K(i:i+3, i:i+3) = K(i:i+3, i:i+3) + k;
    M(i:i+3, i:i+3) = M(i:i+3, i:i+3) + m;
end

% 施加边界条件
K([1, end], :) = 0;
K(:, [1, end]) = 0;
K(1, 1) = 1;
K(end, end) = 1;
M([1, end], :) = 0;
M(:, [1, end]) = 0;
M(1, 1) = 1;
M(end, end) = 1;

% 求解特征值和特征向量
[V, D] = eig(K, M);

% 提取固有频率和振型
f = sqrt(diag(D))/(2*pi);
y = V';

% 绘制前6阶固有频率和振型
for i = 1:6
    figure;
    plot(linspace(0, L, n), y(i, :), '-o');
    title(['第', num2str(i), '阶固有频率为', num2str(f(i)), 'Hz']);
    xlabel('位置（米）');
    ylabel('位移（米）');
end

在这个程序中，我们将梁分成20个节点，并计算每个节点的刚度矩阵和质量矩阵。然后，我们将这些局部矩阵组装成全局刚度矩阵和质量矩阵，并施加边界条件（两端简支）。最后，我们使用MATLAB的eig函数求解特征值和特征向量，提取固有频率和振型，并绘制前6阶固有频率和振型。

请注意，这个程序只是一个简单的示例。在实际应用中，您可能需要根据您的具体问题进行一些修改。