非均布含水率梯度下圆木构件横截面湿度应力解析模型

Analytical Model of Moisture Stress in the Cross-Section of Wood Member Under Non-uniform Moisture Gradient

摘要

建立木构件湿度应力解析模型是研究木材干缩开裂规律及其对承载力影响的基础.首先通过木构件横截面应力平衡分析建立了含水率梯度下横截面应力分布模型；其次借鉴温度应力弹性解析模型和温—热类比关系,建立了非均布含水率梯度下木构件的平衡条件、物理条件、几何条件和协调方程，得到了木构件横截面切向、径向湿度应力的解析解、开裂临界含水率和切向拉压临界点；最后通过数值模拟方法验证了本文解析解的正确性.以杉木圆柱构件为例，分析了不同参数（含水率差、构件直径、初末含水率值）对木构件横截面径向、切向湿度应力分布的影响.结果表明：木构件横截面湿度应力分布与木构件初末含水率值及构件尺寸大小无关，而主要与含水率差和材料性能有关.

Abstract

Establishing the moisture stress analytical model of wood members is the basis for studying the laws of wood shrinkage cracking and their impact on the bearing capacity. Firstly， the stress distribution model of cross-section under moisture content gradient was established by analyzing the stress balance in the cross-section of wood member. Secondly， the equilibrium conditions， physical conditions， geometric conditions and coordination equations for wood members under non-uniform moisture content gradient were established based on the temperature stress elastic analytical model and the thermo-analogy relationship. The analytical solutions of tangential and radial moisture stresses， critical moisture content of cracking and critical point of tangential tension and compression were obtained by solving the equations. Finally， the correctness of the analytical solutions was verified by numerical simulations. Taking Chinese fir members as an example， the influence of different parameters （moisture content difference， member diameter， initial and final moisture content values） on the radial and tangential moisture stress distributions of the cross-section of wood members was analyzed. The results show that the moisture stress distributions in the cross-section of the wood member is not related to the initial and final moisture content of the wood member and the size of the member， but is related to the moisture content difference and material properties.

Graphical abstract

关键词

湿度应力 / 木构件横截面 / 含水率差 / 非均布含水率 / 解析模型

Key words

moisture stress / cross-section of wood member / moisture content difference / non-uniform moisture content / analytical model

由于外界气候条件的变化，木构件在正常使用过程中会相应地发生湿度传递，从而在构件截面上形成含水率梯度，进而产生湿度应力，当木材中的横纹应力达到木材抗拉强度时，木构件会发生开裂.越来越多的古木建筑保护者^［1-2］和现代木结构研究者^［3-4］意识到干缩裂缝对木结构的影响.

目前，少数国内外学者针对木材应力-应变关系进行了定性研究.例如，周宝华^［5］根据在弹性范围内木材的应力与应变成正比，对木材干燥过程中的应力变化进行了探究；Liu等^［6］研究表明：榆树木材干燥过程中，收缩应力逐渐增大，且在达到最大值后趋于稳定；程曦依等^［7］通过研究养生期间木材内部含水率、应变和应力的变化，得到木材养生处理可降低木材内含水率差和内部应力的结论；高璇^［8］推导了落叶松和樟子松温湿度应力计算公式；Arends等^［9］提出了一个基于扩散方程和线弹性的理论模型，以此解决木板单侧暴露于水分时所引起的应力分布问题；Yu等^［10］提出了一种木材在变湿度条件下随时间变化的材料响应改进模型.

还有少部分学者对非原木构件的湿度应力进行了一定研究.陈旭等^［11］采用有限元法分析了胶合木直梁和曲梁的横纹温度与湿度应力，结果表明：胶合木梁因环境温度产生的附加应力很小，而因环境湿度变化产生的附加应力很大；Gereke等^［12］研究了交叉层压实木面板的湿度应力，并进行了干燥试验，结果表明：可以通过施加预应力或控制中间层和外层的初始含水率来减少裂缝；Angst等^［13］通过有限元分析方法模拟了胶合木横截面湿度应力随时间的变化，结果表明，局部应力是否明显大于横截面高度方向上的平均应力，主要取决于胶合木的几何结构；Endo等^［14］研究了云杉在120 °C和不同加热湿度下平衡含水率的变化，得到结论：由于水分的存在加速了热降解过程，较高的加热湿度值会导致更大的质量损失，水热处理导致的吸湿性问题包括可逆和不可逆效应；Franke等^［15］通过实验室试验、现场测量和建模的方法研究了大型胶合木构件的含水率和湿度应力；Autengruber等^［16］提出一种仅根据周围气候的相对湿度来估计木构件边界处含水率的方法.

