Analysis of the Structure and Hydraulic Function of Bordered Pits Using the Lattice Boltzman Method

Abstract

Fluid flow between adjacent tracheids is realized through bordered pits in the xylem of conifers. The pit has an extremely small size and a highly complex structure. This paper presents a mesoscopic analytical method for the relationship between the pit structure and its hydraulic characteristics through mathematical modeling using the lattice Boltzmann method (LBM) and curved boundary treatment. Mongolian Scots pine were selected as the research subject of this study, and the bordered pit structure parameters was collected by scanning electron microscopy (SEM) and transmission electron microscopy (TEM), and the original geometric features were maintained for direct modeling analysis. The model revealed the relationship between various components of the bordered pit and liquid flow velocity/resistance, indicating that margo is the main factor affecting flow resistance. Further anatomical investigation separately analyzed the influence of change in a single factor, including pit diameter, pit aperture diameter, pit depth, torus diameter, and margo thickness, on the overall flow and pressure drop to confirm the importance of various factors in this relationship. Additionally, the influence of pore size and pore location distribution in the margo on the flow rate and pressure drop was further analyzed quantitatively. The results showed that the flow rate through individual pores is the result of the combined effect of pore area and radial position of the pore in the margo. Our study promotes the research and application of the mesoscopic model LBM in simulating flow conditions in the complex flow field of pits, which realizes the numerical analysis of the flow field model based on individualized real bordered pits. In comparison with the classical macroscopic model, the accuracy and effectiveness of the proposed model are proved. This research can provide a promising method for analyzing the physiological and ecological functions of conifer and realizing the efficient utilization of wood resources.

1. Introduction

Bordered pits play an important role in conducting water through the xylem, and act as a bridge to transport water or nutrients between tracheids [1,2]. Bordered pits are a part of the primary wall but do not undergo thickening during secondary wall formation; therefore, bordered pits appear as a depression in the secondary cell wall [3]. Bordered pits of Mongolian Scots pine, a coniferous tree species, were selected as the research subject of this study (Figure 1). The pit border (Figure 1e,f) is located above the pit membrane, which is arched, and the pit membrane separates adjacent pits. A primary thickening, known as torus, is formed in the middle of the pit membrane on the axial tracheid wall (Figure 1c,f), which is generally considered to be non-water-conducting. The diameter of the torus is usually larger than that of the pit aperture, and it acts as a piston to regulate the intertrachedral fluid flow. The pit membrane contains a complex network composed of a large number of microfibrils interweaved around the torus, which is known as the margo (Figure 1c). Most of the microfibrils are radially arranged, while some are arranged in a tangential or oblique direction. The sap flow between adjacent tracheids (Figure 1d,e) in coniferous wood xylem is mainly realized by the pores on the margo [4,5]. These pores enable the liquid to pass through the tracheids in a flexible manner, which is essential for preventing the spread of viral pathogens in tracheids and for preventing embolism in the whole xylem as well as for promoting tree growth [6,7].

Because of the extremely small size and complex structure of the bordered pits, it is extremely challenging to analyze the mechanism of flow through the pits at a microscopic level using physiological tools alone. The early estimation of individual pit resistance was achieved using cellulases, which dissolve the pit membranes; comparing the resistance before and after cellulase treatment revealed a relatively low range of pit resistance (1.0–28.8 MPa s m⁻¹) [8,9,10]. However, cellulase does not affect the pit membranes in some tree species [11,12]. Therefore, some scholars developed a physical model of greatly scaled-up bordered pits to study the hydrodynamic resistance of bordered pits in plant stems [12]. Using this model, the flow resistance was estimated at 1.7 × 10¹⁵ Pa.sm⁻³. This model reveals that the flow resistance of pit varies with the pit shape. However, bordered pits are highly complex and irregular in shape; therefore, determining the flow resistance based on the pit shape is challenging. Many scholars proposed that modeling may be a more feasible method for small structures such as bordered pits. They modeled xylem as a water-conducting medium and estimated the liquid flow inside the plant, revealing that pit membranes can account for 50% or more of the total hydraulic resistance in the xylem [13]. Valli et al. defined the pit margo region as a porous medium in the simulation process, and analyzed the flow in the pit of conifer xylem, showing that 38% of the flow resistance was attributed to the margo [14]. Schultes solved the Navier–Stokes equation, based on the geometry of bordered pits in black spruce, and showed that the margo and torus together constitute more than 80% of the flow resistance [15]. To a certain extent, these studies have achieved a qualitative analysis of many features in the xylem, which are otherwise difficult to observe or measure. However, most models cannot truly estimate the actual morphology and structural characteristics of the overall pit. Therefore, there is a lack of quantitative analysis of the relationship between the structural parameters of bordered pits and the fluid flow resistance characteristics, and no systematic study has been conducted to determine the influence of pit structure on water transmission. In addition, some methods require constant adjustment of the grid density, which leads to a relatively complicated calculation process.

