关于matlab:基于Matlab的有限元网格自动生成算法-Q4Q8Abaqus单元网格

5次阅读

共计 6248 个字符,预计需要花费 16 分钟才能阅读完成。

今日给大家带来的次要内容是二维问题下四边形单元有限元网格如何主动生成?单元网格的造成实际上属于有限元计算中的前解决局部,即确定单元节点信息,当模型较为简单时,用户可在 Abaqus、Ansys 等大型商业有限元软件中进行建模,导出网格信息。当模型较为简单时,如二维立体板模型,用户可基于一些较为根底的网格生成算法,在本人的程序中通过管制模型长、宽等信息,即可生成有限元网格。看似利用无限,然而在一些比较复杂的畛域内,往往须要先在简略的模型中失去实践验证,如此以来,有利于自编程代码的完整性,即前解决、内核计算、后处理于一体。本篇推文,木木就带着大家学习一下 Q4、Q8 单元网格的主动生成以及 Abaqus 网格节点程序解读。代码获取:基于 Matlab 的有限元网格主动生成算法 | Q4、Q8、Abaqus 单元网格 Q4 单元网格单元主动网格划分如下图所示,为 4 节点四边形单元网格生成示意图,图中 NXE 和 NYE 别离是模型横向和纵向单元个数,dhx 和 dhy 别离是单元的横向、纵向长度。

单元主动网格划分立体板模型被划分为若干个小矩形区域,共有 4 个节点,别离是、、、,一个矩形中是 1 个四边形单元。该模型总的单元数目和节点数目别离为,。网格生成算法代码(纵向排序)global nnd nel nne  nodof eldof  n global geom connec dee nf Nodal_loadsglobal Length Width NXE NYE X_origin Y_origin dhx dhy %%nnd = 0;k = 0;% 纵向排序 for i = 1:NXE   for j=1:NYE        k = k + 1;        n1 = j + (i-1)(NYE + 1);        geom(n1,:) = [(i-1)dhx – X_origin    (j-1)dhy – Y_origin ];        n2 = j + i(NYE+1);        geom(n2,:) = [idhx – X_origin       (j-1)dhy – Y_origin  ];        n3 = n1 + 1;        geom(n3,:) = [(i-1)dhx – X_origin       jdhy – Y_origin  ];        n4 = n2 + 1;        geom(n4,:) = [idhx- X_origin       jdhy – Y_origin       ];        nel = k;        connec(nel,:) = [n1  n2  n4  n3];        nnd = n4;        endend 代码解读采纳全局变量 global 的模式,进行变量的传递;从两层循环构造上看,最外层是 for i = 1:NXE…end 阐明网格划分的过程中,x 不动,遍历每一个 y,节点 纵向排序;先确定节点号 和,n3 = n1 + 1、n4 = n2 + 1 阐明 和 在 和 的根底上,编码加 1;n1 = j + (i-1)(NYE + 1) 行不动,每次依照列减少,阐明 依照纵向排序;n2 = j + i(NYE+1) 比 多了一列的节点,阐明 与 同行;、、、的坐标 geom 由 绝对地位 - 坐标轴原点 (X_origin,Y_origin),该数值由主程序中给出;nel = k 指的是每次循环中矩形个数,当两层循环完结时,nel 指的是全副四边形单元的个数;connec(nel,:) 指的是四边形单元节点编码程序;connec(nel,:) = [n1  n2  n4 n3] 指的是该单元由、、、节点组成,并不是固定的,依照肯定的程序即可(顺时针或者逆时针);nnd = n4 当循环完结时,的数值就是节点的最大值,也就是节点的个数。网格绘制 figure(‘Name’,’Q4 单元有限元网格模型 ’);patch(‘Faces’, connec, ‘Vertices’, geom,  ‘Facecolor’,’none’,’EdgeColor’,’b’,…   ‘LineWidth’,2,’Marker’,’pentagram’,’MarkerSize’,8,’MarkerFaceColor’,’r’,…    ‘MarkerEdgeColor’,’r’);axis off% 节点编号显示 for i=1:nnd    txt =num2str(i);    text(geom(i,1)+dhx/10,geom(i,2)+dhy/10,txt);end% 单元编号显示 (与单元节点编码程序无关)for j=1:nel    txt1 = num2str(j);    text(geom(connec(j,1),1)+dhx/2,geom(connec(j,1),2)+dhy/2,txt1,’Color’,’red’,’FontWeight’, ‘bold’);end

