2022年偏微分方程的有限元法求解 .pdf

《2022年偏微分方程的有限元法求解 .pdf》由会员分享，可在线阅读，更多相关《2022年偏微分方程的有限元法求解 .pdf（22页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 1 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 2 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 3 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 4 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结

2、- - - - - - -第 5 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 6 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 7 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 8 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 9 页，共 22 页精品资料欢迎下载精选学习资料

3、- - - - - - - - - 名师归纳总结 - - - - - - -第 10 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 11 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 12 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 13 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - -

4、 -第 14 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 15 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 16 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 17 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 18 页，共 22 页精品资料欢迎下载精选学习资料 - - -

5、- - - - - - 名师归纳总结 - - - - - - -第 19 页，共 22 页精品资料欢迎下载精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 20 页，共 22 页精品资料欢迎下载%对于 d2u/dx2=f 的 FEM 解算器，其中f=x*(1-x) % %边界条件u(0)=0, u(1)=0. %精确解用以比对xx=linspace(0,1,101);% 产生 0-1 之间的均分指令，101 为元素个数uex=(1/6).*xx.3-(1/12).*xx.4-(1/12).*xx; %对力项设置高斯点的数目NGf=2; if (NGf

6、=2) xiGf=-1/sqrt(3);1/sqrt(3);%1、2的值aGf=1 1; else, NGf=1; xiGf=0.0; aGf=2.0; end %单元数目Ne=5; %建立网格节点x=linspace(0,1,Ne+1); %零刚性矩阵K=zeros(Ne+1,Ne+1); b=zeros(Ne+1,1); %对所有单元循环计算刚性和残差for ii=1:Ne, kn1=ii; kn2=ii+1; x1=x(kn1); x2=x(kn2); dx=x2-x1;% 每一个单元的长度dxidx=2/dx;%d /dx dxdxi=1/dxidx;%dx/d dN1dxi=-1/2

7、;%d 1/ddN2dxi=1/2;%d 2/ddN1dx=dN1dxi*dxidx;%-1/(xj-xj-1) dN2dx=dN2dxi*dxidx;%1/(xj-xj-1) K(kn1,kn1)=K(kn1,kn1)-2*dN1dx*dN1dx*dxdxi;%Rj的第二项K(kn1,kn2)=K(kn1,kn2)-2*dN1dx*dN2dx*dxdxi; K(kn2,kn1)=K(kn2,kn1)-2*dN2dx*dN1dx*dxdxi; K(kn2,kn2)=K(kn2,kn2)-2*dN2dx*dN2dx*dxdxi; %用高斯积分估计力项的积分for nn=1:NGf%NGf=2 x

8、iG=xiGf(nn);% 得到高斯点的 N1=0.5*(1-xiG);% 求 N1 和 N2(即在 xiG 的权重 /插值 ) 形状函数在 的值精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 21 页，共 22 页精品资料欢迎下载N2=0.5*(1+xiG);% 值fG=xiG*(1-xiG);%对点求 f gG1=N1*fG*dxdxi;%在节点处估计权函数在高斯点的被积函数gG2=N2*fG*dxdxi;%估计是个积分值b(kn1)=b(kn1)+aGf(nn)*gG1;% aGf为 1 b(kn2)=b(kn1)+aGf(nn)*gG2;

9、end end %在 x=0 处设置 Dirichlet 条件kn1=1; K(kn1,:)=zeros(size(1,Ne+1); K(kn1,kn1)=1; b(kn1)=0; %在 x=1 处设置 Dirichlet 条件kn1=1; K(kn1,:)=zeros(size(1,Ne+1); K(kn1,kn1)=1; b(kn1)=0; %求解方程v=Kb;%v 为 Kx=b 的解plot(x,v,*-);% 画图并比较hold on; plot(xx,uex); hold off; xlabel(x); ylabel(u); 精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 22 页，共 22 页