Loading... 已知一个平面4节点单元体的顶点坐标分别是 1(-2,1)、2(-1,-1)、3(2,1)、4(0,2) 若单元节点速度向量是{2.0,0,2.1,0,2.2,0,2.05,0}T,试计算(ξ,η)=(0,0)点的速度分量及其作用点位置,应变速率分量、等效应变速率、体积应变速率。 ```MATLAB % 已知任一平面四边形单元的顶点坐标和各点的速度分量 % 计算母单元(xi,eta)=(,)的速度分量、对应的作用点、应变速率、等效应变速率、体积应变速率 % 母单元节点编号按照逆时针顺序编写,本题1号节点在第二象限。 clear; u0 = [-2 1 -1 -1 2 1 0 2]'; %实体单元的节点坐标 v0 = [2 0 2.1 0 2.2 0 2.05 0]'; %实体单元的节点速度 a = [-1 -1 1 1]; b = [1 -1 -1 1]; syms xi eta; for i = 1:4 Ni(i) = 0.25*(1+a(i)*xi)*(1+b(i)*eta); % Ni为节点的插值函数 end for j = 1:2:7 N(1,j) = Ni((j+1)/2); end for k = 2:2:8 N(2,k) = Ni(k/2); end v = N*v0; xi_s = 0; eta_s = 0; %所求母单元坐标点 disp('该点的速度分量'); disp(subs(v,{xi,eta},{xi_s,eta_s})); u = N*u0; disp('该点对应的实体单元坐标'); disp(subs(u,{xi,eta},{xi_s,eta_s})); % 获得形函数对xi和eta的导数 for i = 1:4 Ni_xieta(1,i) = diff(Ni(i),xi); Ni_xieta(2,i) = diff(Ni(i),eta); end for j = 1:4 A(j,:) = [u0(2*j-1),u0(2*j)]; B(j,:) = [v0(2*j-1),v0(2*j)]; end % 生成雅可比变换矩阵 J = Ni_xieta*A; Ni_xy = J\Ni_xieta; % Ni对x和y的偏导数,形成的2×4的矩阵 C = Ni_xy*B; %中间矩阵 epsilon_(1,1) = C(1,1); epsilon_(2,1) = C(2,2); epsilon_(3,1) = (sqrt(2)/2)*(C(1,2)+C(2,1)); disp('应变速率'); disp(subs(epsilon_,{xi,eta},{xi_s,eta_s})); epsilon__ = sqrt((2/3)*epsilon_.'*epsilon_); disp('等效应变速率'); disp(subs(epsilon__,{xi,eta},{xi_s,eta_s})); epsilon_V = epsilon_(1,1)+epsilon_(2,1); disp('体积应变速率'); disp(subs(epsilon_V,{xi,eta},{xi_s,eta_s})); ``` 最后修改:2021 年 03 月 03 日 © 允许规范转载 打赏 赞赏作者 支付宝微信 赞 如果觉得我的文章对你有用,请随意赞赏