Q4 单元网格生成(纵向排序)网格生成算法代码(横向排序)书中所给的代码一律纵向排序,我略微批改了一下,使得单元依照 方向横向排序,具体如下:global nnd nel nne  nodof eldof  n global geom connec dee nf Nodal_loadsglobal Length Width NXE NYE X_origin Y_origin dhx dhy %%nnd = 0;k = 0;% 横向排序 for i = 1:NYE   for j=1:NXE        k = k + 1;        n1 = j + (i-1)(NXE + 1);        geom(n1,:) = [(j-1)dhx – X_origin    (i-1)dhy – Y_origin ];        n2 = n1 + 1;        geom(n2,:) = [jdhx – X_origin       (i-1)dhy – Y_origin  ];                n3 = j + i(NXE+1);        geom(n3,:) = [(j-1)dhx – X_origin       idhy – Y_origin  ];        n4 = n3 + 1;        geom(n4,:) = [jdhx- X_origin       idhy – Y_origin       ];        nel = k;        connec(nel,:) = [n1  n2  n4  n3];        nnd = n4;        endend 代码解读从两层循环构造上看,最外层是 for i = 1:NYE…end 阐明网格划分的过程中,y 不动,遍历每一个 x,节点 横向排序;先确定节点号 和,n2 = n1 + 1、n4 = n3 + 1 阐明 和 在 和 的根底上,编码加 1;n1 = j + (i-1)(NXE + 1) 列不动,每次依照行减少,阐明 依照横向排序;n3 = j + i(NXE+1) 比 多了一行的节点,阐明 与 同列;、、、的坐标 geom 由 绝对地位 - 坐标轴原点 (X_origin,Y_origin),该数值由主程序中给出;nel = k 指的是每次循环中矩形个数,当两层循环完结时,nel 指的是全副四边形单元的个数;connec(nel,:) 指的是四边形单元节点编码程序;connec(nel,:) = [n1  n2  n4 n3] 指的是该单元由、、节点组成,并不是固定的,依照肯定的程序即可(顺时针或者逆时针);nnd = n4 当循环完结时,的数值就是节点的最大值,也就是节点的个数。网格绘制

Q4 单元网格生成(横向排序)绘图润饰家喻户晓,Matlab 的可视化能力强的一批,接下来木木略微批改一下 patch 函数外面的参数,即可更改填充面的色彩以及标记的形态:

批改绘图细节 Abaqus-Q4 单元在这个大节,木木想要根据本人的了解,解释一下:Abaqus-Q4 单元节点编码程序为什么是 1 -2-4-3,而不是依照顺时针或逆时针的程序编排。如下图所示,是在 Abaqus 中建设一个 CPS4(Q4) 单元的节点程序:

Abaqus-CPS4 单元节点程序模型节点读取在 Abaqus 建设立体板模型,划分为 个单元,节点、单元编码如下图所示:

Abaqus 划分网格的节点单元编码由上图可知,Abaqus 在进行 CPS4 单元节点编码时,时一排一排地排序,所以说,在一个单元中节点编码先是同一方向上的 1 -2,再是第二行的 3 -4。在 Matlab 中应用 Readmesh.m 函数,将.inp 文件的节点、单元信息导入进 Matlab 中,命令如下:[geom,connec] = Readmesh(‘Q4_abaqus.inp’); 网格绘制 figure(‘Name’,’Q4 单元有限元网格模型 ’);patch(‘Faces’, connec, ‘Vertices’, geom,  ‘Facecolor’,’c’,’Marker’,’o’,’MarkerSize’,4,’MarkerFaceColor’,’k’);axis off% 节点编号显示 [nnd,ntemp] = size(geom);[nel,netemp] = size(connec);for i=1:nnd    txt =num2str(i);    text(geom(i,1)+0.5,geom(i,2)+0.5,txt);end% 单元编号显示 (与单元节点编码程序无关)for j=1:nel    txt1 = num2str(j);    text(geom(connec(j,1),1)-2.5,geom(connec(j,1),2)-2,txt1,’Color’,’red’,’FontWeight’, ‘bold’);end

