付子义;王晨旭;王立国;长谷川弘治
【摘 要】结合声波与电磁波的物理特性,建立光子晶体与声子晶体物理量之间的联系,通过类比光子晶体能带结构的求解过程,对声子晶体能带结构进行仿真计算.以二维声子晶体为研究对象,通过重新建立数学模型,利用Born-von-karman周期性边界条件和布洛赫态,推导出声子晶体偏微分形式的本征方程,基于有限元仿真软件COMSOL Multiphysics系数偏微分模块(coefficient form PDE),求解出相应的本征频率从而求解出能带.与已有计算方法相比,该方法在适用性,难易程度等方面具有明显的优越性.
【期刊名称】《软件》
【年(卷),期】2018(039)012
【总页数】4页(P6-9)
【关键词】声子晶体;能带结构;有限元;PDE;本征方程
【作 者】付子义;王晨旭;王立国;长谷川弘治
【作者单位】河南理工大学电气工程与自动化学院,河南焦作 454150;河南理工大学电气工程与自动化学院,河南焦作 454150;河南理工大学电气工程与自动化学院,河南焦作
454150;室兰工业大学情报电子工学系,日本室兰 0500071
【正文语种】中 文
声子晶体是由嵌入均匀主体材料中的声波散射体组成的有限尺寸周期阵列的人造晶体。声波与周期性排列的散射体之间的相互作用导致形成声波不能通过该结构的频带。这种频带被称为带隙。基于这一特性使得声子晶体具有广阔的应用前景。
声子晶体的概念是由光子晶体的概念演绎而来的，因此与光子晶体相类似的是，能带结构也是声子晶体最重要的研究内容之一。文献[1]基于TMM方法，对声子晶体能带结构进行求解计算，但是TMM方法主要用于一维声子晶体的能带计算。虽然结合PWE方法可用于二维情况，但计算量较大。文献[2]基于PWE方法，但是PWE方法对于组分材料声学参数相差较大的声子晶体需要更多的计算时间和内存，而且不易给出正确结果。文献[3]采用MST方法，此方法理论推导复杂，并且只限于处理圆柱及球形散射体单元结构的声子晶体，应用上存在较大的局限性。文献[4]采用的FDTD方法可以模拟各种复杂的声子晶体结构但当组分材料的声学参数差异较大时，需要更细的网格来达到更高的精度，计算量较大。
本文主要的研究对象是声子晶体的二维结构，不仅因为它们比三维结构更简单和容易地构造，而且因为它们具有实际有用的性质。首先总结并比较了声波和电磁波的基本特性，构建出光子晶体与声子晶体之间物理量的联系。其次，在众多方法中，有限元方法凭借其良好的运算速度、精度、和收敛性，在各领域获得了广泛的应用，并已开发了许多商业化的有限元软件，因此基于有限元仿真软件COMSOL的系数偏微分模块，由声波方程出发，
结合布洛赫态，推导出偏微分形式的标量方程，开发出快速，准确的声子晶体能带结构计算方法，简化了数学建模过程，并且为进一步实现声子晶体的优化提供一个简单易行的办法。
声子晶体实际上是光子晶体的声波版本[5-7]，其主体材料可以是固体，在固体基质材料中，同时存在纵波和横向剪切波并且彼此耦合。相反，当主体材料为流体时，声波在流体中传播是只存在纵波，横波不必考虑[8]。作为本文理论模型的验证，我们采用文献[2]所用参数，以便对比，因此本文采用水银/和四氯化碳组成的二维双组分液相体系，在此只需考虑纵波即可。
本文总结并比较了声波和电磁波的基本特性。首先，我们回顾一下描述声波和电磁波行为的运动方程的归一化版本。
式中表示粒子速度，表示声压，它们是声波引起的静压的压力偏差。介质的材料参数是体模量和密度。由于通常表示粒子的位置，因此粒子速度由其时间导数给出。通过以下归一化过程简化了该表示法。这里不包括耗散，因为耗散项对人造晶体没有重要作用。为了将这些方程与电磁方程进行比较，我们通过简单均匀介质的材料参数对它们进行归一化，该介质将是声波晶体的主体材料。密度用表示，体积模数用表示，声速由表示，声学特性阻抗由 表示。粒子速度的归一化物理量为，其具有与声压相同的维度，时间的归一化物理量为。非均匀结构的材料参数由归一化密度或密度的对比度，和归一化体积模量或体积模量对比度来描述。
声波归一化后的运动方程：
用于分析声子晶体的声波运动方程，可以由等式（2）导出以下形式：
无损介质中的电磁波运动方程也以与声波运动方程相同的方式进行归一化。我们用简单均匀介质的材料参数表示光子晶体的主体材料，即介电常数，磁导率，光速并且特征阻抗由。在此对声波和电磁波使用相同的符号c0和Z0。磁场归一化后的物理量为G=Z0H，其具有与电场相同的维度。时间归一化后的物理量为。对于归一化的介电常数和磁导率表示为和。
电磁波归一化后的运动方程如下：
用于分析光子晶体的电磁波运动方程，由等式（4）可以导出以下形式：
式中表示波数，通常=1。
声波的运动方程，等式（2）仅包括梯度或散度的操作，而电磁波的运动方程，等式（3）通过旋转的矢量操作来描述。这些方程本质上看起来完全不同，没有对应关系。然而，在二维情况下，在笛卡尔坐标中发现它们之间完全对应，如表1所示。TE波的磁场H以及TM波的电场E表现得像标量的声压P。TE波的横向分量E以及TM波的H垂直于传播方向，而声波的粒子速度与传播方向平行。
光子晶体能带结构的求解过程是直接由空间部分的亥姆霍兹方程出发，分析电磁波的传播特性，将矢量形式的布洛赫态转化为标量形式的布洛赫态，推导出偏微分形式的本征方程，然后基于有限元仿真软件COMSOL的系数偏微分模块求解本征值。类比光子晶体能带的求解，以下对声子晶体建立数学模型。
本文研究对象为水银/四氯化碳体系的简单立方晶格结构的声子晶体，每一个声子晶体原胞只含有一个散射体。所采用的文献[2]的材料特性参数r_ = 0.1, K¯ = 20。假设声子晶体是由圆柱体A以正方点阵形式排列与基体B中组成的二维双组分复合材料体系。晶格常数为a，圆柱体横截面半径为r=0.3a。圆柱体轴线与z轴平行，二维周期性结构分布在x-y坐标平面内，如图1所示。
上一节已经推导出用于分析声子晶体的声波运动方程，假设声压处于时谐模式：
则等式（3）可转化为以下形式：
则等式（7）可转化为：
周期结构中的波传播可以通过使用布洛赫定理限制为单个单元，在这种周期性结构中传播的波由以下形式所描述：
式中是周期为a的周期函数，即 ，布洛赫波是由周期函数调制的平面波，其具有与晶格结构相同的周期性。
将上式布洛赫态代入声波运动方程（9）中，可得将等式转化为如下标量形式：
采用文献[2]所用参数，需要对上式进行归一化处理，频率由无量纲量（）所表示：
利用上一部分建立的数学模型，结合COMSOL系数偏微分（coefficient form PDE）模块，可以开展对二维声子晶体能带结构的仿真计算。使用者可以根据实际问题在数学方
程上进行修改，具有极大的灵活性，因此某些特殊的不便于用现成专题模块描述的问题可通过数学模块获得合理的解决。
在二维声子晶体中，此系统沿着Z轴方向是匀质，所以波矢没有限制，而该系统在xy方向具有离散平移对称性。仅沿着布里渊边缘绘制波矢k的原因是给定带状（决定带隙）的最大值和最小值总是出现在布里渊区域的边缘。在第一布里渊区域的边界会产生带隙，所以只需要沿着第一布里渊区域的边界绘制就可以表现出全部需要研究的能带情况。将所得结果与文献[2]进行了比较，讨论了其收敛性和计算时间。
计算得到能带结构如图2所示，低频处带隙的宽度远大于高频处带隙，因此本文主要在低频较宽的带隙频率范围内进行研究。在二维声子晶体结构中，布洛赫波矢被限制在第一布里渊区，并沿着第一布里渊区边界运动（），在此，将第一布里渊区分为三部分来进行仿真计算（）,（），（）由以上结果可以看出，两种方法数值结果符合得较好。图3给出了两种方法的收敛性比较。选取图2中A点，比较了分别采用文献[2]中的平面波法和有限元法对此点上频率计算中的算法收敛性，从中可以看出，平面波方法计算结果从高到低逐渐收敛，有限元法从低到高逐渐收敛，但是可以看出，有限元法的收敛速度快于平面波方法。并且还比较了两种方法的计算时间，当自由度数目相同时，有限元法的计算时间少于平面波方法。
通过总结并比较声波与电磁波之间的物理特性和联系，从光子晶体能带结构求解出发，建立求解声子晶体能带结构的数学模型，结合有限元软件COMSOL Multiphysics系数偏微分模块（coefficient form PDE）对此数学模型进行仿真计算，数值结果与传统的平面波展开法吻合较好，能大大提高收敛性和计算速度。不同于以往需要大量复杂的数学推导，
本文从声波的运动方程出发，结合布洛赫态，推导出偏微分形式的本征方程，通过COMSOL系数偏微分模块内置有限元算法直接对归一化后的本征方程进行数值计算，具有极大的灵活性。此方法对于周期性人造晶体结构能带求解具有较好的借鉴意义，为今后研究更复杂声子晶体模型提供了一种更为便捷的数值仿真手段。
【相关文献】
[1] Hou Z, Kuang W, Liu Y. Transmission property analysis of a two-dimensional phononic crystal[J]. Physics Letters A, 2004, 333(1): 172-180.
[2] 齐共金, 杨盛良, 白书欣, 等. 基于平面波算法的二维声子晶体带结构的研究[J]. 物理学报, 2003, 52(3): 000668-671.
[3] Mei J, Liu Z, Shi J, et al. Theory for elastic wave scattering by a two-dimensional periodical array of cylinders: An ideal approach for band-structure calculations[J]. Physical Review B, 2003, 67(24): 841-845.
[4] Tanaka Y, Tamura S I. Band structures of acoustic waves in phononic lattices[J]. Physica B Physics of Condensed Matter, 2002, 316(1): 237-239.
[5] Chen A, Wang Y, Yu G, et al. Elastic wave localization in two-dimensional phononic crystals with one-dimensional quasi-periodicity and random disorder[J]. 固体力学学报(英文版), 2008, 21(6): 517-528.
[6] 温激鸿, 王刚, 刘耀宗, 等. 基于集中质量法的一维声子晶体弹性波带隙计算[J]. 物理学报, 2004, 53(10): 3384-3388.
[7] 刘启能. 固-流结构声子晶体中弹性波能带的色散研究[J]. 人工晶体学报, 2009, 38(1): 112-117.
[8] Miyashita T. Sonic crystals and sonic wave-guides[J]. Measurement Science & Technology, 2005, 16(5): 47-63.
因篇幅问题不能全部显示，请点此查看更多更全内容