基于MATLAB的有限元法求解EIT技术的正问题

有限元法（Finite Element Method）是基于!近代计算机的快速发展而发展起来的一种近似数值方法, 用来解决力学，数学中的带有特定边界条件的偏微分方程问题（PDE）。而这些偏微分方程是工程实践中常见的固体力学和流体力学问题的基础。有限元和计算机发展共同构成了现代计算力学 （Computational Mechanics）的基础。有限元法的核心思想是“数值近似”和“离散化”， 所以它在历史上的发展也是围绕着这两个点进行的。

作者：菲兹

链接：https://www.zhihu.com/question/27696855/answer/43989189

常常需要求解各类微分方程，而许多微分方程的解析解一般很难得到，使用有限元法将微分方程离散化后，可以编制程序，使用计算机辅助求解。有限元法在早期是以变分原理为基础发展起来的，所以它广泛地应用于以拉普拉斯方程和泊松方程所描述的各类物理场中（这类场与泛函的极值问题有着紧密的联系）。

求解EIT技术的正问题背景：圆形被测物，16个电极排列。



%利用网格划分工具箱将圆形分为多个三角形，并标记节点坐标与三角形单元编号。保存为txt文件为tri_in1.txt。

网格划分工具箱请看链接

这里写链接内容

然后读取文本的信息，%总结点数,单元数,受约束点数%COORD结点坐标数组 LNODS 单元定义数组 通过ASSEMBLE函数ASTIF=ASSEMBLE(NPOIN,NELEM,YOUNG,POISS,THICK,COORD,LNODS)调用单元刚度矩阵函数[ESTIF,SMATX]=FORMESTIF(i,COORD,LNODS)%i为当前所计算的单元号获得刚度矩阵的集成。

%变量说明2008.11.10
%NPOIN NELEM NVFIX
%总结点数,单元数,受约束点数
%COORD结点坐标数组 LNODS 单元定义数组

%FORCE?结点力数组(n,3) n:受力结点数,(n,1):作用结点号,(n,2):x方向,(n,3):y方向
%FIXED?约束信息数组(n,3) n:受约束结点数, (n,1):约束结点号，

%(n,2)与(n,3):约束点x方向和y方向的约束情况,受约束为1,否则为0
%ASTIF总刚度矩阵,ASLOD总荷载向量
% format short e %设定输出类型
clear
FP1=fopen('tri_in1.txt','rt');
%打开数据文件 FP1数据文件指针
%读入初始数据
NPOIN=fscanf(FP1,'%d',1);
NELEM=fscanf(FP1,'%d',1);NFORCE=fscanf(FP1,'%d',1);NVFIX=fscanf(FP1,'%d',1);
YOUNG=fscanf(FP1,'%e',1);POISS=fscanf(FP1,'%f',1);THICK=fscanf(FP1,'%f',1);
LNODS=fscanf(FP1,'%f',[3,NELEM])';%单元定义数组
COORD=fscanf(FP1,'%f',[2,NPOIN])';%结点坐标数组
FORCE=fscanf(FP1,'%d',[3,NFORCE])';%结点力数组
FIXED=fscanf(FP1,'%d',[3,inf])';%?约束数组
ASTIF=ASSEMBLE(NPOIN,NELEM,YOUNG,POISS,THICK,COORD,LNODS);%?生成总刚
m=fix(NPOIN/4);%m输入电流节点
n=fix(NPOIN/3);c=fix(NPOIN/2);%n输出电流节点
ASTIF(1:c-1,c)=0;
ASTIF(c+1,c)=0;
ASTIF(c,c)=1;%c参考电势为0
b=zeros(NPOIN,1);
b(m,1)=1;
b(n,1)=-1;
[Q1,Q2,Q3,Q4]=gaus(ASTIF,b);
scatter(COORD(:,1),COORD(:,2),5,Q4);%散点图
[X,Y,Z]=griddata(COORD(:,1),COORD(:,2),Q4,linspace(-1,1,100)',linspace(-1,1,100),'v4');
figure(1)
contourf(X,Y,Z,20); %等值线图
% figure(2)
% pcolor(X,Y,Z);shading interp%伪彩色图
figure(3)
surf(X,Y,Z);%三维曲面
figure(4)
scatter(COORD(:,1),COORD(:,2),5,Q4);%散点图
fclose(FP1);%?关闭文件

tri_in1.txt文档可以自己通过网格划分获取。

图像显示：







具体代码请留言！互相交流，或者点击链接

这里写链接内容