#### Abstract

Natural fiber bundle like hemp fiber bundle usually includes many small lumens embedded in solid region; thus, it can present lower thermal conduction than that of conventional fibers. In the paper, characteristic of anisotropic transverse thermal conductivity of unidirectional natural hemp fiber bundle was numerically studied to determine the dependence of overall thermal property of the fiber bundle on that of the solid region phase. In order to efficiently predict its thermal property, the fiber bundle was embedded into an imaginary matrix to form a unit composite cell consisting of the matrix and the fiber bundle. Equally, another unit composite cell including an equivalent solid fiber was established to present the homogenization of the fiber bundle. Next, finite element thermal analysis implemented by ABAQUS was conducted in the two established composite cells by applying proper thermal boundary conditions along the boundary of unit cell, and influences of the solid region phase and the equivalent solid fiber on the composites were investigated, respectively. Subsequently, an optional relationship of thermal conductivities of the natural fiber bundle and the solid region was obtained by curve fitting technique. Finally, numerical results from the obtained fitted curves were compared with the analytic Hasselman-Johnson’s results and others to verify the present numerical model.

#### 1. Introduction

Natural fibers like kenaf fiber [1], hemp fiber [2], sisal fiber [3], date palm fiber [4], wood fiber [5], and bamboo fiber [6] have unique advantages of low density, high specific properties, biodegradable nature, and low cost; thus, composites filled with natural fibers, such as natural fiber reinforced polymer/cement composites, are usually viewed as green and environmentally friendly composites and have attracted much attention of researchers for potential engineering application. As one of inherent material properties of natural fibers, thermal property of natural fibers is of great importance in natural fiber reinforced composites, due to inherent hollow microstructure of natural fibers. Recent researches have shown that natural fibers consisting of cellulose or lumens can present extremely lower thermal conduction than conventional fibers like glass fibers and carbon fibers [7]; thus, natural fibers reinforced composites can be considered to be thermal insulator in such engineering as building and furniture. In addition, it is viable to achieve the aim of light-weight and proenvironment composite materials by considering natures of hollow microstructure and biodegradability of natural fibers.

In the past few years, many researchers investigated thermal properties of natural fibers and composites filled with them. For example, El-Shekeil et al. experimentally investigated the influence of fiber content on the mechanical and thermal properties of kenaf fiber reinforced thermoplastic polyurethane composites [8]. Liu et al. evaluated the transverse thermal conductivity of Manila hemp fiber in solid region by the finite element method and analytical Hasselman-Johnson’s model [2]. Also, they studied the effect of the microstructure of natural fiber on the transverse thermal conductivity of unidirectional composite with abaca and bamboo fibers, by experiment and finite element simulation [9]. Behzad and Sain predicted the thermal conductivity for hemp fiber reinforced composites by experimental measurement [10], and subsequently they developed a finite element simulation procedure to predict the temperature profile and the curing behavior of the hemp fiber/thermoset composite during the molding process [11]. Mangal experimentally measured the thermal properties of pineapple leaf fiber reinforced composites [12]. Takagi and cooperators analyzed thermal conductivity of PLA-bamboo fiber composites using a hot-wire method [13]. All these works mentioned above have been beneficial in understanding of thermal transfer mechanism in the natural fiber and design of natural fiber filled composites with desirable thermal properties.

As important fillers of green composites, it is necessary to establish comprehensive understanding of thermal properties of natural fibers or fiber bundles. In this paper, the emphasis is put on the study of thermal properties of the natural hemp fiber bundle. From the viewpoint of composite, the natural hemp fiber bundle can be viewed as composite material with outstanding thermal properties, because, in the natural fiber bundle of interest, there are a large number of lumens filled with air in transverse direction of it (see Figure 1(a) for the cross-section morphology of the hemp fiber bundle [2]). The thermal properties of the natural fiber bundle vary considerably depending on lumen volume and size and also the thermal property of solid region phase, which encloses the lumens in the natural fiber bundle. Figure 1(b) shows a schematic illustration of the natural fiber bundle consisting of lumens and solid region. It is observed in Figure 1 that the large-scale fiber bundle is filled with many small-scale lumens in the solid region; thus, thermal properties of the lumen and the solid region are important parameters of the natural fiber bundle. In practice, the lumen filler is filled with air; thus, its thermal conductivity is usually specified with very small value, for example, 0.026 W/(mK), which is normally the thermal conductivity of air measured at the standard atmosphere. Therefore, the thermal conductivity of the fiber solid region more significantly affects the whole thermal performance of the fiber bundle than that of lumen. Here, the main purpose of this study is to investigate the effect of material thermal property of the solid region phase on the equivalent anisotropic thermal property of the natural fiber bundle by finite element simulation [14, 15] of composite microstructure [16] and then establish an optional interrelationship between them by curve fitting technique [17] to give a rapid and highly accurate prediction of material thermal properties for both of them.

