
FEPG有限元的自动生成与应用.ppt
89页2021/8/61一、通用有限元软件的局限性一、通用有限元软件的局限性 通用软件并不通用,只能求解很少的一部分有限元问题需要有一支训练有素、组织严密、长期稳定的有限元软件队伍,不断维护和改进其通用有限元程序只能求解今天已有的微分方程,不能求解明天的有限元问题,只能采用今天已出现的算法,不可能采用明天才有的算法不公开源代码,对用户来说程序是不可改动的,极大地限制了用户的创造性2021/8/62二、自己编程的种种困难二、自己编程的种种困难有限元程序复杂,代码庞大编程者缺乏软件专业训练,代码可读性差,再用性差,程序难以维护,程序缺乏严密组织与发展大多数工程师、教授与研究人员,他们主要是研究有限元方法或把有限元方法应用于工程设计与科学研究,而不是致力于有限元软件的开发,因此他们不可能花很多时间编写程序2021/8/63三、开发三、开发FEPGFEPG的动机与目的的动机与目的许多专家希望有自己的有限元程序,但不愿花许多时间编程FEPG应该不受领域的限制,各种专业的有限元问题都可采用本系统求解本系统不可能象通用有限元软件那样,通过不断增加单元模块和算法模块来逐步扩大其求解范围不管什么样的有限元问题,不管采用什么样的算法,都能采用本系统求解。
本系统是一种公开源代码的自由软件,任何非本系统开发者都可维护和修改2021/8/64四、由有限元公式自动产生有限元程序四、由有限元公式自动产生有限元程序有限元是求解偏微分方程的数值方法有限元公式包括: (1)弱形式的微分方程表达式 (2)形函数表达式 (3)算法表达式 a. 形成代数方程组的矩阵表达式 b. 计算步骤与迭代控制2021/8/65五、软件实现思想五、软件实现思想1、元件化程序设计方法、元件化程序设计方法 把有限元程序分解为五个元件程序2021/8/66数据转换程序可视化程序后 处 理程序前处理程序后处理计算程序单元计算程序单元子程序代数方程组求解程序有限元计算部分调用有限元程序结构有限元程序结构2021/8/67有限元程序流程框图有限元程序流程框图• 线性稳态问题流程前处理程序START元件程序E元件程序SOLV求解器U元件程序DEFER元件程序计算结果图形显示2021/8/68有限元程序流程框图有限元程序流程框图• 线性动态问题流程前处理程序START元件程序BFT元件程序E元件程序SOLV求解器U元件程序DEFER元件程序时间计算是否结束计算结果图形显示否 是2021/8/69有限元程序流程框图有限元程序流程框图• 非线性稳态问题流程前处理程序START元件程序E元件程序SOLV求解器U元件程序迭代是否结束计算结果图形显示否 是DEFER元件程序2021/8/610有限元程序流程框图有限元程序流程框图• 非线性动态问题流程前处理程序START元件程序BFT元件程序E元件程序SOLV求解器U元件程序DEFER元件程序迭代是否结束计算结果图形显示否 是时间计算是否结束 是否2021/8/611自动生成系统由近百个元件程序组成。
元件程序是以文件名为输入输出参数的可独立运行的程序由脚本(批命令)文件把元件程序组成一个“程序系统”大大降低程序的复杂性,大大提高程序的可读性和再用性实现软件自由化(即自由软件)的一种软件设计方法2021/8/612提供一种用以描述有限元问题的语言,以下称为有限元(描述)语言有限元语言是一种有层次的脚本语言有限元语言应包括微分方程表达式,形函数表达式,算法表达式等有限元公式的书写形式有限元语言应与有限元的专业语言尽可能一致,以便广大科技工作者易于理解与接受有限元语言要适用于各种领域和各种有限元方法与有限元问题关键问题是如何描述有限元问题,采用什么样的方式准确表达有限元问题及其求解算法2、人工智能技术、人工智能技术2021/8/613用用FEPG生成求解二维线弹性生成求解二维线弹性平面应变问题的有限元程序平面应变问题的有限元程序六、程序生成实例六、程序生成实例2021/8/614问题的方程描述平衡方程几何方程本构方程2021/8/615虚功方程虚功方程虚功方程弱形式虚功方程弱形式将本构关系代入上式即可得以位移为基本未知量的虚功方程2021/8/616需填写的需填写的FEPG公式文件公式文件GIO文件文件VDE文件文件GCN文件文件2021/8/617Strn.gcn文件文件defidefia ell &a ell &b str ab str astartsin astartsin asolvsin asolvsin astress bstress b2021/8/618Strn.gio文件文件dispdispstrsstrs#elemtype t3#elemtype t32dxy2dxy2021/8/619Disp.vde文件文件disp u vdisp u vcoor x ycoor x yfunc ex ey exyfunc ex ey exyshap %1 %2shap %1 %2gaus %3gaus %3mate pe pv fu fv 1.0e10;0.3;0.0;-100.0;mate pe pv fu fv 1.0e10;0.3;0.0;-100.0;$c6 fact = pe/(1.+pv)/(1.-2.*pv)$c6 fact = pe/(1.+pv)/(1.-2.*pv)vect u u vvect u u vvect ed ex eyvect ed ex eyvect fs fu fvvect fs fu fvmatr de 2 2matr de 2 2(1.-pv) pv(1.-pv) pvpv (1.-pv)pv (1.-pv)2021/8/620Disp.vde文件文件funcfuncex=+[u/x]ex=+[u/x]ey=+[v/y]ey=+[v/y]exy=+[u/y]+[v/x]exy=+[u/y]+[v/x]stifstifdist=+[ed_i;ed_j]*de_i_j*fact+[exy;exy]*(0.5-pv)*factdist=+[ed_i;ed_j]*de_i_j*fact+[exy;exy]*(0.5-pv)*factload=+[u_i]*fs_i load=+[u_i]*fs_i endend2021/8/621Disp.vde文件文件funcfuncex=+[u/x]ex=+[u/x]ey=+[v/y]ey=+[v/y]exy=+[u/y]+[v/x]exy=+[u/y]+[v/x]stifstifdist=+[ed_i;ed_j]*de_i_j*fact+[exy;exy]*(0.5-pv)*factdist=+[ed_i;ed_j]*de_i_j*fact+[exy;exy]*(0.5-pv)*factload=+[u_i]*fs_i load=+[u_i]*fs_i endend2021/8/622Strs.vde文件文件disp sa,sb,sc,disp sa,sb,sc,coef u vcoef u vcoor x ycoor x yshap %1 %2shap %1 %2gaus %3gaus %3mass %1mass %1mate pe pv fu fv 1.0e10;0.3;0.0;-100.0;mate pe pv fu fv 1.0e10;0.3;0.0;-100.0;$c6 fact = pe/(1.+pv)/(1.-2.*pv)$c6 fact = pe/(1.+pv)/(1.-2.*pv)vect u u vvect u u vvect ed ex eyvect ed ex eyvect fs fsa fsbvect fs fsa fsbvect ss sa sb scvect ss sa sb scvect ff fsa fsb fscvect ff fsa fsb fscmatr de 2 2matr de 2 2(1.-pv) pv(1.-pv) pvpv (1.-pv)pv (1.-pv)2021/8/623Strs.vde文件文件stifstif$cv ex=+{u/x}$cv ex=+{u/x}$cv ey=+{v/y}$cv ey=+{v/y}$cv fs_i=+de_i_j*ed_j*fact$cv fs_i=+de_i_j*ed_j*fact$cv func=+{u/y}+{v/x}$cv func=+{u/y}+{v/x}$c6 fsc=func*(0.5-pv)*fact$c6 fsc=func*(0.5-pv)*factdist=+[sa;sa]*0.0dist=+[sa;sa]*0.0load=+[ss_i]*ff_iload=+[ss_i]*ff_iendend2021/8/624三维正交各向异性材料计算三维正交各向异性材料计算——FEPG系统在复合材料中的应用2021/8/625三维弹性本构关系三维弹性本构关系2021/8/626其中其中2021/8/627其中其中2021/8/628需填写的需填写的FEPG公式文件公式文件GIO文件文件VDE文件文件GCN文件文件2021/8/629GCN文件文件defidefia ell &a ell &b str ab str astartsin astartsin asolvsin asolvsin astress bstress b2021/8/630Gio文件文件dispdispstrsstrs#elemtype w4#elemtype w43dxyz3dxyz2021/8/631Disp.fde文件文件disp u v wdisp u v wcoor x y zcoor x y zfunc gex gey gez geyz gexz gexyfunc gex gey gez geyz gexz gexyshap %1 %2shap %1 %2gaus %3gaus %3array de(6,6)array de(6,6)mate e1 e2 e3 v12 v13 v23 g23 g13 g12 cita fx fy fz \mate e1 e2 e3 v12 v13 v23 g23 g13 g12 cita fx fy fz \127.6d3;10.3d3;10.3d3;0.32;0.32;0.32;6.0d3;6.0d3;6.0d3;45127.6d3;10.3d3;10.3d3;0.32;0.32;0.32;6.0d3;6.0d3;6.0d3;45;0.0;0.0;5.e3;;0.0;0.0;5.e3;$c6 pi=3.14159$c6 pi=3.14159$c6 s=sin(cita*pi/180.0)$c6 s=sin(cita*pi/180.0)$c6 c=cos(cita*pi/180.0)$c6 c=cos(cita*pi/180.0)$c6 c2=cos(2.0*cita*pi/180.0)$c6 c2=cos(2.0*cita*pi/180.0)$c6 s2=sin(2.0*cita*pi/180.0)$c6 s2=sin(2.0*cita*pi/180.0)fvect le 6fvect le 6fvect ge 6fvect ge 6vect gev gex gey gez geyz gexz gexyvect gev gex gey gez geyz gexz gexyvect u u v wvect u u v wvect fu fx fy fzvect fu fx fy fz2021/8/632Disp.fde文件文件matrix te 6 6matrix te 6 6c*c s*s 0 0 0 s2c*c s*s 0 0 0 s2s*s c*c 0 0 0 (-s2/2)s*s c*c 0 0 0 (-s2/2)0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 c (-s) 00 0 0 c (-s) 00 0 0 s c 00 0 0 s c 0(-s2/2) s2/2 0 0 0 c2(-s2/2) s2/2 0 0 0 c2matr Re 6 6matr Re 6 61/e1 (-v12/e1) (-v13/e1) 0 0 01/e1 (-v12/e1) (-v13/e1) 0 0 0(-v12/e1) 1/e2 (-v23/e2) 0 0 0(-v12/e1) 1/e2 (-v23/e2) 0 0 0(-v13/e1) (-v23/e2) 1/e3 0 0 0(-v13/e1) (-v23/e2) 1/e3 0 0 00 0 0 1/g23 0 00 0 0 1/g23 0 00 0 0 0 1/g13 00 0 0 0 1/g13 00 0 0 0 0 1/g120 0 0 0 0 1/g122021/8/633Disp.fde文件文件funcfunc$cv de_i_j=re_i_j$cv de_i_j=re_i_j@l revmatri n 6 de@l revmatri n 6 de@l deform.xyz f le x y z u v w@l deform.xyz f le x y z u v w@a ge_i=+[le_j]*te_i_j@a ge_i=+[le_j]*te_i_j@w gev ge@w gev ge stifstifdist=+[gev_i;gev_j]*de_i_jdist=+[gev_i;gev_j]*de_i_jload=+[u_i]*fu_iload=+[u_i]*fu_iendend2021/8/634平面二维水沙运动问题平面二维水沙运动问题——FEPG系统在水力学中的应用2021/8/635问题描述问题描述 如右图所示,水在1/4环形水渠内流动,进口处x方向单宽流量为4, y方向单宽流量为0,水深为2,出口处水面高度为4。
求水渠内的单宽流量和水边高度 2021/8/636该问题的偏微分方程该问题的偏微分方程水流连续性方程:水流连续性方程:水流动量方程:水流动量方程:2021/8/637该问题的偏微分方程该问题的偏微分方程悬移质不平衡输移方程悬移质不平衡输移方程 悬移质引起的河床变形方程悬移质引起的河床变形方程 2021/8/638需填写的需填写的FEPG公式文件公式文件GIO文件文件VDE文件文件GCN文件文件2021/8/639GCN文件文件defidefia shola &a shola &b sholb b sholb c sholcc sholcd sholdd sholde sholee shole startsin astartsin astartsin bstartsin bstartsin cstartsin cstartsin dstartsin dstartsin estartsin ecall transcall trans2021/8/640if exist stop del stopif exist stop del stop:1:1bftbftsolvsin asolvsin acopy unod unodacopy unod unodaif exist end del endif exist end del end:2:2solvsin bsolvsin bif not exist end goto 2if not exist end goto 2solvsin csolvsin csolvsin dsolvsin dsolvsin esolvsin ecall postcall postif not exist stop goto 1if not exist stop goto 12021/8/641GIO文件文件fldfldflcflcsddsddsdcsdczbhzbh #elemtype q4#elemtype q42dxy 2dxy 2021/8/642算子分裂法求解该问题算子分裂法求解该问题采用算子分裂法分五步骤求解上述平面二维水沙运动方程 1.解水流扩散方程和水位.解水流扩散方程和水位 2021/8/643FLD.VDE文件文件disp hu,hv,zdisp hu,hv,zcoor x,ycoor x,ycoef hun,hvn,zn,hn,coef hun,hvn,zn,hn,func divfunc divvect hu hu hvvect hu hu hvvect x x yvect x x yshap %1 %2shap %1 %2gaus %3gaus %3mass %1 volmass %1 volmate epsilon eg en fx fy 0.1d-1 9.8d0 2.0d-2 0.0d0 0.0d0mate epsilon eg en fx fy 0.1d-1 9.8d0 2.0d-2 0.0d0 0.0d0vect fl fx fyvect fl fx fy2021/8/644FLD.VDE文件文件funcfunc$c6 vol = 1.0d0$c6 vol = 1.0d0$c6 const=1.0$c6 const=1.0div = +[hu/x]+[hv/y]div = +[hu/x]+[hv/y] stifstif$cv u2v2=dsqrt(hun*hun+hvn*hvn)$cv u2v2=dsqrt(hun*hun+hvn*hvn)dist=+[hu_i/x_j;hu_i/x_j]*vol*epsilondist=+[hu_i/x_j;hu_i/x_j]*vol*epsilon+[div;z]*vol+[z/x_i;hu_i]*vol*eg*hn+[div;z]*vol+[z/x_i;hu_i]*vol*eg*hn+[hu_i;hu_i]*eg*vol*en*en*u2v2/hn**(7.0d0/3.0d0)+[hu_i;hu_i]*eg*vol*en*en*u2v2/hn**(7.0d0/3.0d0)+[z/x_i;z/x_i]*vol*const+[z/x_i;z/x_i]*vol*const load=+[hu_i]*vol*fl_iload=+[hu_i]*vol*fl_i end end 2021/8/645算子分裂法求解该问题算子分裂法求解该问题2.用基于牛顿法的最小二乘法求解水流对流方程.用基于牛顿法的最小二乘法求解水流对流方程 2021/8/646FLC.VDE文件文件disp hu hvdisp hu hvcoor x y coor x y func fhu fhvfunc fhu fhvcoef hun1 hvn1 hun hvn un1 vn1 un vn hn1coef hun1 hvn1 hun hvn un1 vn1 un vn hn1mate rou 1.0 mate rou 1.0 shap %1 %2shap %1 %2gaus %3gaus %3mass %1 1.0mass %1 1.0vect hun hun hvn vect hun hun hvn vect x x y vect x x y vect fhun1 fhun1 fhvn1vect fhun1 fhun1 fhvn1vect un un vnvect un un vnvect un1 un1 vn1vect un1 un1 vn1vect hu hu hvvect hu hu hvvect fhu fhu fhvvect fhu fhu fhv2021/8/647funcfunc$cv fhun1=+{un_j/x_j}*hun+{hun/x_j}*un_j$cv fhun1=+{un_j/x_j}*hun+{hun/x_j}*un_j$cv fhun1=+fhun1*dt+hun1$cv fhun1=+fhun1*dt+hun1$cv fhvn1=+{un_j/x_j}*hvn+{hvn/x_j}*un_j$cv fhvn1=+{un_j/x_j}*hvn+{hvn/x_j}*un_j$cv fhvn1=+fhvn1*dt+hvn1$cv fhvn1=+fhvn1*dt+hvn1fhu=+[hu]+[hu_j]*{un/x_j}*dt+[hu_j/x_j]*un*dtfhu=+[hu]+[hu_j]*{un/x_j}*dt+[hu_j/x_j]*un*dt +[hu]*{un_j/x_j}*dt+[hu/x_j]*un_j*dt +[hu]*{un_j/x_j}*dt+[hu/x_j]*un_j*dt fhv=+[hv]+[hu_j]*{vn/x_j}*dt+[hu_j/x_j]*vn*dtfhv=+[hv]+[hu_j]*{vn/x_j}*dt+[hu_j/x_j]*vn*dt +[hv]*{un_j/x_j}*dt+[hv/x_j]*un_j*dt +[hv]*{un_j/x_j}*dt+[hv/x_j]*un_j*dt stifstifdist=+[fhu_i;fhu_i]dist=+[fhu_i;fhu_i] load=+[fhu_i]*fhun1_iload=+[fhu_i]*fhun1_i endend2021/8/648算子分裂法求解该问题算子分裂法求解该问题3.解挟沙方程扩散项.解挟沙方程扩散项 2021/8/649SDD.VDE文件文件disp hsdisp hscoor x,ycoor x,ycoef hn1,u,vcoef hn1,u,vfunc gradhsx gradhsyfunc gradhsx gradhsyvect gradhs gradhsx gradhsyvect gradhs gradhsx gradhsyvect x x yvect x x yshap %1 %2shap %1 %2gaus %3gaus %3mass %1 volmass %1 volmate edt rho rous rfai omeg ek eg \mate edt rho rous rfai omeg ek eg \0.1d-1 2650 1.0 1.0 1.5d-4 2.2d-1 9.8 0.1d-1 2650 1.0 1.0 1.5d-4 2.2d-1 9.8 2021/8/650funcfunc$c6 vol = 1.0d0$c6 vol = 1.0d0gradhsx = +[hs/x]gradhsx = +[hs/x] gradhsy = +[hs/y]gradhsy = +[hs/y] stifstif$c6 u2v2=dsqrt(u*u+v*v)$c6 u2v2=dsqrt(u*u+v*v)$c6 ao=rfai*omega$c6 ao=rfai*omega$c6 sstar=ek*(u2v2**3.0d0/(eg*hn1*omega))**em$c6 sstar=ek*(u2v2**3.0d0/(eg*hn1*omega))**emdist=+[gradhs_i;gradhs_i]*vol*edt+[hs;hs]*rfai*omeg*rho/hn1/rdist=+[gradhs_i;gradhs_i]*vol*edt+[hs;hs]*rfai*omeg*rho/hn1/rousous load=+[hs]*sstar*rfai*omeg*rho/rousload=+[hs]*sstar*rfai*omeg*rho/rous end end 2021/8/651算子分裂法求解该问题算子分裂法求解该问题4.解挟沙方程对流项.解挟沙方程对流项 2021/8/652SDC.VDE文件文件disp hsdisp hscoor x y coor x y func fhs func fhs coef hs1 un1 vn1 hn1 hncoef hs1 un1 vn1 hn1 hnshap %1 %2shap %1 %2gaus %3gaus %3mate rho 1.0 mate rho 1.0 vect x x y vect x x y vect un1 un1 vn1 vect un1 un1 vn1 2021/8/653funcfuncfhs=+[hs]+[hs/x_j]*un1_j*dtfhs=+[hs]+[hs/x_j]*un1_j*dt +[hs]*{un1_j/x_j}*dt +[hs]*{un1_j/x_j}*dt stifstifdist=+[fhs;fhs]dist=+[fhs;fhs] load=+[fhs]*(hs1-rho*(hn-hn1))load=+[fhs]*(hs1-rho*(hn-hn1)) end end 2021/8/654算子分裂法求解该问题算子分裂法求解该问题5.解河床变形方程.解河床变形方程 2021/8/655ZBH.VDE文件文件disp zbdisp zbcoor x,ycoor x,ycoef hs,u,v,hcoef hs,u,v,hshap %1 %2shap %1 %2gaus %3gaus %3mass %1 rhosmass %1 rhosmate rhos rfai omeg ek eg em \mate rhos rfai omeg ek eg em \1.4d1 1.0d0 1.5d-4 2.2d-1 9.8d0 7.6d-11.4d1 1.0d0 1.5d-4 2.2d-1 9.8d0 7.6d-1 2021/8/656stifstif$c6 u2v2=dsqrt(u*u+v*v)$c6 u2v2=dsqrt(u*u+v*v)$c6 ao=rfai*omeg$c6 ao=rfai*omeg$cv sstar=ek*(u2v2**3.0d0/(eg*h*omeg))**em$cv sstar=ek*(u2v2**3.0d0/(eg*h*omeg))**emdist=+[zb;zb]*0.0d0dist=+[zb;zb]*0.0d0 load=+[zb]*ao*(hs/h-sstar)load=+[zb]*ao*(hs/h-sstar) end end 2021/8/657结果图结果图2021/8/658线性正弦稳态涡流问题线性正弦稳态涡流问题((TEAM Workshop TEAM Workshop 问题问题7 7)) —— FEPG系统在电磁学中的应用2021/8/659问题描述问题描述 TEAM TEAM WorkshopWorkshop 问题7是用来检验三维线性正弦涡流场分析方法计算精度的实验模型,其结构如右图所示。
该模型由带有方孔的铝板和载流线圈组成,线圈中通有随时间按正弦规律变化的电流,其幅值为2742安匝铝板的电导率为 求涡流和磁通密度分布 2021/8/660求解方程及边界条件求解方程及边界条件2021/8/661需填写的需填写的FEPG公式文件公式文件GIO文件文件VDE文件文件GCN文件文件2021/8/662GCN文件文件defidefia ell &a ell &b str a &b str a &c str a &c str a &startc astartc asolvc asolvc astress bstress bstress cstress c2021/8/663GIO文件文件afiafimagemagemagbmagb2021/8/664方程的方程的afi.vde文件文件disp Ax,Ay,Az,udisp Ax,Ay,Az,ucoor x,y,zcoor x,y,zfunc curlx,curly,curlz,divfunc curlx,curly,curlz,divshap %1 %2shap %1 %2gaus %3gaus %3mass %1 1.0d0mass %1 1.0d0mate sigma epsilon omega fmu fjx fjy fjz mate sigma epsilon omega fmu fjx fjy fjz \0.0d0;0.0d0;50.d0;1.256d-6;0.0d0;0.0d0;0.0d0;\0.0d0;0.0d0;50.d0;1.256d-6;0.0d0;0.0d0;0.0d0;vect A Ax Ay Azvect A Ax Ay Azvect x x y zvect x x y zvect fj fjx fjy fjzvect fj fjx fjy fjzvect curl curlx curly curlzvect curl curlx curly curlz2021/8/665方程的方程的afi.vdeafi.vde文件文件disp Ax,Ay,Az,udisp Ax,Ay,Az,ucoor x,y,zcoor x,y,zfunc curlx,curly,curlz,divfunc curlx,curly,curlz,divshap %1 %2shap %1 %2gaus %3gaus %3mass %1 1.0d0mass %1 1.0d0mate sigma epsilon omega fmu fjx fjy fjz mate sigma epsilon omega fmu fjx fjy fjz \0.0d0;0.0d0;50.d0;1.256d-6;0.0d0;0.0d0;0.0d0;\0.0d0;0.0d0;50.d0;1.256d-6;0.0d0;0.0d0;0.0d0;vect A Ax Ay Azvect A Ax Ay Azvect x x y zvect x x y zvect fj fjx fjy fjzvect fj fjx fjy fjzvect curl curlx curly curlzvect curl curlx curly curlz2021/8/666afi.vde文件文件funcfunccurlx=+[Az/y]-[Ay/z]curlx=+[Az/y]-[Ay/z]curly=+[Ax/z]-[Az/x]curly=+[Ax/z]-[Az/x]curlz=+[Ay/x]-[Ax/y]curlz=+[Ay/x]-[Ax/y]div=+[Ax/x]+[Ay/y]+[Az/z]div=+[Ax/x]+[Ay/y]+[Az/z]2021/8/667afi.vde文件文件stifstif$c6 eu1=omega*sigma$c6 eu1=omega*sigma$c6 eu2=1/fmu$c6 eu2=1/fmudist=+[curl_i;curl_i]*|eu2;0.0|dist=+[curl_i;curl_i]*|eu2;0.0|+[div;div]*|eu2;0.0|+[div;div]*|eu2;0.0|+[A_i;A_i]*|0.0;eu1|+[A_i;A_i]*|0.0;eu1|+[u/x_i;A_i]*|sigma;0.0|+[u/x_i;A_i]*|sigma;0.0|+[A_i;u/x_i]*|0.0;eu1|+[A_i;u/x_i]*|0.0;eu1|+[u/x_i;u/x_i]*|sigma;0.0|+[u/x_i;u/x_i]*|sigma;0.0|load=+[A_i]*fj_i*|1.0;0.0|load=+[A_i]*fj_i*|1.0;0.0|endend2021/8/668计算涡流的mage.pde文件disp exr exi eyr eyi ezr ezidisp exr exi eyr eyi ezr ezicoef axr axi ayr ayi azr azi ur uicoef axr axi ayr ayi azr azi ur uicoor x y zcoor x y zshap %1 %2shap %1 %2gaus %3gaus %3mass %1mass %1load =fexr fexi feyr feyi fezr feziload =fexr fexi feyr feyi fezr fezimate sigma epsilon omega fmu fjx fjy fjz mate sigma epsilon omega fmu fjx fjy fjz \0.0d0;0.0d0;50.d0;1.256d-6;0.0d0;0.0d0;0.0d0;\0.0d0;0.0d0;50.d0;1.256d-6;0.0d0;0.0d0;0.0d0;2021/8/669mage.pde文件文件stifstif$cv fexr = +axi*omega*sigma-{ur/x}*sigma$cv fexr = +axi*omega*sigma-{ur/x}*sigma$cv fexi = -axr*omega*sigma-{ui/x}*sigma$cv fexi = -axr*omega*sigma-{ui/x}*sigma$cv feyr = +ayi*omega*sigma-{ur/y}*sigma$cv feyr = +ayi*omega*sigma-{ur/y}*sigma$cv feyi = -ayr*omega*sigma-{ui/y}*sigma$cv feyi = -ayr*omega*sigma-{ui/y}*sigma$cv fezr = +azi*omega*sigma-{ur/z}*sigma$cv fezr = +azi*omega*sigma-{ur/z}*sigma$cv fezi = -azr*omega*sigma-{ui/z}*sigma$cv fezi = -azr*omega*sigma-{ui/z}*sigmadist=+[exr;exr]*0.0dist=+[exr;exr]*0.0endend2021/8/670计算磁通密度的计算磁通密度的magb.pde文件文件disp bxr bxi byr byi bzr bzidisp bxr bxi byr byi bzr bzicoef axr axi ayr ayi azr azi ur uicoef axr axi ayr ayi azr azi ur uicoor x y zcoor x y zshap %1 %2shap %1 %2gaus %3gaus %3mass %1mass %1load =fbxr fbxi fbyr fbyi fbzr fbziload =fbxr fbxi fbyr fbyi fbzr fbzimate sigma epsilon omega fmu fjx fjy fjz mate sigma epsilon omega fmu fjx fjy fjz \0.0d0;0.0d0;50.d0;1.256d-6;0.0d0;0.0d0;0.0d0;\0.0d0;0.0d0;50.d0;1.256d-6;0.0d0;0.0d0;0.0d0;2021/8/671magb.pde文件文件stifstif$cv fbxr = +{azr/y}-{ayr/z}$cv fbxr = +{azr/y}-{ayr/z}$cv fbxi = +{azi/y}-{ayi/z}$cv fbxi = +{azi/y}-{ayi/z}$cv fbyr = +{axr/z}-{azr/x}$cv fbyr = +{axr/z}-{azr/x}$cv fbyi = +{axi/z}-{azi/x}$cv fbyi = +{axi/z}-{azi/x}$cv fbzr = +{ayr/x}-{axr/y}$cv fbzr = +{ayr/x}-{axr/y}$cv fbzi = +{ayi/x}-{axi/y}$cv fbzi = +{ayi/x}-{axi/y}dist=+[bxr;bxr]*0.0dist=+[bxr;bxr]*0.0endend 2021/8/672结果图结果图导体表面导体表面z=0.019mz=0.019m处涡流的矢量图处涡流的矢量图 2021/8/673结果图结果图BzBz的计算值与测量值的比较图(的计算值与测量值的比较图( y=0.072m,z=0.034m y=0.072m,z=0.034m )) 2021/8/674地下厂房洞室群的地下厂房洞室群的三维围岩稳定分析三维围岩稳定分析 —— FEPGFEPG系统在水电工程中的应用系统在水电工程中的应用2021/8/675工程背景工程背景进行索风营电站的建设和进一步优化索风营电站地下工程的设计和施工。
通过模拟计算揭示出索风营围岩稳定状态,进行围岩稳定评价,为设计的优化和施工程序的安排提供指导 2021/8/676地下厂房结构图地下厂房结构图 2021/8/677模拟计算的主要困难模拟计算的主要困难 1.考虑横观各向同性的弹塑性; 2.主变室及尾水隧洞所处岩层大部分为薄层; 3.大量支护锚杆的存在如何处理; 4.岩体渗流对围岩稳定的影响;2021/8/678岩石的横观各向同性特征岩石的横观各向同性特征 利用FEPG生成横观各向同性弹塑性计算程序2021/8/679软弱结构面三维造型软弱结构面三维造型 利用组合网格技术处理软弱面2021/8/680加锚杆岩体单元加锚杆岩体单元 利用等效方法构造加锚杆的岩体单元2021/8/681渗流计算网格图渗流计算网格图 利用生(死)节点法处理渗流自由面2021/8/682 Step15_厂房位移云图厂房位移云图 Step15_厂房变形示意图厂房变形示意图2021/8/683Step15_厂房第一主应力厂房第一主应力 Step15_厂房屈服区厂房屈服区 2021/8/684Step15_衬砌支护上的变形衬砌支护上的变形2021/8/685Step15_衬砌支护上的第一主应力衬砌支护上的第一主应力2021/8/686七、结论七、结论采用元件化程序设计方法,大大降低了有限元程序的复杂性。
采用脚本语言形式的有限元语言,免去了大量繁琐的,费力耗时的有限元编程由有限元公式自动产生有限元代码,不仅节省了编程时间,并且大大提高了程序的可读性和再用性,为有限元程序的维护和发展创造了前所未有的前景FEPG完全基于有限元方法的基本原理(虚位移原理),不受专业领域的限制,各种有限元问题和有限元方法均可采用FEPG是一个开放源码的,任何人都可参与开发的自由软件2021/8/687联系我们联系我们fem@2021/8/688谢谢大家!2021/8/689。