Abaqus 网格信息导入进 Matlab 由上图可看到,节点、单元程序与 Abaqus 显示统一。Q8 单元网格单元主动网格划分如下图所示,为 8 节点三角形单元网格生成示意图,图中 NXE 和 NYE 别离是模型横向和纵向单元个数,dhx 和 dhy 别离是单元的横向、纵向长度。

单元主动网格划分立体板模型被划分为若干个小矩形区域,共有 8 个节点,别离是、、、、、、、,一个矩形就是一个四边形单元。网格生成算法代码 global nnd nel nne  nodof eldof  n global geom connec dee nf Nodal_loadsglobal Length Width NXE NYE X_origin Y_origin dhx dhy %%nnd = 0;k = 0;for i = 1:NXE   for j=1:NYE    k = k + 1;%     n1 = (i-1)(3NYE+2)+2j – 1;     n2 = i(3NYE+2)+j – NYE – 1;     n3 = i(3NYE+2)+2j-1;      n4 = n3 + 1;     n5 = n3 + 2;     n6 = n2 + 1;     n7 = n1 + 2;     n8 = n1 + 1; %    geom(n1,:) = [(i-1)dhx – X_origin    (j-1)dhy – Y_origin ];    geom(n3,:) = [idhx – X_origin       (j-1)dhy – Y_origin  ];     geom(n2,:) = [(geom(n1,1)+geom(n3,1))/2  (geom(n1,2)+geom(n3,2))/2];    geom(n5,:) = [idhx- X_origin       jdhy – Y_origin ];     geom(n4,:) = [(geom(n3,1)+ geom(n5,1))/2  (geom(n3,2)+ geom(n5,2))/2];    geom(n7,:) = [(i-1)dhx – X_origin       jdhy – Y_origin  ];    geom(n6,:) = [(geom(n5,1)+ geom(n7,1))/2  (geom(n5,2)+ geom(n7,2))/2];        geom(n8,:) = [(geom(n1,1)+ geom(n7,1))/2  (geom(n1,2)+ geom(n7,2))/2];   %        nel = k;           nnd = n5;    connec(k,:) = [n1 n2 n3 n4 n5 n6 n7 n8];   endend 代码解读先确定节点号、、,n4 = n3 + 1、n5 = n3 + 2、n6 = n2 + 1、n7 = n1 + 2、n8 = n1 + 1 阐明、、、、在、、的根底上,编码加 1、加 2;留神 n1、n2、n3 的编码构造,与下面示意图相符;的坐标 geom 由 绝对地位 - 坐标轴原点 (X_origin,Y_origin),该数值由主程序中给出;nel = k 指的是每次循环中单元个数,当两层循环完结时,nel 指的是全副四边形单元的个数;connec(k,:) 指的是 8 节点四边形单元节点编码程序;connec(k,:) = [n1 n2 n3 n4 n5 n6 n7 n8] 指的是该单元由、、、、、、、节点组成,并不是固定的,依照肯定的程序即可(顺时针或者逆时针);nnd = n5 循环完结后,节点总个数等于最初一次循环的 n5 编码号。网格绘制 figure(‘Name’,’Q8 单元有限元网格模型 ’);patch(‘Faces’, connec, ‘Vertices’, geom,  ‘Facecolor’,’c’,’Marker’,’o’,’MarkerSize’,4,’MarkerFaceColor’,’k’);axis off% 节点编号显示 for i=1:nnd    txt =num2str(i);    text(geom(i,1)+dhx/10,geom(i,2)+dhy/10,txt);end% 单元编号显示 (与单元节点编码程序无关)for j=1:nel    txt1 = num2str(j);    text(geom(connec(j,1),1)+dhx/2,geom(connec(j,1),2)+dhy/2,txt1,’Color’,’red’,’FontWeight’, ‘bold’);end

Q8 单元网格主动生成本文的次要参考内容及 Matlab 代码均来自 Amar Khennane 编写的《Introduction to Finite Element Analysis Using MATLAB and Abaqus》,想要进一步理解有限元编程的小伙伴能够动手,强烈推荐!往期举荐

文章起源:技术邻 – 易木木响叮当

原文链接:https://www.jishulink.com/pos…

正文完
 0