**(a)**

**(b)**

#### 2. Finite Element Model for Anisotropic Natural Fiber Bundle Reinforced Composites

In this paper, the natural fiber bundle is assumed to be embedded into a polymer matrix with constant thermal conductivity to form a square representative volume element (RVE) (or unit cell) (see Figure 2), as was done by many researchers in the analysis of heterogeneous materials [18–20]. Then, the finite element model of natural fiber bundle reinforced composites [11, 15, 16] will be established to investigate the influence of solid region phase on the fiber bundle. The established unit cell consists of three different regions, that is, matrix, solid region, and lumen. Each region has isotropic thermal conductivity.

**(a)**

**(b)**

Theoretically, the distribution of lumens in the natural fiber bundle will cause anisotropy of the composite under consideration. Therefore, in this study, the anisotropic thermal conductivities of the composite will be investigated. In steady-state heat conduction problem, the temperature field within the anisotropic representative volume element satisfies the quasiharmonic differential equation: where the thermal conductivities and are piecewise constant. Because the matrix, the solid region, and the lumen are assumed to be locally isotropic and homogeneous, in the matrix, in the solid region and equal to in the lumen, respectively.

Besides, in the heat conduction system, the heat flux components and are, respectively, defined by temperature gradient: In the present composite computational model, the proper thermal boundary conditions should be applied along the boundary of the cell shown in Figure 2 to construct a complete composite heat transfer system being solved by finite element technique [15, 21, 22], which has been successfully employed by many researchers for the analysis of effective thermal properties of unidirectional fiber composites [23–25]. According to the work of Islam and Pramila [24], the prescribed temperature boundary conditions , on the vertical or horizontal boundaries of the cell can produce the most accurate results up to a relatively high fiber volume fraction, and the remaining boundaries of the cell are assumed to be insulating, as shown in Figure 2. It is assumed that ; thus, the average heat flux components and loaded on the data collection face, for example, the left side face for the case (a) and the bottom side face for the case (b), are positive. As a result, based on Fourier’s law in heat transfer (see (2)), the effective thermal conductivities of the composite can be given as [24, 25] where is the side length of the square cell, for the case in Figure 2(a), and for the case in Figure 2(b). The integrals in (4) and (5) can be evaluated by trapezoidal numerical integration.

In the practical computation, the side length of the square cell is set to be 1, which is a normalized length. If the volume fraction of the fiber bundle to the cell is assumed to be a moderate value of 50%, the normalized radius of the fiber bundle is 0.3990. Furthermore, if the volume fraction of lumen to the fiber bundle keeps constant, that is, 30.87%, which is the experimental result [2], the normalized radius of the lumen can be evaluated by , where is the number of lumens in the fiber bundle. For example, if the number of lumens is 106, which is close to the actual disperse of lumens in the fiber bundle (see Figure 1), then the normalized radius of each lumen is 0.0215. Besides, the specified temperature boundry conditions along the two vertical edges of the unit cell are set to be 20 and 0, respectively. Moreover, the thermal conductivities of matrix, solid region, and lumen are, respectively, normalized with the reference value 0.026 W/(mK), which is the thermal conductivity of lumen, in the finite element procedure below. In the paper, the symbols , , and , respectively, indicate the thermal condctivities of matrix, solid region, and lumen.

#### 3. Numerical Results and Discussions

##### 3.1. Convergence Investigation

