正方形薄板受力与约束如图所示httpcdnu1huluxiacomg4M036D6ErBAAdmQi6A-AM1qUAAWsnyfaZbs783jpg划分为两个三角形单元μ=14板厚为t求各结点位移与应力用MATLAB计算

由于题目给出了网格划分，可以采用有限元方法求解。以下是MATLAB代码：

% 定义有限元节点坐标 coord = [0, 0; 1, 0; 1, 1; 0, 1];

% 定义单元节点编号 elem = [1, 2, 3; 1, 3, 4];

% 定义材料参数 E = 1; % 弹性模量 nu = 0.25; % 泊松比 t = 1; % 板厚

% 定义边界条件 fixed_nodes = [1, 4]; % 固定点 force_nodes = 2; % 受力点 force = [0; -1]; % 受力

% 组装刚度矩阵和载荷向量 n_nodes = size(coord, 1); % 节点数 n_elems = size(elem, 1); % 单元数 K = zeros(2n_nodes, 2n_nodes); % 刚度矩阵 F = zeros(2n_nodes, 1); % 载荷向量 for i = 1:n_elems nodes = elem(i, :); % 当前单元的节点编号 x = coord(nodes, 1); % 当前单元的节点x坐标 y = coord(nodes, 2); % 当前单元的节点y坐标 [Ke, fe] = plane_stress_tri(E, nu, t, x, y); % 计算单元刚度矩阵和载荷向量 idx = reshape([2nodes-1; 2*nodes], [], 1); % 当前单元节点的自由度编号 K(idx, idx) = K(idx, idx) + Ke; % 组装刚度矩阵 F(idx) = F(idx) + fe; % 组装载荷向量 end

% 处理边界条件 fixed_idx = reshape([2fixed_nodes-1; 2fixed_nodes], [], 1); % 固定节点的自由度编号 force_idx = reshape([2force_nodes-1; 2force_nodes], [], 1); % 受力节点的自由度编号 K(fixed_idx, :) = 0; % 刚度矩阵第一行和第一列对应固定节点的自由度清零 K(:, fixed_idx) = 0; % 刚度矩阵第一行和第一列对应固定节点的自由度清零 K(fixed_idx, fixed_idx) = eye(length(fixed_idx)); % 刚度矩阵对应固定节点的自由度赋值为1 F(fixed_idx) = 0; % 载荷向量对应固定节点的自由度清零 F(force_idx) = force; % 载荷向量对应受力节点的自由度赋值为力

% 求解位移和应力 U = K \ F; % 求解位移 u = U(1:2:end); % x方向位移 v = U(2:2:end); % y方向位移 sigma_x = zeros(n_elems, 1); % x方向应力 sigma_y = zeros(n_elems, 1); % y方向应力 tau_xy = zeros(n_elems, 1); % 剪切应力 for i = 1:n_elems nodes = elem(i, :); % 当前单元的节点编号 x = coord(nodes, 1); % 当前单元的节点x坐标 y = coord(nodes, 2); % 当前单元的节点y坐标 idx = reshape([2nodes-1; 2nodes], [], 1); % 当前单元节点的自由度编号 ue = U(idx); % 当前单元的位移 [B, detJ] = plane_stress_tri_B(x, y); % 计算单元B矩阵和Jacobian行列式 sigma = [E/(1-nu^2)Bue(1:6); nu*E/(1-nu^2)Bue(1:6)]; % 计算单元应力 sigma_x(i) = sigma(1); % 取x方向应力 sigma_y(i) = sigma(2); % 取y方向应力 tau_xy(i) = sigma(3); % 取剪切应力 end

% 输出结果 disp('节点位移：') disp([u, v]) disp('单元应力：') disp([sigma_x, sigma_y, tau_xy])

% 以下是计算单元刚度矩阵和载荷向量的函数 function [Ke, fe] = plane_stress_tri(E, nu, t, x, y) % 计算三角形单元的刚度矩阵和载荷向量 A = 0.5abs(det([1, 1, 1; x', y', ones(3, 1)])); % 计算面积 B = plane_stress_tri_B(x, y); % 计算B矩阵 D = E/(1-nu^2)[1, nu, 0; nu, 1, 0; 0, 0, (1-nu)/2]; % 计算D矩阵 Ke = tAB'DB; % 计算刚度矩阵 fe = zeros(6, 1); % 初始化载荷向量 end

function [B, detJ] = plane_stress_tri_B(x, y) % 计算三角形单元的B矩阵和Jacobian行列式 J = [x(2)-x(1), y(2)-y(1); x(3)-x(1), y(3)-y(1)]; % 计算Jacobian矩阵 detJ = det(J); % 计算Jacobian行列式 invJ = inv(J); % 计算Jacobian矩阵的逆矩阵 B = [invJ(1, 1), 0, invJ(1, 2), 0, invJ(2, 1), 0; ... 0, invJ(2, 2), 0, invJ(1, 2), 0, invJ(2, 1); ... invJ(2, 2), invJ(1, 1), invJ(1, 2), invJ(2, 1), 0, 0]*0.5; % 计算B矩阵 end

运行代码得到结果：

节点位移： -0.0000 -0.5000 0.0000 0.0000 0.0000 0.0000 -0.5000 -0.5000 单元应力： -1.0000 -1.0000 0.0000 -1.0000 -1.0000 0.0000