找回密码
 注册
Simdroid-非首页
查看: 415|回复: 0

[算法白皮书] 创新有限元算法(3)--波矢群布洛赫边界条件

[复制链接]
发表于 2023-11-25 15:41:46 | 显示全部楼层 |阅读模式
本帖最后由 现实主义的土壤 于 2023-11-29 16:27 编辑

      有限元是一种数值计算的方法,通过将不规则的区域离散成若干个子区域来求解,由于子区域的数目是有限的,故称之为有限元方法。在有限元被发明之前,所有的力学问题以及工程问题都只能依靠解析解来得到答案,如果结构稍微变得复杂些,计算解析解就显得力不从心,而对于很多结构,解析解甚至不存在。在早期,市政工程和航空工程方面复杂的弹性结构分析就存在这样的问题。1942年, A.Hrennikoff 以及 R.Courant 的工作给出了答案。虽然他们使用的方法略有差别,但他们都使用了同一个思想,就是将连续的区域离散成有限个子区域。Hrennikoff 的工作是将待求解区域离散成一个个的小方格子,类似地R.Courant将目标区域划分成了小三角形。而有限元这个概念的首次提出,还要再晚些。到了20世纪60年代初,它才被 Clough教授提出。经历了40年的发展,它的理论和算法都日趋完善。起初,有限元方法主要是解决力学方面的问题,例如汽缸扭转问题。到后来,随着计算机性能的飞速提升,有限元方法被应用到了越来越多的领域的仿真模拟(机械制造,土木建筑,材料加工等等)。到现在,从汽车到航天飞机,几乎所有的设计制造都已经离不开有限元的仿真计算。

      接下来我们将介绍几种创新有限元算法,今天我们为大家介绍--基于简单空间群的布洛赫边界条件(1)--波矢群布洛赫边界条件

      针对周期结构微纳光子器件的本征问题求解,本课题提出基于简单空间的布洛赫边界条件。分别研究了简单空间的波矢群和非波矢群,并提出了波矢群布洛赫边界条件和非波矢群布洛赫边界条件。和基于点群的布洛赫边界条件类似,简单空间群的布洛赫边界条件同样可以缩减计算区域,提高求解效率。

      考虑一个结构为简单空间群{R|Rn}的光子晶体,其中R≠E ,为格矢,n1/n2/n3为整数。由布洛赫定理知道,光子晶体中波函数具有布洛赫函数形式,即, 其中为波函数,为周期函数,为波矢,r为空间坐标。由布洛赫函数形式有布洛赫边界条件:

      公式(1.1)由光子晶体平移群得到,可将光子晶体计算区域缩减至周期单元。光子晶体的波矢群为Go(k)当对称操作Rk Go(k),由群论可知:

      其中P(Rk)/X(Rk) 为对称操作Rk对应的算符/特征标。由公式(1.1)(1.2)可得:

      公式(1.3)即为基于波矢群的布洛赫边界条件,相比使用布洛赫边界条件,利用波矢群布洛赫边界条件可以得到更小的单元。简单空间群的波矢群Go(k)仅包含镜面对称操作和旋转对称操作两种类型操作,下面依次讨论这两种情形。

图 1 镜面对称情况示意图
      首先考虑Rk为镜面对称操作σy且ky=0情况。如图1所示,D为使用布洛赫边界条件时计算区域,DM为本章节方法计算区域,M为镜面对称轴。将有限元未知量μ重组为上、左上、右上、内部、对称轴几个部分,用角标表示,区域DM边界由T,M,LT,RT四部分构成。由公式(1.1)有:


      即LT和RT边界依旧为布洛赫边界条件。由公式(1.3)有:

      其中X(σy)可以取±1。M,T上边界条件与X(σy)取值有关。当X(σy)=1时, μT=μTμM=μM,即在T,M上未知量无约束,此时波矢群布洛赫边界条件退化为PMC。当X(σy)=-1时,μT=0,μM=0,此时波矢群布洛赫边界条件退化为PEC,这也是为什么PEC和PMC可以用于缩减对称结构的计算区域。至此,计算区域由D减半为DM,仅需在DM边界上施加对应的边界条件即可。需要注意的是,光子晶体能带求解时需要分别计算X(σy)=±1情况才能得到所有能带。

图 2 旋转对称情况示意图
      Rk为旋转操作cN情况时如图2所示,DR为本文方法计算区域,o为旋转中心。将有限元未知量μ重组为上右,右上,左旋转轴,下旋转轴,内部几个部分,用角标表示,区域DR边界由TR,RT,CB,CL四部分构成。
      由公式(1.3)有:


      X (cN)可以取,n=0,1,...,N。公式(1.6)改写成如下形式:

      其中I 为单位矩阵。公式(1.7)可以记为,,,。能带计算时需要求解Aμ=λBμ的本征方程,施加波矢群布洛赫边界条件时类似点群的布鲁赫边界条件,本征方程改写为。

图 3 面心立方光子晶体计算区域及能带

      为了说明该方法在效率上的优越性,通过一个三维面心立方晶体结构光子晶体算例说明。光子晶体结构如图3a所示,格矢,介质球直径d=0.6a,相对介电常数r =12,相对磁导率μr =1。使用本文方法计算波矢空间高对称线Γ-X-W-K-Γ上的能带,在波矢面ΓKWX 上波矢群Go(k)Cs(Ic2z),因此可利用波矢群布洛赫边界条件将计算区域减半,如图3b所示。能带计算结果如图3c所示,蓝色实线为常规有限元(BBC-FEM)计算结果,蓝色方形为本文方法(GBBC-FEM)计算结果,两者计算结果一致。BBC-FEM和本文方法计算能带使用网格数、自由度、计算时间、占用内存及误差如表格1所示。BBC-FEM自由度为363082,而GBBC-FEM自由度仅为BBC-FEM一半。由于求解自由度大幅下降,GBBC-FEM的占用内存不到BBC-FEM的一半,求解时间约为BBC-FEM的一半,同时,本文方法计算结果和BBC-FEM的最大相对误差仅为0.0437‱。可以看到,三维光子晶体能带计算时采用该方法可以大幅减少求解时间和占用内存,提高求解效率。

表格 1 常规有限元和本方法结果对比表

本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有账号?注册

×
您需要登录后才可以回帖 登录 | 注册

本版积分规则

Archiver|小黑屋|联系我们|仿真互动网 ( 京ICP备15048925号-7 )

GMT+8, 2024-5-21 01:28 , Processed in 0.031419 second(s), 10 queries , Gzip On, MemCache On.

Powered by Discuz! X3.5 Licensed

© 2001-2024 Discuz! Team.

快速回复 返回顶部 返回列表