Generally, the finite element (FE) solution will be more accurate as the model is subdivided into smaller elements. The only sure way to know if we have sufficiently dense mesh is to make several models with different grids of elements and check the convergence of the solution. In order to investigate the convergence of the FE solution, the composite model in Figure 2(a) is studied by ABAQUS and the FE size is changed from very coarse to very fine. The element type employed in ABAQUS is DC2D8. In each of FE size levels, the average horizontal heat flux component at the left wall of the square unit cell is calculated. Table 1 gives a summary of the output of these size levels indicating the number of elements and nodes used in the computational domain corresponding to each element size. In this table, Size 1 stands for very coarse elements and Size 4 means very fine elements. The table also indicates the deviation between the average heat flux component at the left wall of the square unit cell calculated using various element sizes and that calculated using the finest element size of 4. The summarized results in Table 1 indicate that the maximum deviation between the solution using the finest element of 4 that corresponds to 51680 elements and the coarsest element of 1 that corresponds to 2124 elements is 0.0302%. This reflects clearly that the numerical solution obtained via this FE simulation is mesh size independence. Additionally, looking for high accuracy, the authors decided to use a fine element size of 4 in the following computation.

##### 3.2. Anisotropy Investigation

It is known that the distribution mode of lumen may cause anisotropy of both fiber bundle and composite. To investigate this effect, let us consider the composite model involving polymer matrix, solid region, and lumens, as displayed in Figure 3(a), in which 106 lumens are regularly distributed in the fiber bundle to approximate the real distribution of lumens in the practical natural fiber bundle (see Figure 1). Figure 3(b) presents the computational mesh of element Size 4. Results in Table 2 display the change of anisotropic thermal conductivities of the composites for various thermal conductivities of the solid region phase in the natural fiber bundle. It is obvious that the thermal conductivities of the composite along two directions are extremely similar, so it is concluded that the approximated practical distribution of lumen in the fiber bundle causes the composite and fiber bundle isotropy. The similar conclusion was drawn by Liu et al., who predicted that the anisotropy of the composite became smaller with the number of lumens increasing [2]. Therefore, it is reasonable to assume the isotropy of the composite and the fiber bundle in the following analysis.

**(a)**

**(b)**

##### 3.3. Effect of the Thermal Conductivity of the Solid Region in the Natural Fiber Bundle

To estimate the effect of the solid region on the composite, it is assumed that the normalized thermal conductivity of the solid region changes in the interval [1, 10]. By finite element computation, the distribution of horizontal heat flux component and the corresponding average value of it can firstly be obtained on the data-collection surface for each specific value of the normalized thermal conductivity of the solid region. Then, the normalized effective thermal property of the composite can be evaluated by (3). The variation of the effective thermal conductivity of the composite is given in Figure 4, from which it is observed that the simulated effective thermal conductivity of composite increases with the increasing thermal conductivity of solid region, as we expect. Simultaneously, it is found that the variation shown in Figure 4 shows slight nonlinearity, instead of linearity. Thus, to describe the nonlinear variation shown in Figure 4, the following quadratic and cubic polynomial expressions from curve fitting technology are, respectively, employed.(i)Quadratic polynomial fitting: (ii)Cubic polynomial fitting:

##### 3.4. Effect of the Thermal Conductivity of the Homogenized Fiber Bundle

In this section, the composite model shown in Figure 5(a) is taken into consideration to investigate the effect of the homogenized fiber bundle on the composite. In the model, the homogenized fiber bundle is represented by a solid fiber with the same size. Also, the same thermal boundary conditions as those used in the composite model above are applied along the outer boundaries of the cell. To conduct the finite element analysis, a total of 12247 quadratic quadrilateral elements of type DC2D8 with 37142 nodes are generated by ABAQUS to discretize the computational domain (see Figure 5(b)).

**(a)**

**(b)**

It is assumed that the normalized thermal conductivity of the solid fiber changes in the interval [1, 10]; thus, the effective thermal conductivity of the composite can be evaluated by finite element simulation for any specific value of thermal conductivity of the equivalent solid fiber. The variation of the effective thermal conductivity of the composite against the equivalent solid fiber is displayed in Figure 6, which clearly shows that the simulated effective thermal conductivity of the composite nonlinearly increases with the increasing value of the thermal conductivity of the equivalent solid fiber. To accurately capture the nonlinear variation shown in Figure 6, the following quadratic and cubic polynomial curves are, respectively, employed by means of curve fitting technology.(i)Quadratic polynomial fitting: (ii)Cubic polynomial fitting:

