matlab实现根据梁模型的长度，截面积，杨氏模量，密度，惯性矩计算质量矩阵和刚度矩阵（Ｋｅ ＝ＫｅＬ ＋ＫｅN，刚度矩阵要计算这两个之和）的计算（要求最终输出矩阵为6*6的矩阵）

根据梁模型的长度，截面积，杨氏模量，密度，惯性矩计算质量矩阵和刚度矩阵的步骤如下：

1. 定义梁的几何形状和边界条件，例如梁的长度、截面积、杨氏模量、密度等参数；

2. 根据梁的几何形状和参数计算梁的惯性矩；

3. 计算每个单元的局部坐标系下的刚度矩阵 KeL 和旋转矩阵 R；

4. 将局部坐标系下的刚度矩阵 KeL 通过旋转矩阵 R 转换到全局坐标系下的刚度矩阵 KeN；

5. 将局部坐标系下的质量矩阵 Me 和全局坐标系下的刚度矩阵 KeN 相加得到单元的刚度矩阵 Ke；

6. 将所有单元的质量矩阵和刚度矩阵组合成整体质量矩阵和整体刚度矩阵；

7. 将整体质量矩阵和整体刚度矩阵组合成一个 6 x 6 的矩阵。

下面是一个MATLAB代码示例：

% 定义梁的几何形状和边界条件
L = 5; % 梁的长度
A = 0.01; % 梁的截面积
E = 2e11; % 梁的杨氏模量
rho = 7800; % 梁的密度
I = 1/12*A*(L^3); % 梁的截面惯性矩
n = 10; % 将梁分成n个单元
dx = L/n; % 单元长度
BCL = 0; % 左端边界条件
BCR = 0; % 右端边界条件

% 计算梁的质量矩阵和刚度矩阵
K_global = zeros(6, 6);
M_global = zeros(6, 6);
for i = 1:n
    % 计算单元的惯性矩
    I_elem = 1/12*A*(dx^3) + (dx/2)^2*A*(dx/2);
    % 计算单元的局部坐标系下的刚度矩阵 KeL 和旋转矩阵 R
    KeL = E*I_elem/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];
    theta = atan(2*dx/L);
    cos_theta = cos(theta);
    sin_theta = sin(theta);
    R = [cos_theta, sin_theta, 0, 0, 0, 0;
         -sin_theta, cos_theta, 0, 0, 0, 0;
         0, 0, 1, 0, 0, 0;
         0, 0, 0, cos_theta, sin_theta, 0;
         0, 0, 0, -sin_theta, cos_theta, 0;
         0, 0, 0, 0, 0, 1];
    % 计算单元的全局坐标系下的刚度矩阵 KeN
    KeN = R' * KeL * R;
    % 计算单元的刚度矩阵 Ke
    Ke = KeL + KeN;
    % 计算单元的质量