目前关于木构件湿度应力的研究方法基本采用试验或数值模拟，但这2种方法未能进一步揭示木构件应力的分布机理和应力间的变换关系，为了进一步定量分析湿度应力，需建立湿度应力解析解方程，只有将木材或木构件的湿度应力分布问题解决后，才能进一步研究木构件干缩裂缝预测模型.

1 木构件横截面湿度应力分布机理分析

木构件含水率若均匀下降Δw，径向干缩湿胀系数α_R、弦向干缩湿胀系数α_T，横截面从虚线位置均匀收缩到实线位置（图1a）.

取夹角为θ的扇形块，截面外边δ_R（R）=Rα_RΔw；径向收缩时弦向相应的变形值（扇形范围内虚线与实线差）δ{}_{\mathrm{T}}^{\mathrm{\text{'}}}（R） =Rθ-（R-δ_R）θ=δRθ=Rθα_RΔw；弦向收缩值δ_T（R）=Rθα_TΔw.

根据对称性，在通过圆心的任意纵向横截面上没有切应力；弦向外部受拉、内部受压，并且拉、压合力等于零；在ρ=ρ₀处，弦向应力为零；由于裂缝首先在外周出现，可以推断外周的弦向拉应力最大，见图1b.

当弦向存在正应力时，径向一定也存在正应力，才能满足径向平衡条件，见图1c.周边（ρ=R）径向压应力应为0；ρ越小，压应力越大（应力增量dσ_R>0）、压缩变形越大.径向受压变形将减小径向总变形（收缩变形和受力变形）与弦向收缩变形之间的差额，从而导致弦向拉应力减小.当ρ<ρ₀时，径向压应力增量dσ_R<0，弦向转为压应力.径向始终是压应力，其数值从外表的0逐渐增大.

在ρ=ρ₀位置，弦向收缩变形刚好等于径向收缩变形与压缩变形之和.当ρ<ρ₀时，径向收缩变形与压缩变形之和大于弦向收缩变形，弦向必须有压缩变形才能满足变形协调条件.

取图1c所示的微元，建立平衡微分方程和变形协调方程，再利用物理条件（应力-应变关系和干缩-含水率关系，两者都是线性关系）建立微元力和变形之间的关系.需要注意的是，ρ位置处的径向收缩变形可表示为δ_R（ρ）=ρα_RΔw（在含水率均匀分布的情况下），但压缩变形需要从0到ρ积分.弦向因轴对称，无论是收缩变形还是拉伸变形，都是均匀的.

2 湿度应力方程建立及求解

Jia等^［17］和战剑锋等^［18］在木材干燥研究中将含水率在截面上的分布近似为二次多项式，陈孔阳^［19］通过实验室木构件含水率试验测试，也将含水率在截面上的分布回归为二次多项式.本文假定含水率在截面上的分布为Δw（t，ρ）=Δw［2ρ/R-（ρ/R）²］.

2.1 平衡条件

根据微元径向力平衡条件（以拉应力为正）：

\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }=\frac{{\sigma }_{\mathrm{T}}-{\sigma }_{\mathrm{R}}}{\rho }，{\sigma }_{\mathrm{T}}=\rho \frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }+{\sigma }_{\mathrm{R}} .（1）

2.2 物理条件

对于正交各向异性材料，径向受力应变ε_R和弦向受力应变ε_T为

\left.\begin{array}{l}{\epsilon }_{\mathrm{R}}=\frac{{\sigma }_{\mathrm{R}}-{\nu }_{\mathrm{R}\mathrm{T}}{\sigma }_{\mathrm{T}}}{E_{\mathrm{R}}}=\frac{{\sigma }_{\mathrm{R}}}{E_{\mathrm{R}}}-\frac{{\nu }_{\mathrm{R}\mathrm{T}}}{E_{\mathrm{R}}}\left[\rho \frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }+{\sigma }_{\mathrm{R}}\right]\Rightarrow {\epsilon }_{\mathrm{R}}=\frac{(\mathrm{1}-{\nu }_{\mathrm{R}\mathrm{T}}){\sigma }_{\mathrm{R}}}{E_{\mathrm{R}}}-\frac{{\nu }_{\mathrm{R}\mathrm{T}}}{E_{\mathrm{R}}}\rho \frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }，\\ {\epsilon }_{\mathrm{T}}=\frac{{\sigma }_{\mathrm{T}}-{\nu }_{\mathrm{T}\mathrm{R}}{\sigma }_{\mathrm{R}}}{E_{\mathrm{T}}}=\frac{\rho }{E_{\mathrm{T}}}\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }+\frac{{\sigma }_{\mathrm{R}}}{E_{\mathrm{T}}}-\frac{{\nu }_{\mathrm{T}\mathrm{R}}{\sigma }_{\mathrm{R}}}{E_{\mathrm{T}}}\Rightarrow {\epsilon }_{\mathrm{T}}=\frac{\rho }{E_{\mathrm{T}}}\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }+\frac{{\sigma }_{\mathrm{R}}(\mathrm{1}-{\nu }_{\mathrm{T}\mathrm{R}})}{E_{\mathrm{T}}}.\end{array}\right\}（2）

