《有限单元法》--王勖成，习题2.12 MATLAB 程序
节点、单元编号
只计算到节点位移，应力没计算
%% 《有限单元法》--王勖成，习题2.12 MATLAB 程序% 仅供参考function fem()P=[-10 0 10 0 10 0 -10 0]';%结构节点载荷列向量[K,P]=xQiYi(P,make_K());solu(K,P);%% 存储各单元的节点 x、y 值function [a,b]=get_xy(m)x=[0 12 0;0 12 12;];y=[4 0 0;4 4 0];a=x(m,:);b=y(m,:);end%% 求解并显示function solu(K,P)n=1000;JieDianWeiYi=inv(K)*P% YingLi=WeiYi(3,1)/12*2e5 %x向拉、压时可用plot([0 12 12 0],[4 4 0 0],'b');hold on;xlim([-1 14]);ylim([-4 8]);plot(JieDianWeiYi([1 3 5 7],1)'*n+[0 12 12 0],JieDianWeiYi([2 4 6 8],1)'*n+[4 4 0 0],'r');end%% 用“对角元素改 1 法” ，引入位移约束并使 总的K 不奇异function [b,a]=xQiYi(P,K)%---------------------------------------------------------------------------------%d=[0 0 1 1 1 1 0 0]';%左侧两结点铰支% 共 4 个节点，顺序为 ：左上角开始1 ，顺时针排序% 共 8个位移分量，[dx1 dy1 dx2 dy2 ... ]% 位移为 0的分量，对应项填 0，其他填非 0 值% 例：右侧两点铰支，则为 d=[ 1 1 0 0 0 0 1 1]'for i=1:1:8if d(i,1)==0for j=1:1:8K(i,j)=0;K(j,i)=0;endP(i,1)=0;K(i,i)=1;endenda=P;b=K;end%% 组装总的结构刚度矩阵function K=make_K()K=[su([1 1],[1 1]) su(0,[1 2]) su([1 2],[1 3]) su([1 3],0);su(0,[2 1]) su(0,[2 2]) su(0,[2 3]) su(0,0);su([2 1],[3 1]) su(0,[3 2]) su([2 2],[3 3]) su([2 3],0);su([3 1],0) su(0,0) su([3 2],0) su([3 3],0)];end%% 分块计算矩阵的和function b=su(a,b)K1=get_K(1);K2=get_K(2);if a==0K1=zeros(6,6);a=[1,1];endif b==0K2=zeros(6,6);b=[1,1];endh11=a(1)*2-1; h12=a(1)*2; l11=a(2)*2-1; l12=a(2)*2;h21=b(1)*2-1; h22=b(1)*2; l21=b(2)*2-1; l22=b(2)*2;b=[K1(h11,l11)+K2(h21,l21) K1(h11,l12)+K2(h21,l22);K1(h12,l11)+K2(h22,l21) K1(h12,l12)+K2(h22,l22)];end%% 计算 m 单元的 刚度矩阵function K=get_K(m)%-------------------------------------------------------%v=0.3;D=2e5/(1-v^2)*[ 1 ,v ,0 ;v ,1 , 0 ;0 ,0 ,(1-v)/2];%弹性矩阵[x,y]=get_xy(m);B=get_B(x,y);A=get_A(x,y);K=B'*D*B*A;end%% 计算单元应变矩阵 Bfunction B=get_B(x,y)[~,b,c]=get_abc(x,y);B=zeros(3,6);[x,y]=get_xy(1);A=get_A(x,y);for i=1:1:3B(:,i*2-1:i*2)=[b(i) 0;0 c(i);c(i) b(i)];endB=B/(2*A);end%% 计算 ai bi cifunction [a,b,c]=get_abc(x , y)a=zeros(2,3);b=a;c=a;for l=0:1:2i=l+1;j=rem(l+1,3)+1;m=rem(l+2,3)+1;a(i)=det([x(j) y(j);x(m) y(m)]);b(i)= y(j) - y(m);c(i)= -x(j) + x(m);endend%% 计算单元面积function A=get_A(x,y)A=-det([1 x(1) y(1); 1 x(2) y(2) ; 1 x(3) y(3)])/2;endend《有限单元法》--王勖成，习题2.12 MATLAB 程序相关推荐
- 二阶偏微分方程组 龙格库塔法_有限单元法(Finite Element Method)实现声波方程模拟（Part 2）...
2.1 前言 承接上一篇文章,前面我们已经介绍了一维声波方程有限元求解: 蓝不是蓝:有限单元法(Finite Element Method)实现声波方程模拟(Part 1)zhuanlan.zhih ...
- 用Matlab求解一维非稳态周期性导热问题（有限单元法+隐式离散+高斯赛德尔迭代法）
本次求解不一定对,请先看最后说明 一.问题描述与分析 本次问题条件如下: 计算模拟如下一维常物性无内热源非稳态导热的温度场,以及内外壁面的热流密度,并进行温度场和热流的特点分析,相关参数如下. 室内温 ...
- 岩土工程渗流问题之有限单元法：理论、模块化编程实现、开源程序手把手实操技术
有限单元法在岩土工程问题中应用非常广泛,很多商业软件如Plaxis/Abaqus/Comsol等都采用有限单元解法.尽管各类商业软件使用方便,但其使用对用户来说往往是一个"黑箱子" ...
- Plaxis Python 命令流自动化处理、岩土工程渗流问题之有限单元法
目录 岩土工程渗流问题之有限单元法:理论.模块化编程实现.开源程序手把手实操应用 基于python命令流及代码的Plaxis自动化建模与典型案例实践应用 岩土工程渗流问题之有限单元法:理论.模块化编程 ...
- 现代控制理论习题解答与Matlab程序示例
现代控制理论习题解答与Matlab程序示例 大部分现代控制理论习题都可以通过计算机辅助解决,如Matlab或Octave Online. 2019更新版课程要求和Matlab简明教程: https:/ ...
- 有限单元法基本原理和数值方法_SPH法介绍
SPH法介绍 Smoothed Particle Hydrodynamics 基于网格的数值方法虽然已经有广泛的应用,但是在很多方面仍存在不足之处,比如在计算流体动力学的大变形.运动物质交界面.自由表 ...
- 有限单元法基本原理和数值方法_有限元法分析结果的四类误差，你知道吗？
本文指出了有限元法分析结果的误差影响存在于其每一操作步骤,并对这些误差进行了归类分析.随后,结合工程实例,通过改变单元类型(形状和精度).调整单元尺寸大小和应用多种分网方式,显示理想化误差和离散化误差 ...
- 有限单元法基础 -- ING
基本概念 虚位移原理 / 最小势能原理 / 卡氏第一定理(principle of virtual displacment / principle of minimum potential energ ...
- 离散方法（一）——有限单元方法(FEM)
有限元法,也叫有限单元法,它的基本思想是将一个结构或连续体的求解域离散为若干个子域(单元),并通过它们边界上的结点相互联结成为组合体. 有限元法用每一个单元内所假设的近似函数来分片地表示全求解域内待求 ...
- 有限元法、有限差分法和有限体积法的区别（转载）
有限元法.有限差分法和有限体积法的区别(转载) (2011-06-12 21:52:50) 有限差分方法(FDM)是计算机数值模拟最早采用的方法,至今仍被广泛运用.该方法将求解域划分为差分网格,用有限 ...
最新文章
- [译] Bob，函数式编程是什么鬼？
- 微软算法100题26 左旋转字符串
- 电脑故障扫描修复软件_253个电脑故障修复工具
- clion卸载+clion连接docker编写trt程序
- 开红字信息表的时候服务器为空,对方开了红字信息表但系统没有显示怎么办
- 有限自动机与有限状态机
- 我的装机软件清单-备忘
- VB 几种打开文本文件速度对比
- dedecms自定义表单提交成功后提示信息修改和跳转链接修改
- 转载--gulp入门
- linux显示防火墙端口命令,Linux查询端口是否被防火墙屏蔽 firewall-cmd命令用法
- NRPE: Unable to read output 问题处理总结
- Java-接口第一篇认识Interface
- Win8.1激活方法
- 自动化/控制工程专业英语01——什么是控制[考研/保研面试]
- C语言阿拉伯数字转大写汉字,将阿拉伯数字转大写汉字程序的类
- Python学习笔记（5）
- 关于Android短信拦截
- redis数据库(一)
- 九龙证券|这一刻，资本市场进入全新时代！