Based on these considerations, the mesoscopic dynamic lattice Boltzman method (LBM) and curve boundary treatment are proposed to simulate the relationship between complex boundary of the pit structures and flow characteristics. The curved boundary treatment is modified according to the characteristics of the bordered pit structure. Given that LBM has the advantage of not only the few assumptions of the microscopic method but also the macroscopic method, which does not require details of molecular motion, it has more advantages and potential than the traditional numerical method in dealing with the description of complex flow phenomena such as multi-scale, multiphase system, interface dynamics, and seepage, among others. When analyzing the pit structure and hydraulic characteristics, including the influencing factors of fluid migration process under its internal complex structure, the LBM-C model can be directly analyzed with the original geometric characteristics of the pit. Compared with the traditional method, this method shows great potential in the efficient realization of boundary conditions.

2. Materials and Methods

Fresh material from small twigs and branches of Mongolian Scots pine (2–3 years) were collected during September and October from the forest farm of Northeast Forestry University. The wood samples were cut into 1-cm long sections and initially split in half. To prevent the deviation of torus caused by the relatively high interfacial tension between air and water during drying, organic solvents with low surface tension (alcohol and acetone) were used to replace the water with high surface tension, so as to keep the pit state of raw wood as much as possible. The wood sample was split along the radial direction, and the sapwood part was selected to collect the bordered pit image. It can be seen from the TEM and SEM images that the curvature of the interface at the edge of the pit borders changed greatly; pores formed by the irregular interweaving of microfibrils on the membrane have different shapes and sizes. The boundary conditions for the passage of fluid are very complicated, and the boundary movement causes the extrusion and expansion of the fluid area, and the space between the borders of two adjacent pits is an important feature that affects the distribution of water flowing through the margo. In addition, liquids in plants flow with a low Reynolds number (mainly adhesion) [16]. If the flow space is small, pores will generate high shear force and limit the liquid flow. It can be seen that the space of two adjacent pits borders is directly related to the inter-facial curvature of the pit borders. If the boundary of the analysis object is handled improperly, the simulation results will not be able to support the investigation of the mechanism and factors influencing the water transport through bordered pits.

Based on these considerations, the application of the curved boundary treatment for the fluid flow analysis of the complex pit boundary conditions is proposed in this paper. In the LBM model, the function of boundary conditions is to transfer the information of various physical properties of the boundary to the fluid. This transmission is completed by solving the discrete Boltzmann equation of particle distribution on a regular lattice. Particles first propagate to adjacent lattice points at each time step and then interact through local collisions, where the momentum is redistributed. Hydrodynamic variables such as velocity u are obtained from the velocity moment of the distribution function.

2.1. Lattice Boltzmann Method

LBM describes the motion of a fluid particle distribution function with discrete velocity on a fixed lattice. According to non-equilibrium statistical mechanics, the Boltzmann equation describing the contribution of collisions to the distribution function can be written as:

(1)$f_i(x+c_i\Delta t,t+\Delta t)-f_i(x,t)=-\frac{1}{\tau }\left(f_i(x,t)-f_i{}^{eq}(x,t)\right)$

where x is a lattice point on the computational lattice L, τ is the dimensionless relaxation time, f_i(x,t) is the density of fluid particles moving in the c_i direction, and c_i is the vector pointing to the neighboring lattice node. The LBM used in this study is the 9-velocity incompressible lattice D2Q9 model for the 2D study [17] and the 19-velocity incompressible lattice D3Q19 model for the 3D study [18,19]. The D2Q9 speed configuration can be written as follows:

(2)${\mathrm{c}}_i=c\{\begin{array}{ll}\left(0,0\right)& i=0\\ \left(\pm 1,0\right)\left(0,\pm 1\right)& i=1,2,3,4\\ \left(\pm 1,\pm 1\right)& i=5,6,7,8\end{array}$

Moreover, for the D3Q19 lattice:

(3)${\mathrm{c}}_i=c\{\begin{array}{ll}\left(0,0,0\right),& i=0\\ \left(\pm 1,0,0\right),\left(0,\pm 1,0\right),\left(0,0,\pm 1\right)& i=1,2,\dots ,6\\ \left(\pm 1,\pm 1,0\right),\left(\pm 1,0,\pm 1\right),\left(0,\pm 1,\pm 1\right)& i=7,8,\dots ,18\end{array}$

Lastly, f_i^eq is the equilibrium distribution that can be selected in different ways, and the common choice can be written as follows:

(4)${f_i}^{eq}={\xi }_i\rho \left(1+\frac{1}{c_s^2}(c_i\cdot u(x,t)+\frac{1}{2c_s^4}{\left((c_i\cdot u(x,t)\right)}^2-\frac{1}{2c_s^2}{u}^2(x,t)\right)$

where ξ_i is a weighting factor, which depends on the length of the connection vector $c_i$ in our simulation. For the D2Q9 model, ξ₀ = 4/9; i is the nearest neighbor lattice displacement ξ_1–4 = 4/9 and the diagonal ξ_5–8 = 1/36. For the D3Q19 model, ξ₀ = 1/3; ξ_1–6 = 1/18; ξ_7–18 = 1/36. Furthermore, c_s is the sound velocity of this model, which is calculated as c_s = c/$\sqrt{3}$, where c is the flow velocity, which is calculated as c = ∆x⁄∆t, where ∆x and ∆t are the lattice step and discrete time step, respectively. The kinematic viscosity of the simulated fluid (v) can be calculated as follows:

(5)$v=c^2\left(2\tau -1\right)\Delta t/6$

2.2. Curved Boundary Treatment for LBM of Micro-Scale Liquid Flow in Pits

Many authors have explored the curve boundary processing condition [20,21,22], which is used to modify the theoretical results according to the structural characteristics of the pit. The lattice diagram of the flow field at the pit boundary using the curved boundary treatment is shown in Figure 2. The red solid line refers to the physical curved wall between the lattice points, and squares filled with crisscross lines represent individual lattice points. If the square intersects the boundary, the point is considered a boundary point.

However, if the square does not intersect the boundary, the fluid lattice points on one side are represented by hollow circles, and the wall points on the other side are represented by filled diamonds. Then, the lattice points of the solid and fluid regions are represented by x_w and x_f, respectively, as their position vectors. In the fluid evolution process, the distribution function f_i (x + t) after the collision at the x_f lattice point is known. The solid point of the green filled circle (x_b) is the intersection of the physical boundary and the line connecting x_w and x_f. To reflect the distance, the location of the physical boundary line is indicated by the distance ratio between the fluid point and the boundary lattice point:

(6)$\lambda =\left|x_f-x_w\right|/\left|x_f-x_b\right|(0\le \lambda <1)$

The boundary condition denotes the distribution function from the boundary point to the fluid point. The distribution function is divided into the equilibrium and non-equilibrium parts, similar to the c_i direction in the Figure 2:

(7)$f_i(x)=f_i{}^{eq}(x)+f_i{}^{neq}(x)$

The non-equilibrium part of the w-point distribution function is extrapolated from that of the two-point distribution function of f and ff as follows:

(8)$f_i{}^{neq}(x_b)=2f_i{}^{neq}(x_f)-f_i{}^{neq}(x_{ff})$

The density of point w is extrapolated from that of f and ff as follows:

(9)${\rho }_w=2{\rho }_f-{\rho }_{ff}$

The speed of point w is extrapolated from that of f and n as follows:

(10)$v_w=\frac{v_n+(\lambda -1)v_f}{\lambda },0.5\le \lambda <1.5$

Then, f_i^eq(x) is calculated by substituting Equation (10) into the equilibrium distribution expression (Equation (4)).

This boundary condition is used to simulate the liquid flow in the pits. The width of the corresponding position of the pit border is ω(x) and the pressure at x is p(x). Assuming that the pressure p(x) has a linear relationship with the width ω(x), Equation (11) is established:

(11)$p(x)-p_0=\sigma \left(\omega (x)-{\omega }_0\right)$

where ω₀(x) is the width of the pressure inlet, i.e., the width of the pit channel formed by the tip of the pit border, and σ is almost constant. Because the Reynolds number (Re) of the fluid in the xylem is small (Re << 1), the flow is stable, and the flow state is generally laminar [15], which can be approximated as a poiseuille flow with radial changes, and the flow velocity along the pipe direction at (x, y) is calculated as:

(12)$u(x,y)=u_0(x)\left(1-\frac{4y^2}{{\omega }^2\left(x\right)}\right)$

where ω₀(x) is the center linear velocity. Assuming that the entrance velocity is a parabola, its peak velocity (entry center) is calculated as:

(13)$u_0(x)=\frac{{\omega }^2\left(x\right)}{4\eta }\frac{\partial p(x)}{\partial x}$

The pressure gradient here is arbitrarily selected as a typical biological characteristic. The selected value is 0.02–0.03 MPa m⁻¹ [23], and the peak velocity obtained using Equation (12) is between 0.001 and 0.003 ms⁻¹. η is the viscosity coefficient. In this study, assuming that the inner wall is a non-slip boundary condition, the definition of flow Q can be calculated as:

(14)$Q(x)={\displaystyle {\int }_{-R}^{R}u(x,y)}dy=\frac{2}{3}{u}_0(x)\omega (x)=\frac{{\omega }^3\left(x\right)}{6\eta }\frac{\partial p(x)}{\partial x}$

Based on Equation (11):

(15)$\frac{\partial p(x)}{\partial x}=\frac{\sigma d\omega (x)}{dx}$

Consequently:

(16)$Q(x)=\frac{\sigma {\omega }^3\left(x\right)}{6\eta }\frac{d\omega (x)}{dx}$

Integrating this equation, for the steady-state flow through the pits, the relationship among the boundary of the inner curve of the pit, the pressure, and x is expressed as:

(17)${\omega }^4\left(x\right)-{\omega }^4\left(0\right)=Bx,B=-48vQ/\lambda$

(18)${\left({\omega }_0+\frac{p\left(x\right)-p_0}{\lambda }\right)}^4-{\left(a_0+\frac{p(0)-p_0}{\lambda }\right)}^4=Bx$

2.3. Materials and Numerical Parameters

Structural parameters on the bordered pits of Mongolian Scots pine branches were collected via SEM and TEM. The detailed structure of the pit membrane obtained by SEM is shown in Figure 3. Intact and undamaged images of the margo were used for subsequent analysis. However, in the actual image acquisition process, margo microfibrils were partially damaged in many instances because of the pre-processing technology, cutting angle, cutting position, and other factors. More than 87% of the 52 pit structures photographed showed varying degrees of damage. To utilize the damaged part of the margo microfibril, the intact area was copied or repaired according to the microfibril direction connection, and the original microfibril arrangement form was retained.

The overall structure and size of the pit, including pit border and pit chamber, were obtained by TEM (Figure 4a). It should be emphasized that the cross-sectional view, based on TEM, showing the structural shapes and dimensions of the pit border, and the front view of the pit membrane captured by SEM could not be obtained from the same pit. However, in a large number of images, the shape and relevant size of the pit border showed obvious changes. Therefore, it was necessary to integrate the information of multiple groups of pit images at similar positions to establish a consensus model, and then to build a three-dimensional model as shown in Figure 4c. Then, image processing was used to perform noise reduction, edge detection, segmentation, and measurement of the image to obtain the structure and size of each component of the pit, as shown in Figure 4b.

The model was discretized with several grid spacings to obtain a sufficient lattice size. On the real wall of traicheid, as shown in Figure 4a, the pits are very close to each other. Therefore, to ensure that the pressure drop calculation results are sufficiently accurate, the simulation volume chosen in this study was only slightly larger than the size of the pit chamber [14]. As shown in Figure 4b, the inner diameter of the pit is 10.23 μm and the depth is 3.38 μm. In the simulation, the pit is discretized as 600 × 600 lattice sites, with a channel width comprising 276 sites. Initial density ρ = 1.0, pressure boundary condition is implemented at the inlet and outlet. Inlet density ρ_in = 1.001, outlet density ρ_out = 0.999. The flow direction is perpendicular to the pit membrane, the relaxation parameter is set as τ = 1.0, the periodic boundary is applied in all directions, and the density and momentum of the fluid are averaged at the inlet and outlet surfaces of the system.

3. Results and Discussion

3.1. Simulations of the Water Flow through Different Components of the Bordered Pits

Figure 5 shows the flow field in the pit of Mongolian Scots pine, and the color indicates the velocity distribution value. To analyze the influence of the structure of the pit components, such as border, torus, and margo, on the flow field in the pit, only the channel was established in the initial calculation model without the pit border and pit membrane (Figure 5a). The curvature of these velocity contours passing through the channel was small, indicating that the pressure drop of the channel and the flow resistance were very small. In Figure 5b, the pit border has been added to the model, focusing on the analysis of the influence of the pit border structure on the liquid flow. The flow at the narrow pit aperture formed by the pit borders was greatly increased because of the reduction of the cross-section, and the flow velocity was also increased. In this process, without the participation of the pit membrane, the water in the pit chamber flows only in the horizontal basin of the aperture size, and most of the water between the adjacent pits on the same side of the channel does not participate in the flow.

Then, the pit torus was added to the model (Figure 5c). In this model, when the fluid flows through the space near the torus, under the perturbation of the pit border and torus, the water flow will separate and there may be turbulence locally. This is because of the flow channel shrinkage and changes in the flow area increase the velocity gradient, resulting in a greater strain rate. Further analysis showed that most of the flow occurs in the area close to the torus, and the inner boundary of the pit border is simulated according to the real structure of the pit itself. Compared with the boundary conditions of Kvist et al. [24] and Schulte [15], which were approximately a straight line, the flow space formed by the curved boundary is relatively large, which has an important influence on the distribution of water flow. By integrating velocity in this area, more than 98% of the flow rate was generated through the 40% area inside the adjacent pit chambers close to the torus. The model was developed further by adding the part of the margo connected between the pit border and the torus (Figure 5d), which showed that the flow almost filled the narrow space between the adjacent pit borders, and the path dispersion also started to change as the circulation area decreased, increasing the tortuosity.

Assuming that the area between pit borders, torus, and margo pores are connected in series, this method roughly estimates the proportion that each component of the pit occupies in the total flow resistance of the pit. The analysis showed that the pressure drop in the four models increased gradually, and the effect of each component of the pit on the flow resistance was calculated separately. The loss along the way of the no-border model was very small or even negligible. The pressure loss caused by the pit border can be directly calculated from the inlet and outlet pressure drop, and the pressure loss caused by the torus can be obtained by subtracting the pressure drop of the pit border model from that of the no-border model. The pressure loss caused by the margo was obtained by subtracting the pressure drop of the model without margo from that of the complete model. The resistance was obtained by dividing the model pressure loss by the flow rate through the model (3.28 × 10⁻¹⁴ m³ s⁻¹). The pressure drop of each component was divided by the total pressure drop of the pit to obtain the proportion of the total flow resistance in the pit (Table 1).

3.2. Influence of Pit Structure on Flow Rate and Pressure Drop

Analysis of several groups of bordered pit samples by TEM and SEM revealed variability in the size of the pit. Under normal circumstances, liquids passed through large pores easily [6]. The diameter of the torus was slightly larger than half of the diameter of the pit membrane, but it varied greatly from less than half to two-third of the diameter of the pit membrane. In addition, the thickness of microfibril interweaving and the coverage of amorphous materials also varied, resulting in variable margo thickness. Within the water-conducting network, the hydraulic resistance was associated with pit membrane thickness [25]. In order to explore the relationship between the pit structure and the internal flow field, the single factor analysis method is used. The size of each part of a pit structure is adjusted and calculated separately, so as to analyze the influence of the change of this factor on the flow rate and pressure drop. The size change is based on the size range of 52 groups of pit images obtained by SEM and TEM. The total pressure drop and flow rate are calculated respectively, and the influence of the size of each part of the bordered pit on the flow characteristics is analyzed (Figure 6).

Simulation results of the flow field in the pits with different diameters are shown in Figure 6a. After changing the diameter of the pit, the dimensions of each part of the pit were scaled proportionally. The calculation results showed that as the pit diameter increased, the flow rate through the pits also increased, whereas the total pressure drop decreased significantly. This is because the flow velocity was unchanged. Thus, when the pit diameter increases, the flow area of the liquid through the pit aperture and pit membrane increases, resulting in increased flow rate and decreased pressure drop.

Figure 6b shows the effect of different aperture diameters on the flow of liquid through the pits. Only the size of the pit aperture diameter was adjusted, while other structural parameters remained unchanged. The simulation results showed that the flow rate through the pit remained constant, despite changes in the pit diameter. This is because when the flow velocity and pit diameter are constant, the overall flow rate remains unchanged. When the pit diameter increased, the total pressure drop decreased. This is because the increase in pit aperture size means that the distance between the tips of the pit border increases and the local resistance loss decreases, so that the total pressure drop flowing through the pit decreases. Under the condition of maximum or minimum diameter, the proportion of flow resistance at the pit border to the total flow resistance in the pit differs by more than 40%. Figure 6c shows the influence of the torus diameter on the pressure drop and flow rate. As the torus diameter increased and other structural parameters remained unchanged, the total pressure drop of the pit increased significantly. The torus is generally considered to exhibit extremely low permeability [26,27], and the flow of liquid through the torus is negligible. If the size of the pit membrane is constant, the torus becomes larger, and the area of the pit margo inevitably decreases; the circulation area also decreases accordingly. Additionally, because the flow rate is constant, the velocity of flowing through the membrane will increase, and the local resistance loss will also increase, thereby increasing the pressure drop. When the torus diameter was small, the flow resistance of the torus accounted for 22.12% of the total flow resistance of the pit. As the torus diameter increased, the flow resistance of the torus accounted for up to 33.6% of the total flow resistance of the pit.

Pit depth is an important feature of the pit structure, and is the distance between two adjacent pit borders (through the pit membrane in the middle). The influence of pit depth on the pressure drop and flow rate is shown in Figure 6d. The larger the pit depth, the larger is the capacity of the chamber. With the increase in pit depth, the flow efficiency of water in the pit was increased, and thus, the total pressure drop of the pit was reduced. This is because when the liquid normally flows through the pits, the pit membrane is in the middle of the pit cavity. As the pit depth becomes larger, the space between the pit border and the pit membrane increases, and the entire pit cavity becomes larger. Consequently, the local resistance loss is reduced, and the flow in the pit becomes smooth, thus decreasing the pressure drop.

Figure 6e shows the influence of margo thickness on the pressure drop and flow rate. The model only adjusted the margo thickness, while other structural parameters remained unchanged. As the margo thickness increased, the flow rate of the pit remained unchanged but the total pressure drop of the pit increased. When the thickness of the margo increases, the pores in the margo are either broken down by the interwoven microfibrils into smaller pores or covered by amorphous materials to lengthen the channel, and the resistance of the fluid passing through the margo increases, thus increasing the total pressure drop of the pit.

3.3. Flow Analysis of Margo Pores

The flow of liquid in the pit structure is realized by the pores woven by microfibrils on the margo. The margo model of Mongolian Scots pine contains 1197 pores of different sizes and diverse shapes, including triangles, diamonds, and circles. Calculations showed that the average diameter and area of each pore were 0.38 μm and 0.116 μm², respectively, and the Reynolds number was very small, indicating that the laminar flow state was dominated by viscous force. The model solution was used to plot the flow velocity of the liquid through the margo pore (Figure 7). The flow rate is obtained by calculating the area of each pore and integrating the flow velocity [15]. Moreover, the pores are gradually added to the margo region of the model so that it is possible to see how additional pores affect the overall model solution. The approach we took was to gradually add the pores to the margo region of the model so that we could see how the additional pores affected the solution of the entire model. In the model development process, the larger pores near the pit torus were first added, after which the smallest pores were added after about 5–6 steps. Calculations were stopped until further addition of the remaining very small pores did not result in any further significant reduction in pressure drop. The speed scale on the left side of the image represents the maximum flow velocity. The maximum velocity of a single pore seems to depend largely on the size of the pore. The area of the largest hole in the middle position (near the left side) of the margo accounted for 2.0% of the total pore area, and the flow rate accounted for 4.7% of the total flow rate through the pit membrane.

By analyzing the relationship between pressure drop and pore area, we found that micropores on the outer edge of the margo have little effect on the pressure drop. Therefore, when determining the flow through a single pore in the margo, only the relationship between the position distributions of 136 relatively large pores in the membrane and the flow conditions is discussed (Figure 8). The position of the circle represents the flow rate through a single pore and the size of the circle indicates the radial position of the centroid of each pore in the margo; a small circle means close to the center and a big circle means close to the outer edge of the margo. Different color scales were selected to make the radial position of each pore visible. The flow rate through individual pores showed a nonlinear relationship with the pore area (R² = 0.917). However, the radial position of the pores was also important. The pores in the middle of the margo showed a higher flow rate than those on the outer edge. The largest pores in the middle almost dominated the flow, which is supported by the velocity shown in Figure 7.

A small number of macropores in the margo may have a great influence on the overall flow and membrane resistance [28]. To specifically analyze the relationship between the flow rate of individual pores and the pressure drop of the margo model under the combined influence of pore area and location, we considered the centroid of the pit membrane as the center of the circle, and divided the margo into five concentric zones (green dashed lines) starting from the porous area on the edge of the torus (radial range is 2–4.5 μm), each ring scale is 0.5 μm), as shown in Figure 9a. Since most of these pores were irregular in shape, the pore area was classified according to the annular area where the centroid of the pore was located, and the total flow rate of all margo pores in each annular area was calculated. Figure 9b shows the cumulative pore area and flow rate in each annular area, indicating that the pore size distribution and the pores position have a considerable impact on the pressure drop.

The area of margo pores in the 2.5–3 μm annular zone was the largest, accounting for 35.7% of the total area, and represented 52.3% of the total flow rate. This is because the number of big pores was greater than that of small pores in this area, and the distance between the pit border and pit membrane (depth of the pit chamber) was relatively large. The space was wide, and the flow resistance was relatively small, leading to increased flow rate. In the outer ring area of the margo (3.5–4.5 μm), the number of small holes was much greater than that of big holes. The outer ring area is located in a narrow position of the pit cavity, and the flow shear force and flow resistance are large, resulting in a small flow rate. Therefore, the velocity and flow rate in each pore of margo are the result of a combined effect of the pore area and the radial position of the pore.

3.4. Verification and Analysis of Flow Resistance

To verify the effectiveness of the LBM-C mesoscopic method proposed in this paper, based on the same six pits samples of the bordered pits of Mongolian Scots pine branches, we applied the CFD method [15,16,27,29] based on the macro-continuous model (represented by MA-C) for analysis and comparison. The flow of liquid through each component of the pits was calculation with two methods respectively, and the calculated results of the two methods for the flow resistance of each part in the pits were roughly identical (Figure 10a–f), which verifies the availability of the LBM-C model proposed in this paper for the characteristics of the pit structure. The average value of the total flow resistance of the six pit samples were 114.75 (LBM-C) and 124.22 (MA-C), and the value obtained by MA-C method 9.47% larger than the LBM-C method. Additionally, the flow resistance values of aperture and torus were not significantly different under the two analysis methods, whereas the flow resistance of the margo was always the highest, the flow resistance values caused by margo were lower in LBM-C method than that of MA-C method in comparison of six samples. As shown in Figure 10a, the resistance of LBM-C method to analyze margo is 90.9 × 10¹⁴ Pas m⁻³ and the result of MA-C method is 98.9 × 10¹⁴ Pas m⁻³, 8.8% more than that of LBM-C method. The average value was calculated of the flow resistance ratio of each component of these samples, as shown in Table 2. The resistance distribution trends of the components were roughly the same; compared with LBM-C method, the ratio of margo resistance by MA-C method is 3.88% larger.

The possible reasons for the difference were analyzed; on the one hand, the size of the bordered pit is very small (in the order of micron and sub-micron) and the structure of the bordered pit is highly complex [30]. The LBM-C method provided a more accurate estimation of the pit morphology and analyzed the adjustment and barrier mechanism of each part of the structure boundary to the hydraulic network on a mesoscopic scale. Moreover, the LBM-C method analyzed the relationship between the structure and flow in more detail, especially for the internal curved boundary of the pit borders forming the pit chamber and the margo structure with an irregular pore shape and size. Moreover, the LBM-C method could also perform a detailed analysis of the flow through pores in the 4–4.5 μm ring area away from the pit aperture. In addition, the algorithm process was relatively simple. It was not necessary to adjust the mesh density step-by-step according to the complexity of shape profile or pore size variation, and the deviation was relatively small and stable. We know that the micro-particle background of the LBM-C method makes the method more intuitive and convenient to analyze the interaction between the fluid and surrounding environment. As a result, it overcomes the disadvantages of traditional numerical methods in dealing with complex numerical boundaries, such as mesh reconstruction or using a porous plate instead of actual margo analysis [14,31]. With respect to the capability of analyzing the interaction between the complex structure inside the pit and fluid particles, the LBM-C method is intermediate between the micro-fluid and macro-continuous models, demonstrating that it is more applicable than the macro-continuous model. This may explain the ability of the LBM-C method to decrease the resistance value to a greater extent compared to the macro-continuous method while analyzing the flow in the margo. Because of the different shapes and sizes of the pores on the margo, the LBM-C method can capture the interface flow field and flow shear force more accurately and analyze the relationship between structure and flow in greater detail, thus improving the accuracy of the analysis of flow resistance and flow and enabling the numerical study of the flow field model based on individualized, actual bordered pits.

The LBM-C model proposed in this study aids in an exploratory analysis of the internal flow field of perforated structures and a feasibility analysis in applying this method to simulate flow problems in complex biological systems. In further studies, we will explore the application of LBM-C to investigate seepage in the xylem flow network composed of tracheids, rays, bordered pits, and cross-field pits of early wood and late wood. For bordered pits, the two-phase flow state of air on one side of the pit membrane [32,33], which is the same as that in case of embolized tracheids [34,35], should also be considered. The LBM-C model can more accurately calculate the interface dynamics in vapor–liquid two-phase flow, dynamically simulate the complex flow field inside the xylem, and, at the same time, comprehensively analyze the axial movement inside the tracheid and lateral movement inside the pit. Alternatively, the dynamic relationship of the interaction between the pit structure and fluid is simulated from the perspective of fluid–solid coupling, and the tensile tension and stress distribution of the pit border and the pit membrane are described. The selection of parameters is used to establish the liquid circulation pattern of the torus at different positions in the pit chamber, in order to infer more information about the mechanical properties of the conifer pit membrane. Our goal is to better understand the relationship between the tracheid–pit structure and permeability characteristics of the coniferous wood xylem through flow simulation, and elucidate how the specific structural details of individual pits facilitate the fluid flow and prevent embolism from spreading. This study could deepen our understanding and grasp of the liquid conductivity of the xylem, on which conducting experimental research on the cellular scale is not easy, and we can fully understand the factors influencing the accessibility and permeability of the porous network of xylem and improve our understanding of the gas and liquid conductivity of conifer xylem.

4. Conclusions

The structure of the bordered pit is essential to maintain hydraulic conductivity and in wood treatment. Due to its mesoscopic dynamic characteristics, LBM is highly suitable to simulate flow problems in complex biological systems. In this study, by modifying the curved boundary processing scheme, we proposed the use of LBM-C method to analyze the complex flow field inside the pit, and obtained the relationship of the bordered pit structure components with velocity, flow rate, pressure drop, and resistance. The results show that when the diameter of the pit increased, flow rate increased, whereas total pressure drop decreased. At a constant pit diameter, flow rate remained unchanged. Total pressure drop decreased with increase in pit diameter or depth; however, it increased with the increase in torus diameter or margo thickness. The main structure affecting flow resistance is the margo. The flow through the margo is the result of the joint action of the pore area and position on the margo, a few large holes in the inner ring region of the margo contribute to most of the flow. Finally, the classical macro-continuous method was used to carry out comparative verification. Through the comparative analysis of resistance data generated by the components in six bordered pits samples, the average resistance ratio difference of each part obtained by the two methods is less than 5%, which verifies the effectiveness of LBM-C method. Moreover, through the difference analysis of the calculation results of the two methods and the comparison of the operation process, the efficiency of LBM-C method is revealed. This study provides a promising method for analyzing the physiological and ecological functions of conifers and realizing the efficient utilization of wood resources, which has a positive impact on the selection, cultivation, processing, and utilization of wood.

Author Contributions

Conceptualization, W.Q.; Data curation, W.Q.; Investigation, Y.J.; Methodology, Y.M.; software, G.Q.; Writing—original draft, W.Q. and C.Y.; Writing—review and editing, W.Y. validation, G.Q. and Y.J. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by “Fundamental Research Funds for the Central Universities (2572020DR12)” and “Major special research and development projects in Guangdong Province (2020B020216001).”

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Acknowledgments

We thank Yuan-yuan Miao and Wheeler for fruitful discussion on xylem anatomy and technical aspects of the measurements.

Conflicts of Interest

The authors declare no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures and Tables

Figure 1. Schematic diagram of the pit structure in Mongolian Scots pine. (a) Sectioning of wood segments and setting of coordinate axes. (b) SEM image of wood cross section with some early and late tracheids (RT plane). (c) The SEM image of the bordered pit membrane (LR plane). (d) The SEM images of tracheids of the earlywood from radial plane (LR plane). (e) The SEM images of the pits on the tracheids of early wood from the radial plane (LR plane). (f) The TEM images of the pit from the RT plane.

Figure 2. Schematic diagram of flow field lattice analysis of the pit boundary by the curved boundary treatment.

Figure 3. The scanning electron microscopy (SEM) image of the pit membrane.

Figure 4. Pit structure of the Mongolian Scots pine. (a) Transmission electron microscope (TEM) image of Mongolian Scots pine bordered pits. (b) Approximate structure and size of extracted the pit. (c) The geometry of the 3D pit model.

Figure 5. Velocity (m s−1) through different parts of the pit model in the section perpendicular to the pit membrane (LT plane). (a) No border pit. (b) Pit with borders but no membrane. (c) Pit with borders and a torus but no margo. (d) Full pit model.

Figure 6. Relationship of the structural parameters of the pit with the total pressure drop and flow rate. (a) The relationship of the different pit diameter parameters with the pressure drop and flow rate. (b) The influence of different aperture diameters on the pressure drop and flow rate. (c) The influence of torus diameters on the pressure drop and flow rate. (d) The influence of pit depth on the pressure drop and flow rate. (e) The influence of margo thickness on the pressure drop and flow rate.

Figure 7. Velocity (m s−1) of flow through the margo pores.

Figure 8. Relationship of the area and position of individual pores in margo with the flow rate.

Figure 9. Cumulative flow values of different radial zones in the margo. (a) The center of the pit membrane is the center of the circle, and the margo is divided into five concentric regions (green dashed line) from the ring edge to the boundary between margo and torus (radial range is 2~2.5 μm, 2.5~3 μm, 4~4.5 μm). (b) The cumulative pore area and flow rate in each annular area.

Figure 10. Comparison of the hydraulic resistance attributable to each component in six bordered pits analyzed using the LBM-C and MA-C methods. (a–f) are the resistance of aperture, torus, margo and the total from 6 different samples of Mongolian Scots pine. The blue column represents the calculation result of the LBM-C method, and the gray column represents the result of the MA-C method.

Table 1

Pressure drop (ΔP) along the model and calculated resistance (R) and proportion of resistance (%) for various pit components.

Component	ΔP (Pa)	Resistance (Pa s m−3)	Fraction of Resistance (%)
No Border	0.573	1.747 × 1013	0.116
Border	31.21	9.515 × 1014	6.36
Torus	127.4	3.884 × 1015	25.97
Margo	331.3	1.009 × 1016	67.55

Table 2

Fraction of average resistance (%) attributable to the three components of the bordered pits analyzed by the LBM-C and MA-C methods.

Sections	Resistance LBM-C (%)	Resistance MA-C (%)
Aperture	15.07	12.02
Torus	20.42	19.59
Margo	64.51	68.39