其中：E_R，E_T为径向弹性模量和弦向弹性模量； ν_RT，ν_TR分别为径向-弦向泊松比和弦向-径向泊松比.

径向收缩变形 δ_R，sh（ρ）和弦向收缩变形δ_T，sh（ρ） 为

\left.\begin{array}{l}{\delta }_{\mathrm{R},\mathrm{s}\mathrm{h}}(\rho )=\frac{{\alpha }_{\mathrm{R}}\mathrm{\Delta }w}{R}{\int }_{\mathrm{0}}^{\rho }r\mathrm{d}r=(\mathrm{2}{\alpha }_{\mathrm{T}}-\mathrm{0.5}{\alpha }_{\mathrm{R}})E_{\mathrm{T}}\mathrm{\Delta }w{\rho }^{\mathrm{2}}/R-({\alpha }_{\mathrm{T}}-{\alpha }_{\mathrm{R}}/\mathrm{3})E_{\mathrm{T}}\mathrm{\Delta }w{\rho }^{\mathrm{3}}/R^{\mathrm{2}}，\\ {\delta }_{\mathrm{T},\mathrm{s}\mathrm{h}}(\rho )=\mathrm{d}\theta {\alpha }_T\mathrm{\Delta }w{\rho }^{\mathrm{2}}/R.\end{array}\right\}（3）

径向受力变形δ_R，F（ρ）和弦向受力变形δ_T，F（ρ）为

\left.\begin{array}{l}{\delta }_{\mathrm{R},\mathrm{F}}={\int }_{\mathrm{0}}^{\rho }{\epsilon }_{\mathrm{R}}(r)\mathrm{d}r\Rightarrow {\delta }_{\mathrm{R},\mathrm{F}}=\frac{\mathrm{1}}{E_{\mathrm{R}}}{\int }_{\mathrm{0}}^{\rho }{\sigma }_{\mathrm{R}}\mathrm{d}r-\frac{{\nu }_{\mathrm{R}\mathrm{T}}\rho {\sigma }_{\mathrm{R}}}{E_{\mathrm{R}}},\\ {\delta }_{\mathrm{T},\mathrm{F}}(\rho )=\rho \mathrm{d}\theta {\epsilon }_{\mathrm{T}}\Rightarrow {\delta }_{\mathrm{T},\mathrm{F}}=\mathrm{d}\theta [\frac{{\rho }^{\mathrm{2}}}{E_{\mathrm{T}}}\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }+\frac{\rho {\sigma }_{\mathrm{R}}(\mathrm{1}-{\nu }_{\mathrm{T}\mathrm{R}})}{E_{\mathrm{T}}}].\end{array}\right\}（4）

