概述
四叉树方法用来对计算域快速高效的自动离散，同时在边界附近用任意多边形离散来保证计算的准确性。这种方法的优点在于四叉树单元悬挂节点的位移协调条件自动满足，由于类型数量有限，四叉树单元的单元刚度矩阵可以预先计算并存储，用来在计算时直接选取调用，所以分析时只需要计算多边形单元的刚度矩阵，因此计算效率可以大大提高。
四叉树技术由于简单和高效的优点很早之前就吸引了研究者的注意，Yerry和Shephard最先将四叉树用于有限元网格划分。但在普通位移有限元上没有得到广泛运用，主要原因是悬挂节点难以满足位移协调条件。本文描述了一种将杂交应力元和四叉树—多边形混合网格结合的方法用来模拟颗粒增强复合材料。在杂交应力元的框架下采用类似Sukumar和Tabarraei的思想将所有单元视为任意边的多边形单元，这样悬挂节点的位移协调条件将自动满足。这个方法的另一个优点是计算效率高。因为杂交应力元本身的一个特点就是单元刚度矩阵计算比较费时，四叉树网格的单元刚度矩阵可进行预先计算，结合四叉树网格后计算时间会大幅减少。而且对计算域的离散是全自动的，省去了人工划分网格的时间。
生成过程
考虑到四叉树网格单元类型有限的特点，将它用于杂交应力元方法是一个很好的选择。下面将概述四叉树网格的生成过程。
纯四叉树
四叉树网格的生成基于递归和分解原理，图1演示了四叉树网格划分的基本步骤：
图1 四叉树网格划分过程
将表示计算域的如图1所示的数字图像划分为四个等大的正方形区域，一个区域代表一个单元。然后依次判断每一个单元内的像素值的最大差值，如果大于给定的阈值，则继续划分为四个等大正方形。重复此判断过程直到所有单元内的像素值满足阈值的要求或者单元的长度达到给定的最小值则停止划分。
需要注意的是，由于四叉树网格的划分特点，在相邻单元间可能会生成尺寸相差较大的不平衡网格，如图1-3所示的底部的两个大单元。为了使网格过渡平滑，在四叉树划分时还加入了平衡约束条件，即相邻单元边长比最大为2：1。通过平衡条件的限制生成如图1-4所示的平衡四叉树网格。同时，在相邻单元交界边上不可避免的会出现悬挂节点，由于悬挂节点的限制，四叉树网格仅限如图2所示的6种基本单元，考虑到悬挂节点的4个方向的正交旋转位置，共有24种单元，对于各向同性均质材料，4节点单元和8节点单元无需考虑旋转的影响，6节点单元中悬挂节点相对的单元只需考虑2次旋转，所有单元类型可缩减至16种。单元类型可再次进行缩减至6种，即常规的4节点与8节点单元、有一条边存在中间节点的5节点单元、有三条边存在中间节点的7节点单元及临边与对边存在中间节点的2种6节点单元，如图2所示。
图2 6种四叉树单元
对于同一种材料，可预先计算这6种单元的刚度矩阵并存储，计算时只需检索调用即可，这将大大减少刚度矩阵的计算量。
混合四叉树-多边形单元
纯四叉树网格被用来进行计算域离散时，虽然可以节省100%的单元刚度的计算时间，但是由于所有的单元都是规则的正方形形状，在离散夹杂和基体的边界时无法精确的与几何实体拟合。界面处会呈现锯齿状。如图1所示，计算时在这些区域会出现较大的不符合实际的应力集中，导致计算精度下降。虽然可以通过局部加密网格来提高整体计算精度，但是在边界处的锯齿状会一直存在。为了解决这个问题，在离散边界附近区域时，用任意边的多边形单元进行离散，如图3所示。牺牲了一部分计算效率，但是可以保证结果的精确度。
算例
为了验证混合四叉树-多边形网格在杂交应力元框架下的有效性和计算效率。将相同的模型分别用本文的方法和商用有限元软件MARC进行计算，提取结果进行对比。为了同时保证计算效率和计算精度，杂交应力元在构造应力场时选取5阶AIRY应力函数，即应力系数矩阵的项数是18项。
单夹杂模型算例
考虑具有一个夹杂的RVE模型，如图5所示，基体的材料特性如下：
杨氏弹性模量E=72GPa,泊松比V=0.3
夹杂的材料特性如下：弹性模量E=420GPa,泊松比V=0.25.
图4 单夹杂模型
为了更好的捕捉到应力集中，需要在基体和夹杂的界面附近进行局部网格加密，采用四叉树网格方法和普通位移有限元软件MRAC划分的模型如图5所示，MARC模型包含4909个四边形单元，4850个节点，四叉树模型包含1768个多边形单元，2193个节点。杂交应力元计算时选用10阶airy应力函数来表示单元内应力场，展开后含63个多项式系数。
图5 四叉树与MARC网格划分图
图6给出了两个模型的水平和竖直方向的应力分量的应力云图，其分别为杂交应力元模型及MARC程序的一般有限元计算结果。从应力云图上可以看出两种方法的计算结果非常接近，应力分布基本相同，由于边界条件是对称的所以应力分布也是对称的，在夹杂和基体周围都捕捉到了应力集中。
图6 MARC和杂交应力元水平和竖直方向应力云图
算例2
本算例模型为一个含随机夹杂的模型，在20*20mm的正方形基体中包含6个形状随机的夹杂。边界条件，材料参数与算例一相同。如图7所示。
图7 随机夹杂网格划分图
MARC单元模型如图所示，共11146个三角形单元，22130个节点，四叉树模型如图所示，共10106个单元，8182个节点。图8给出了两个模型分别用MARC软件和杂交应力元计算结果的应力云图，图15给出了水平中线路径上的水平和竖直方向的应力曲线图，从图中可以看到两种方法计算的结果基本相同，在形状比较细长的夹杂附近应力比其他夹杂处更高，说明了此方法的准确性。
图8 MARC和杂交应力元水平和竖直方向应力云图
在计算环境完全相同且计算结果基本一致的前提下，四叉树网格模型节省了大约90%的计算时间。这是因为杂交应力元计算时耗时最久的就是单刚的计算，而四叉树网格方法省去了这一步直接进行读取调用，大幅减少了计算时间。