1、收稿日期：修回日期：基金项目：福建省科技计划项目“节能环保装配式立体绿化新技术研究”（）；福建农林大学横向项目“绿化种植容器用可降解耐候型木塑复合材料的研发”（）。第一作者 简 介：吴 圣 众（），男，硕 士 研 究 生，从 事 生 物 质 复 合 材 料 及 其 设 计 应 用 研 究。：。通信作者：饶久平（），男，教授，博士，博士生导师，从事生物质复合材料研究。：。：基于流固耦合分析支撑架对云杉抗风性的影响吴圣众，吴克勤，曾钦志，周建斌，吴美强，饶久平（福建农林大学材料工程学院，福建 福州；福州乐亿生态科技股份有限公司，福建 福州）开放科学标识码（码）摘要：在强风作用下新移栽的树木可能发生
2、歪斜、倒伏等破坏，造成经济损失与人员伤亡，因此，树木的抗风性以及提高树木抗风性研究具有重要意义。以云杉为例，对设立支撑架前、后的云杉分别进行三维简化建模，树高 ，冠幅 ，主干基部直径 ，并赋予模型对应的力学性能，再进行流固耦合数值模拟。结果表明：无支撑时，云杉树梢位移、枝干所受剪切应力和树底弯矩均随着风速的增强而增大。在 风速下，云杉横纹方向剪切应力最大值为 ，大于云杉横纹抗剪强度的，云杉可能发生剪切破坏，通过后处理发现破坏危险点位置在 级分枝与主干的交界处，云杉的枝干易被折断。分别在云杉树高的、处设立支撑架，对比无支撑架时云杉的风致响应，云杉的树梢位移、枝干所受剪切应力和树底弯矩均有所减小，
3、云杉发生横纹剪切破坏的风险降低，发生倒伏的概率大幅下降。当支撑架置于树高 处时云杉所受的剪切应力最小，相对无支撑架时减小了，云杉枝干被折断的风险大幅下降。然而，处于支撑架上方的云杉部分随风摆动的振频较无支撑时有所增大，云杉的稳定性下降，长时间的高频振动亦会增大枝干折断的风险。计算机数值模拟风致树木响应具有很强的自定义性，其模拟结果对树木的抗风性及其支撑架的设计应用具有一定的参考价值。关键词：流固耦合；树木支撑架；抗风性；云杉中图分类号：；文献标识码：文章编号：（），（，；，）：，：；森林与环境学报，（）：第 卷 第 期 年 月在强风作用下，树木可能会发生折枝、歪斜和倒伏等破坏，造成经济损失甚至
4、人员伤亡。年 月，辽宁省沈阳市受强风暴雨影响，致使全市树木倒伏两千多棵，其中在公园、绿地等较为空旷的区域，树木倒伏情况较为严重；年 月，台风“山竹”对深圳城市树木造成严重破坏，其中道路、绿地树木歪斜接近，倒伏高达，其中乔木受损尤为严重。随着城市园林面积的不断扩大，越来越多的树木被移栽到城市中用作行道树或园林景观树。然而，在树木移栽初期，由于其根部尚未扎深扎实，不能提供足够的支撑力，故新移栽的树木在强风作用下更易发生歪斜和倒伏。云杉（）属常绿乔木，树干高大通直，造型美观，故常被用作行道树或园林景观树。为了避免强风对新移栽的云杉造成破坏，应在移栽后对其采取一定的支撑措施。最常见的树木支撑方法有四角
5、支撑法、三角支撑法和扁担式扶架法等，其中四角支撑法对移栽树木的支撑效果最好，且安装简便、造型美观，故四角支撑法在园林绿化保护领域的使用最为广泛。目前国内外研究树木抗风性的方法一般有现场实测法、风洞试验法和数值模拟法。对欧洲赤松、白桦、挪威云杉进行了拔树试验，得出对拉力抵抗程度的排列顺序为白桦欧洲赤松挪威云杉。树木的抗拉能力与胸径大小有关，且当土壤为冻土情况时，树木在拉力作用下主干更容易发生断裂破坏。郑义津等以豆瓣黄杨苗为试验树，在风洞中通过人工落沙的方式模拟风沙流场，控制不同的落沙量形成不同的风沙流场，测量不同树木间距条件下树后风沙流场的情况，并与无树时进行对比。结果表明，树后沙浓度的降低在树
6、冠层最为显著，同时沙浓度也随着树木间距的增大而下降。通过地面激光扫描获得树木的 几何形状，将其转化为 模型后以固定在地面上梁结构的形式构建了有限元模型，并对计算模型进行有限元分析，显示了自重和静态风荷载等荷载工况组合引起树木整体结构的变形，通过结构分析结果可以发现树木潜在的破坏危险点。相较于现场实测法与风洞试验法，数值模拟法拥有很强的自定义性，可以模拟任意风速下不同树况的风致响应。随着计算机科学的快速发展，数值模拟法的精度、速度有了大幅的提升，并且具有成本低、适用范围广等优点，故数值模拟法成为当前研究树木抗风性的热门方法。相较于其他学者单纯研究无支撑保护时树木在风中的响应，本研究运用有限元模拟
7、软件，采用双向流固耦合分析方法，研究树木在设立支撑架前、后的抗风性差异，旨在为园林绿化树木支撑实施方案提供思路。材料与方法 树木及支撑架建模与模型条件树木是一种复杂的空间结构，不同树种之间的差异可能很大，即便是相同的树种，形态也有所区别。但对于同一棵树，由于枝条自身的重力作用，随着着枝深度的增大，枝长随之增大，着枝角度也逐渐增大。城市中新移栽的行道树为了提高存活率，防止在移栽初期消耗过多的养分，一般需要对枝叶进行修剪，故新移栽的树木枝叶并不茂盛。为了确保计算机流固耦合模拟的正常计算，本研究对模型做出如下假设：由于树木根部与周围土壤的结合作用十分复杂，故假设树木底部与地面的连接方式为固结，同理，
8、假设树木支撑杆与地面的连接方式为固结。常见的支撑杆与树干的连接采用钉子、树箍、铁丝等实现紧密、固定连接，因此，本研究只对树木支撑架的支撑杆部分进行建模，并假设支撑杆与树干的连接方式为固结。新移栽的树木枝叶并不茂盛，并且由于树叶与风场的耦合作用十分复杂，若要生成具有大量树叶的有限元模型，并求解所有枝叶的流固耦合是不现实的，因此，本研究未对树叶进行建模，重点研究树木枝干对树木抗风性的影响。本研究模拟短时间（）内树木的抗风性，故假设该期间的风为平均风。参照图，运用 中的 模块对新移栽的云杉进行简化建模（图），云杉树木高 ，冠幅 ，着枝深度 ，云杉主干基部直径 ，锥度。由于树木支撑架的使用方法尚无统一
9、的国家或行业标准，故树木支撑杆建模与支撑方法参考当前市面上常见的木质支撑杆外形与支撑方法。如图 所示，树木支撑杆设置为圆柱体，直径 ，支撑角度，支撑高度分别在树高的、处，支撑方式为四角支撑。云杉模型置于风场中，风场流体域尺寸为 ，云杉模型距左右边界各 ，距风场入口 ，如图 所示。云杉及支撑架模型在风场 第 期吴圣众，等：基于流固耦合分析支撑架对云杉抗风性的影响中的朝向如图 所示。注：移栽的云杉；云杉模型（无支撑架）；云杉模型（有支撑架）；风场模型；云杉及支撑架在风场中的朝向。：；图 云杉、支撑架、风场模型及相对位置 ，有限元模拟树木在强风作用下的运动可以看作由树冠和树干组成的固体域同流体域的空
10、气相互作用过程，而树木与空气之间存在反复的交互作用，即空气的流动造成树木的形变，而树木的形变又反过来影响空气的流动，因此树木在风中的响应采用双向流固耦合分析。先对模型进行前处理，再设定求解条件。模型前处理 模型的前处理主要包括网格划分和材料的定义与赋予。双向流固耦合模拟需要用到、和 模块，其中 模块计算时只考虑风场对树木外形的作用，故在 模块对风场模型进行 运算。分别对风场模型和树木模型进行网格划分，网格划分质量越高，分析计算速度越快，精度越高。由于树木模型不仅为非轴对称结构，而且是完全不规则的，所以划分树木模型网格时采用 自带的 模块进行划分。同理，虽然风场模型外部是规则的长方体，但是由于树
11、木 模 型 与 风 场 模 型 流 固 耦 合 时 进 行 了运算，所以风场内部结构也并不规则，故划分风场模型网格时也采用 自带的 模块进行划分。云杉模型网格及其所处的风场模型网格如图 所示，风场网格的数量为 万个，平均网格质量为，网格类型为四面体网格；云杉模型网格的数量为 万个，平均网格质量为，网格类型为四面体网格。二者的网格质量均满足双向流固耦合的要求。注：风场模型网格划分；云杉模型网格划分。：；图 风场模型网格划分与云杉模型网格划分 材料的定义包括 模块中流体域材料的定义与 模块中云杉与支撑架材料的定义。流体域的材料定义为理想空气，模拟自然风的作用。因为树木是一种正交各向异性材料，所以具
12、有 个弹性常数，本研究模拟的云杉力学性能参考前人测量的树木力学数据，如表 所示。将云杉力学参数通过 模块输入，在 模块中为主干与一级分枝分别建立柱坐标系，再将树木材料赋予树木模型。树木支撑架的材料定义与赋予方法同树木模型，材质为云杉木棒，力学性能参考云杉的力学性能。森 林 与 环 境 学 报第 卷表 云杉的主要力学参数 顺纹弹性模量 水平径向弹性模量 水平弦向弹性模量 顺纹弦面剪切弹性模量 顺纹径面剪切弹性模量 水平 面 剪 切 弹 性 模 量 弦向压力应变 径向延展应 变 纵向压力应变 径向延展应变 弦向压力应变 纵向延展应 变 密度（）模型求解设定 树木及其支撑架在风中的响应为双向流固耦合
13、过程，在 中分别调整 模块与 模块，再将二者数据导入 模块进行双向流固耦合模拟分析求解。模块设定：入口边界设置为速度入口，速度分别设置为、（分别对应、级风，其中 级风为台风）；出口边界设置为压力出口。四周壁面条件为静止无滑移，粗糙程度为。仿真选用瞬态计算和 湍流模型。将树木与流场的交界面创建为系统耦合动网格区域，网格方法设置为光顺与网格重构。在求解参数设置中，由于模拟的风速较大，故求解方法选用耦合显式算法，最后将 混合初始化。模块设定：在完成树木力学模型的建立之后，需要对树木 模型施加约束。设置树木模型底部与支撑架底部为固定支撑，选取树木和支撑架与风接触的面为流固耦合界面，为树木及支撑架增加地
14、球重力。在分析设置中设置定义依据为按时间定义，关闭自动时步，步时设置为 ，步骤结束时间设置为 。模块设定：对 模块创建的树木、流场交界面和 模块创建的流固耦合界面进行数据交互，分别设置模拟步长与总时长为、，最后进行双向流固耦合分析。结果与分析 云杉模型固有频率验证模型的准确性是研究树木风致响应的基础，国内外许多学者针对树木的固有频率与树干参数进行了大量的研究，发现树木的固有频率与胸径成正比，与树高的平方成反比，采用最小二乘法拟合得到如下公式：（）式中：为树木固有频率；为树木胸径（）；为树木高度（）。通过公式（）计算得到云杉的固有频率为，将云杉模型导入模态分析模块，得出云杉模型的固有频率为。计算
15、得出的固有频率与有限元模型的固有频率较为接近，故有限元模型较符合真实情况，验证了建模的合理性。无支撑时云杉的抗风性 风致云杉树梢位移为研究云杉在不同风况条件下的抗风性差异，分别模拟云杉在、风速下的响应。在强风作用下，树木的最大位移可能出现在树梢或 级分枝顶端，而树梢处的位移更能直观地体现树木对风的抵御能力。因此，在模拟的后处理中提取 树梢的位移时程数据。云杉的风致树梢位移如图 所示。在风场作用下，云杉由静止开始摆动，在不同风速下，云杉树梢最大位移均出现在图 无支撑时不同风速下的风致云杉树梢位移 第 期吴圣众，等：基于流固耦合分析支撑架对云杉抗风性的影响 左右处，随着风速的增强，树梢位移非线性增
16、大，而是成倍数增大。树干的振动频率几乎保持不变，树梢来回摆动的幅度随时间逐渐减小趋于缓和。云杉的风振特性与其他学者模拟、实测的树木风振特性具有相似性，但对比实测情况，树木在风中的摆幅均有所减小，这可能是因为树木根部状态简化与忽略树叶等原因共同造成的。风致云杉剪切应力模型的破坏准则应当根据木材的材料特性决定，木材的各项强度指标（顺纹）从大到小排序为：抗拉强度抗弯强度抗压强度抗剪强度。因此，根据木材的剪切强度理论，模型的破坏准则采用最大剪切应力准则。在模拟的后处理中提取 不同风速下横纹、顺纹、截纹方向（模型中分别用、表示）剪切应力随时间变化的数值，如图 所示。木材横纹方向的抗剪强度较顺纹与截纹方向更小，不同风速下，横纹方向剪切应力在 个方向剪切应力中均为最大，故可将横纹抗剪强度与横纹方向剪切应力二者作对比，以此判断树木是否发生剪切破坏。通常木材的横纹抗剪强度约为顺纹径向与弦向抗剪强度平均值的一半。参考云杉的顺纹径面、弦面抗剪强度分别约为、，因此云杉的抗剪强度取 。图 不同风速下云杉在各个方向的剪切应力 由于木材标准强度测定试验选取的是经过标准化处理后的试块，其各项强度指标不同于立木状态的木