2.3 几何条件

位移协调条件：

{\delta }_{\mathrm{T}}=\mathrm{d}\theta {\delta }_{\mathrm{R}},{\delta }_{\mathrm{T},\mathrm{F}}(\rho )+{\delta }_{\mathrm{T},\mathrm{s}\mathrm{h}}(\rho )=\mathrm{d}\theta [{\delta }_{\mathrm{R},\mathrm{F}}(\rho )+{\delta }_{\mathrm{R},\mathrm{s}\mathrm{h}}(\rho )].（5）

2.4 协调方程（基本方程）

建立的协调方程见式（6），其中，α_E=E_T/E_R.

\begin{array}{l}\frac{{\rho }^{\mathrm{2}}}{E_{\mathrm{T}}}\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }+\frac{\rho {\sigma }_{\mathrm{R}}(\mathrm{1}-{\nu }_{\mathrm{T}\mathrm{R}})}{E_{\mathrm{T}}}+{\alpha }_{\mathrm{T}}\mathrm{\Delta }w{\rho }^{\mathrm{2}}/R=\\ \frac{\mathrm{1}}{E_{\mathrm{R}}}{\int }_{\mathrm{0}}^{\rho }{\sigma }_{\mathrm{R}}\mathrm{d}r-\frac{{\nu }_{\mathrm{R}\mathrm{T}}\rho {\sigma }_{\mathrm{R}}}{E_{\mathrm{R}}}+(\mathrm{2}{\alpha }_{\mathrm{T}}-\mathrm{0.5}{\alpha }_{\mathrm{R}})E_{\mathrm{T}}\mathrm{\Delta }w{\rho }^{\mathrm{2}}/R-({\alpha }_{\mathrm{T}}-{\alpha }_{\mathrm{R}}/\mathrm{3})E_{\mathrm{T}}\mathrm{\Delta }w{\rho }^{\mathrm{3}}/R^{\mathrm{2}}\Rightarrow \\ {\rho }^{\mathrm{2}}\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }+(\mathrm{1}-{\nu }_{\mathrm{T}\mathrm{R}})\rho {\sigma }_{\mathrm{R}}+{\alpha }_{\mathrm{E}}{\nu }_{\mathrm{R}\mathrm{T}}\rho {\sigma }_{\mathrm{R}}-\\ {\alpha }_{\mathrm{E}}{\int }_{\mathrm{0}}^{\rho }{\sigma }_{R}dr+\left(\mathrm{2}\right.{\alpha }_T-\mathrm{0.5}{\alpha }_R)E_T\Delta w{\rho }^{\mathrm{2}}/R-({\alpha }_T-{\alpha }_R/\mathrm{3})E_T\Delta w{\rho }^{\mathrm{3}}/R^{\mathrm{2}}=\mathrm{0}.（\mathrm{6}）\end{array}

2.5 方程求解

对基本方程两边关于ρ求导，设ρ=e ^t，则t=lnρ，有

\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}\rho }=\frac{\mathrm{1}}{\rho }\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}t}，\frac{{\mathrm{d}}^{\mathrm{2}}{\sigma }_{\mathrm{R}}}{\mathrm{d}{\rho }^{\mathrm{2}}}=\frac{\mathrm{1}}{{\rho }^{\mathrm{2}}}\frac{{\mathrm{d}}^{\mathrm{2}}{\sigma }_{\mathrm{R}}}{\mathrm{d}{t}^{\mathrm{2}}}-\frac{\mathrm{1}}{{\rho }^{\mathrm{2}}}\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}t}.（7）

微分方程转化为