##### 3.5. Optional Interrelationship between the Thermal Conductivity of the Solid Region and the Fiber Bundle

Finally, the equivalence of the two composite models, respectively, shown in Figures 3(a) and 5(a) requires that the two composite models should have same effective thermal conductivities. Therefore, combining (6)–(9), we have an optional interrelationship between the thermal conductivity of the solid region and the fiber bundle; that is, or from which the variation of in terms of is plotted in Figure 7.

To verify the obtained relation of thermal conductivity between the solid region and fiber bundle, the theoretical Hasselman-Johnson’s model derived from the interface interaction between the circular matrix and circular inclusions embedded in the matrix [2] is taken as reference for the purpose of comparison. Here, an analytical expression from the Hasselman-Johnson’s model can be employed to investigate the thermal conductivity of the solid region with respect to that of the fiber bundle ; that is, where represents the volume fraction of lumen to the fiber bundle.

According to the experiment data in [2], the practical volume content of the lumen to the fiber bundle is about 30.87%; thus, the substitution of into (12) yields from which one can get the variational curve of in terms of , as shown in Figure 7 for the purpose of comparison. Specially, if is taken to be 0.115 W/(mK) [2, 10], the thermal conductivity of the solid region is calculated by (13) as 0.1847 W/(mK).

In Figure 7, it is observed that there is good agreement between the numerical results from either quadratic or cubic curves and the theoretical result of Hasselman-Johnson’s model for the case of moderate change of . For example, for the case of W/(mK), the thermal conductivity of the solid region is calculated as 0.1843 W/(mK) for the quadratic approximation and 0.1857 W/(mK) for the cubic approximation, which has relative derivation of 0.21% and 0.54% of the theoretical solution 0.1847 W/(mK), respectively. Therefore, both quadratic and cubic relations of thermal conductivity of the solid region and the fiber bundle can be used to evaluate thermal properties of the natural fiber bundle or the solid region in the bundle. Also, in contrast to the analytical expression (13), it can be seen from (10) and (11) that either or is given, and one can easily determine another. This means that we can do inverse prediction conveniently by a specified material thermal conductivity of the solid region. This is an advantage of the optional relation presented in the paper over the analytical solution.

Besides, it is obvious in Figure 7 that the existence of lumen significantly weakens the capacity of heat transmission in the fiber bundle. As a result, the thermal conductivity of the fiber bundle is greatly less than that of .

#### 4. Conclusion

In this paper, 2D computational composite model of the natural fiber bundle including numbers of lumens is developed. Due to the geometrical limitation of the fiber bundle, it is not convenient to directly apply thermal boundary conditions to it to perform finite element analysis of composite. To treat this, the fiber bundle is assumed to be embedded into a matrix with known thermal conductivity to construct unit composite cell, which is numerically analyzed by applying proper thermal boundary conditions along the cell boundary. By means of the developed finite element computational composite model, the effect of the solid region in the bundle on the overall thermal property of the composite is studied. Simultaneously, a homogenized composite model is constructed, in which the fiber bundle is replaced by an equivalent solid fiber to investigate the influence of the homogenized fiber bundle. By comparing the two composite models developed in this study, an optional interrelationship between thermal conductivities of the solid region and the homogenized fiber bundle was obtained by curve fitting technique. Finally, the present computational composite model is verified and numerical experiments show that either quadratic or cubic predictions can produce almost similar results for the solid region in the fiber bundle, in contrast to the theoretical Hasselman-Johnson’s model and other numerical results. Moreover, the direct or inverse predictions can be easily performed to evaluate the thermal conductivity of the fiber bundle or the solid region in the practice, if one of them is given. More importantly, the present computational method can be easily extended for the prediction of thermal property of other natural fiber bundles with various lumen patterns.

#### Conflict of Interests

The authors declare no conflict of interests regarding the publication of this paper.