\frac{{\mathrm{d}}^{\mathrm{2}}{\sigma }_{\mathrm{R}}}{\mathrm{d}{t}^{\mathrm{2}}}+\mathrm{2}\frac{\mathrm{d}{\sigma }_{\mathrm{R}}}{\mathrm{d}t}+{\sigma }_{\mathrm{R}}(\mathrm{1}-{\alpha }_{\mathrm{E}})=-K_{\mathrm{1}}{\mathrm{e}}^t/R+K_{\mathrm{2}}{\mathrm{e}}^{\mathrm{2}t}/R^{\mathrm{2}}，（8）

\left.\begin{array}{l}{K}_{\mathrm{1}}=E_{\mathrm{T}}(\mathrm{4}{\alpha }_{\mathrm{T}}-\mathrm{2}{\alpha }_{\mathrm{R}})\mathrm{\Delta }w,\\ K_{\mathrm{2}}=E_{\mathrm{T}}(\mathrm{3}{\alpha }_{\mathrm{T}}-{\alpha }_{\mathrm{R}})\mathrm{\Delta }w.\end{array}\right\}（9）

通过建立特征方程，求特解和通解，再建立边界条件：

\left\{\begin{array}{l}{\int }_{\mathrm{0}}^R{\sigma }_{\mathrm{T}}\mathrm{d}\rho =\mathrm{0},\\ {\sigma }_{\mathrm{R}}\left(R\right)=\mathrm{0}.\end{array}\right.\to \left\{\begin{array}{l}{C}_{\mathrm{2}}=\mathrm{0},\\ C_{\mathrm{1}}=\frac{K}{（\mathrm{4}-{\alpha }_{\mathrm{E}}）}{R}^{(\mathrm{1}-\sqrt[]{{\alpha }_{\mathrm{E}}})}.\end{array}\right.（10）

则在抛物线型含水率条件下圆盘的应力分布为

\left.\begin{array}{l}{\sigma }_{\mathrm{R}}=\frac{E_{\mathrm{T}}(\mathrm{4}{\alpha }_{\mathrm{T}}-\mathrm{2}{\alpha }_{\mathrm{R}})\mathrm{\Delta }w}{（\mathrm{4}-{\alpha }_{\mathrm{E}}）}\left[\frac{\rho }{R}-{(\frac{\rho }{R})}^{(\sqrt[]{{\alpha }_{\mathrm{E}}}-\mathrm{1})}\right]+\frac{E_{\mathrm{T}}(\mathrm{3}{\alpha }_{\mathrm{T}}-{\alpha }_{\mathrm{R}})\mathrm{\Delta }w}{（\mathrm{9}-{\alpha }_{\mathrm{E}}）}\left[{(\frac{\rho }{R})}^{\mathrm{2}}-{(\frac{\rho }{R})}^{(\sqrt[]{{\alpha }_{\mathrm{E}}}-\mathrm{1})}\right],\\ {\sigma }_{\mathrm{T}}=\frac{E_{\mathrm{T}}(\mathrm{4}{\alpha }_{\mathrm{T}}-\mathrm{2}{\alpha }_{\mathrm{R}})\mathrm{\Delta }w}{（\mathrm{4}-{\alpha }_{\mathrm{E}}）}\left[\mathrm{2}\frac{\rho }{R}-\sqrt[]{{\alpha }_{\mathrm{E}}}{(\frac{\rho }{R})}^{(\sqrt[]{{\alpha }_{\mathrm{E}}}-\mathrm{1}）}\right]-\frac{E_{\mathrm{T}}(\mathrm{3}{\alpha }_{\mathrm{T}}-{\alpha }_{\mathrm{R}})\mathrm{\Delta }w}{（\mathrm{9}-{\alpha }_{\mathrm{E}}）}\left[\mathrm{3}{(\frac{\rho }{R})}^{\mathrm{2}}-\sqrt[]{{\alpha }_{\mathrm{E}}}{(\frac{\rho }{R})}^{(\sqrt[]{{\alpha }_{\mathrm{E}}}-\mathrm{1})}\right].\end{array}\right\}（11）

令σ_T（R）等于f_T，可得到临界含水率的变化量Δw_cr，则

\left.\begin{array}{l}{\rho }_{\mathrm{0}}/R=\mathrm{e}\mathrm{x}\mathrm{p}\left(\frac{\mathrm{l}\mathrm{n}\sqrt[]{{\alpha }_{\mathrm{E}}}/\mathrm{2}}{(\mathrm{2}-\sqrt[]{{\alpha }_{\mathrm{E}}})}\right)，\\ \mathrm{\Delta }{w}_{\mathrm{c}\mathrm{r}}=-\frac{(\mathrm{2}+\sqrt[]{{\alpha }_{\mathrm{E}}})(\mathrm{1}-{\nu }_{\mathrm{R}\mathrm{T}}{\nu }_{\mathrm{T}\mathrm{R}})f_{\mathrm{T}}}{E_{\mathrm{T}}(\mathrm{2}{\alpha }_{\mathrm{T}}-{\alpha }_{\mathrm{R}})}.\end{array}\right\}（12）

从式（11）可以看出：弦向应力等于零的位置与Δw无关.在ρ₀处，径向收缩变形与受压变形刚好等于弦向收缩变形；r>ρ₀部分，弦向需有拉应变才能满足位移协调；r<ρ₀部分，径向收缩变形与受压变形超过弦向收缩变形，弦向需有压应变才能满足位移协调.

从式（12）可以看出：径向拉压应力分界点和构件直径的比值是一定值，该定值与材料属性相关，而与含水率无关；临界含水率变化值Δw_cr与构件直径无关，与材料属性有关，随着外界湿度持续下降，开裂会延续.

3 数值算例

3.1 模型建立

从上述湿度应力解析解模型中可以看出，应力分布、开裂条件和材料性能、含水率分布、含水率差有关.为了验证上述圆盘和圆柱类模型的正确性，选取模型中几种代表性工况进行分析，具体参数如表1所示.

通过大型有限元分析软件ABAQUS进行数值模拟，建立了表1所示工况的有限元模型.有限元建模方法借鉴参考文献［20-21］，有限元模型的主要材料参数见参考文献［20-25］，如表2所示.

值得注意的是，在ABAQUS的温度场中，有限元单元选用三维八节点热传导DC3D8单元，网格尺寸分别设为5，2.5和1 mm，综合考虑计算精度和计算速度，选用网格尺寸约为2.5 mm.

3.2 有限元分析结果

建立有限元模型时，1，2，3分别定义为木构件的径向、切向和长向.将上述有限元模型的计算结果列于图2，图中左边为S22（切向）应力云图，右边为S11（径向）应力云图.

3.3 对比结果

将上述模型的径向应力（S11）和切向应力（S22）计算结果绘于图3，分别分析对比圆盘尺寸和初始含水率（图3a）、含水率分布形式（图3b、非线性含水率分布形式（图3c）和圆盘-圆柱构件形式（图3d）等因素对湿度应力的影响.

从图3可以看出，对比不同含水率差、构件尺寸、初末含水率下圆盘的弦向应力和切向应力，理论模型和有限元计算结果基本吻合，说明非均布含水率梯度下木构件横截面湿度应力解析解理论模型和有限元模型是正确的.最大拉应力均发生在外侧，对比结果如下：

1） 不同含水率差情况下，含水率梯度越大，其外侧应力越大；

2） 不同尺寸的圆盘构件，湿度应力分布和木构件尺寸大小无关，和含水率差、相对圆心距离和材料性能有关；

3）圆盘的弦向应力和切向应力与构件初始含水率及最终含水率的数值无关，而与含水率差值有关.

4 不同含水率分布形式下应力对比

基于本文建立的非均布含水率分布下湿度切向和径向应力解析解（式（11）），对均匀含水率分布和线性含水率分布下的木材湿度应力进行分析.在含水率差均为15%的情况下，分别将3种含水率分布下的切向应力和径向应力列于图4.

从图4可以看出，圆木构件弦向应力沿径向呈现出“内部为压力，外部为拉力”的特点；在相同含水率梯度差15%的情况下，构件边缘切向应力在线性含水率分布下可达到11.93 MPa，而在均匀含水率分布下切向应力仅为5.87 MPa；径向应力始终为压力，在圆木构件边缘为0.

5 结 论

1） 本文通过木构件横截面应力平衡分析、温度应力弹性解析模型和温-热类比关系，建立了非均布含水率梯度下木构件切向和径向解析解、开裂临界含水率和切向拉压临界点；最后通过数值模拟方法验证了本文解析解的正确性.

2） 径向拉压应力分界点和构件直径的比值是一定值，与材料属性相关，和含水率无关；临界含水率变化值Δw_cr与构件直径无关，仅与材料属性有关；

3） 木构件湿度应力分布和木构件尺寸大小、初末含水率值无关，和含水率差、相对圆心距离和材料性能有关；

4） 木构件湿度应力分析可借鉴温度应力弹性解析模型、温-热类比关系，大型有限元分析软件ABAQUS可通过传热模块模拟传湿过程的方法是可行的.

5） 在相同含水率梯度差下，构件边缘切向应力在线性含水率分布时最大，而在均匀含水率分布时最小.

参考文献

[1] 袁霄， 陈勇平， 郭文静. 古建筑木结构多发残损的特征及防控分析［J］. 木材科学与技术， 2021， 35（5）： 54-59.

[2] Yuan Xiao， Chen Yong-ping， Guo Wen-jing. Features and prevention of common damages in ancient timber structures［J］. Chinese Journal of Wood Science and Technology， 2021， 35（5）： 54-59.

[3] 李钊，王志涛，郭小东. 残损古建筑木结构力学性能相关研究进展与展望［J］. 林产工业， 2022， 59（12）， 39-46.

[4] Li Zhao， Wang Zhi-tao， Guo Xiao-dong. Research progress and prospect on mechanical properties of damaged ancient timber structures［J］. China Forest Products Industry， 2022， 59（12）： 39-46.

[5] 陈偲，沈肇雨，王正， 等. 正交胶合木平面剪切的开裂形貌及其破坏模式探究［J］. 林产工业， 2022， 59（6），7-13，24.

[6] Chen Si， Shen Zhao-yu， Wang Zheng， et al. Study on plane shear cracking morphology and failure mechanism of cross laminated timber［J］. China Forest Products Industry， 2022， 59（6）： 7-13，24.

[7] Vida C， Lukacevic M， Hochreiner G， et al. Size effect on bending strength of glued laminated timber predicted by a numerical simulation concept including discrete cracking［J］. Materials & Design， 2023， 225： 111550.

[8] 周宝华. 木材干燥过程中内应力的初步研究［J］. 南京林业大学学报（自然科学版）， 1982， 6（2）： 76-90.

[9] Zhou Bao-hua. A preliminary study of the internal stresses of timber［J］. Journal of Nanjing Forestry University， 1982， 6（2）： 76-90.

[10] Liu J X， Wang X M. Effect of drying temperature and relative humidity on contraction stress in wood［J］. BioResources， 2016， 11（3）： 6625-6638.

[11] 程曦依， 李贤军， 黄琼涛， 等. 木材干燥后养生期间含水率及应力变化特点［J］. 林业工程学报， 2016， 1（2）： 38-44.

[12] Cheng Xi-yi， Li Xian-jun， Huang Qiong-tao， et al. The study of moisture content and inner stress variation of kiln-dried wood during equalization treatment［J］. Journal of Forestry Engineering， 2016， 1（2）： 38-44.

[13] 高璇. 湿度循环作用下木材的横纹应力疲劳研究［D］. 哈尔滨： 哈尔滨工业大学， 2016.

[14] Gao Xuan. Study on transverse stress fatigue of wood under humidity cycle［D］. Harbin： Harbin Institute of Technology， 2016.

[15] Arends T， Pel L， Huinink H P. Hygromorphic response dynamics of oak： towards accelerated material characterization［J］. Materials and Structures， 2017， 50（3）： 181.

[16] Yu T Y， Khaloian A， van de Kuilen J W. An improved model for the time-dependent material response of wood under mechanical loading and varying humidity conditions［J］. Engineering Structures， 2022， 259： 114116.

[17] 陈旭， 祝恩淳. 胶合木梁中温度与湿度应力研究［J］. 工业建筑， 2008， 38（sup1）： 862-866.

[18] Chen Xu， Zhu En-chun. Moisture and thermal stresses in glulam beams［J］. Industrial Construction， 2008， 38（sup1）： 862-866.

[19] Gereke T， Niemz P. Moisture-induced stresses in spruce cross-laminates［J］. Engineering Structures， 2010， 32（2）： 600-606.

[20] Angst V， Malo K A. Moisture-induced stresses in glulam cross sections during wetting exposures［J］. Wood Science and Technology， 2013， 47（2）： 227-241.

[21] Endo K， Obataya E， Zeniya N， et al. Effects of heating humidity on the physical properties of hydrothermally treated spruce wood［J］. Wood Science and Technology， 2016， 50（6）： 1161-1179.

[22] Franke B， Franke S， Schiere M， et al. Moisture content and moisture-induced stresses of large glulam members： laboratory tests， in situ measurements and modelling［J］. Wood Material Science & Engineering， 2019， 14（4）： 243-252.

[23] Autengruber M， Lukacevic M， Gröstlinger C， et al. Numerical assessment of wood moisture content-based assignments to service classes in EC5 and a prediction concept for moisture-induced stresses solely using relative humidity data［J］. Engineering Structures， 2021， 245： 112849.

[24] Jia D H， Afzal M T. Modeling of moisture diffusion in microwave drying of hardwood［J］. Drying Technology， 2007， 25（3）： 449-454.

[25] 战剑锋， 顾继友， 艾沐野. 白桦干燥过程的横纹干燥应力［J］. 东北林业大学学报， 2005， 33（4）： 25-28.

[26] Zhan Jian-feng， Gu Ji-you， Ai Mu-ye. Transverse drying stress of white birch wood during drying［J］. Journal of Northeast Forestry University， 2005， 33（4）： 25-28.

[27] 陈孔阳. 木构件腐朽和干缩裂缝的研究［D］. 南京： 东南大学， 2019.

[28] Chen Kong-yang. Study on decay and shrinkage cracks of wood components［D］. Nanjing： Southeast University， 2019.

[29] Chen K Y， Qiu H X， Sun M L， et al. Experimental and numerical study of moisture distribution and shrinkage crack propagation in cross section of timber members［J］. Construction and Building Materials， 2019， 221： 219-231.

[30] Autengruber M， Lukacevic M， Gröstlinger C， et al. Finite-element-based prediction of moisture-induced crack patterns for cross sections of solid wood and glued laminated timber exposed to a realistic climate condition［J］. Construction and Building Materials， 2021， 271： 121775.

[31] 周鸿屹. 木结构钢填板螺栓连接节点力学性能试验研究与分析［D］. 哈尔滨： 哈尔滨工业大学， 2017.

[32] Zhou Hong-yi. Experimental study and analysis on mechanical properties of bolted joints with steel filler plates in wood structures［D］. Harbin： Harbin Institute of Technology， 2017.

[33] Larsen F， Ormarsson S. Numerical and experimental study of moisture-induced stress and strain field developments in timber logs［J］. Wood Science and Technology， 2013， 47（4）： 837-885.

[34] 何朝红. 胶合木—钢夹板螺栓连接力学性能试验研究［D］. 长沙： 中南林业科技大学， 2017.

[35] He Zhao-hong. Experimental study on mechanical properties of bolted connection between glued wood and steel splint［D］. Changsha： Central South University of Forestry & Technology， 2017.

[36] 王娟，许刃文，杨庆山， 等.基于能量分析的唐代殿堂型木构架抗震机理研究［J］.工程力学，2022，39（11）： 73-88.

[37] Wang Juan， Xu Ren-wen， Yang Qing-shan， et al. Energy-based Seismic Performance Analysis of Palace-style Timber Frame of Tang-dynasty［J］. Engineering Mechanics， 2022，39（11）： 73-88.

基金资助

新疆维吾尔自治区自然科学基金资助项目(2022D01C687)