Stochastic Thermal Properties of Fibrous Composites Reinforced with Long, Randomly Distributed Fibers

Material properties are determined during physical tests. The kind of test depends on the kind of properties which are going to be determined. A variety of methods can be used to measure thermal conductivity, descriptions of which can be found, among others, in [1].

One of the most widely used approaches for determining effective thermal conductivity is based on simulation methods, where numerical analyses of heat transfer through a virtual material specimen are performed. These methods most often imitate the physical absolute technique, which belongs to the class of steady-state methods [1]. Performing such simulations requires detailed knowledge of the composite’s internal structure and the material properties of its individual phases. The heat transfer problem in these simulations is typically solved using the finite element method (FEM) in order to obtain relatively accurate results. The finite element method is widely used in solving many engineering problems, such as the static and dynamic behavior of mechanical structures, electromagnetics, fluid mechanics, and heat flow processes. There are thousands of publications where authors use and develop this method in all possible fields of its usage.

Particularly in the fields of thermal and coupled thermo-mechanical problems, the finite element method has been successfully applied numerous times, demonstrating good agreement with experimental data. In [2], the authors used a finite element model to calculate strain and residual stresses from the temperature field. The numerical results were compared with experimental measurements, showing good consistency. The analysis of welded components was presented in [3], where the authors compared experimentally measured deformations with FEM results and achieved good agreement. As another example of successful FEM application, [4] considered multilayered plates produced by electron beam additive manufacturing. In that study, the authors developed a thermo-mechanical model to investigate distortions and residual stresses in the fabricated plates. The FEM results were validated against experimental data to ensure calculation accuracy and to support subsequent optimization processes. In [5], finite element analysis was used to simulate the fabrication of ultra-thin, smooth metal coatings. The finite element model was successfully validated using experimental data for the same materials and process parameters.

Numerical investigations of effective thermal conductivity in composite materials constitute one of the potential applications of the finite element method. In [6], the authors employed FEM to predict the thermal conductivity of particle-filled epoxy composites. They developed a three-dimensional model that took into account the filler geometry and its aspect ratio. In [7], the authors proposed a new method for generating a three-dimensional FEM model with randomly arranged fillers for the calculation of thermal conductivity in composites with randomly dispersed particles. In [8], a finite element model was developed to determine the effective thermal conductivity of asphalt concrete with a random aggregate microstructure. The simulation results were compared with experimental data, demonstrating good agreement and accuracy.

The statistical analysis of material properties in composite materials has been addressed in several studies. A stochastic micromechanical model for predicting the probabilistic characteristics of elastic mechanical properties in isotropic functionally graded materials was considered in [9]. In [10], the probabilistic behavior of hybrid fiber-reinforced concrete was investigated, focusing on its mechanical properties. The stochastic behavior of mechanical properties in composites with inclusions of arbitrary shapes was studied in [11]. The study of statistical parameters of effective thermal conductivity in composites reinforced with long, parallel fibers was presented in [12]. In that work, small perturbations in fiber positions were assumed, limited only to the square region of the representative volume element (RVE) of a single fiber.

A commonly used structural material is carbon-fiber-reinforced epoxy composites. The popularity of these composites results from their favorable strength parameters. The characteristics of these materials include not only elastic and strength coefficients but also coefficients describing their thermal properties. These materials exhibit higher thermal conductivity coefficients than commonly used composites based on epoxy resin reinforced with glass fibers. Such behavior of carbon-fiber composites is a direct consequence of the high thermal conductivity of carbon fibers. Carbon fibers are characterized by a wide variation in thermal conductivity values. In practice, carbon fibers can be found with transverse thermal conductivity values ranging from approximately 10 W/(m·K) for standard-grade fibers up to 160 W/(m·K) for specialized fibers [13]. Extensive studies on the effective thermal conductivity coefficients, in the temperature range of 25–130 °C, of carbon-fiber-reinforced epoxy composites—supported by both physical experiments and modeling analyses—can be found in [14]. This article presents detailed analytical approaches for calculating the effective thermal conductivity coefficients in directions perpendicular and parallel to the fiber orientation, as well as in arbitrary directions. Among other articles dedicated to the study of thermal conductivity coefficients of carbon-fiber-reinforced epoxy composites, examples include [15, 16, 17]. In these articles, the investigations of thermal conductivity coefficients are carried out using various physical measurement methods, often supported by modeling calculations.

In the present paper, the stochastic parameters of the transverse effective thermal conductivity of composites reinforced with randomly distributed fibers of circular cross-section are considered. It is assumed that these fibers are mutually parallel and parallel to the thermally loaded surface. The required analysis of thermal fields was conducted on a virtual experimental setup. Models of the setup and the composite specimen tested were developed using the finite element method. This article is a continuation of the previous work published in [12].

All necessary procedures and functions were written in the FORTRAN programming language (GNU FORTRAN Compiler 6.3.0-1).

The planned study requires performing virtual measurements of the thermal conductivity coefficient (i.e., selecting the type of test and the method of its execution) for multiple virtual specimens of a material containing parallel but randomly distributed fibers. The results of these simulations will then be processed statistically. At this stage, appropriate probability density functions will be determined, based on which the mode and standard deviation will be calculated.

2.1.

Virtual experimental setup

To accomplish the stated objective, a series of virtual heat conduction tests must be conducted in order to determine the thermal conductivity coefficient of the materials investigated. The virtual test employed in this study is based on physical measurement methods for determining such coefficients, which rely on a steady-state heat flow [1]. A schematic diagram of the virtual experimental setup for measuring thermal conductivity coefficients is shown in Figure 1, which was originally proposed in [18].

Fig. 1.

Schematic diagram of the virtual experimental setup

The configuration of the virtual measurement setup consists of two main components - the specimen and the heat diffuser. In this case, the specimen represents a cross-section of a fibrous composite segment. The purpose of the heat diffuser is to equalize the temperature within the contact region between the diffuser and the specimen. The system is subjected to a constant heat flux q along boundary Γ_q and a prescribed temperature T₂=0 along boundary Γ_T2, while boundaries Γ_i1 and Γ_i2 are thermally insulated. For a relatively high thermal conductivity of the diffuser material (several orders of magnitude greater than the thermal conductivities of the composite phases), the temperature T₁ on boundary Γ_T1 can be considered uniform and must be determined through numerical simulation. Assuming a constant temperature T₁, the effective thermal conductivity in the direction normal to boundary Γ_q can be expressed as: (1) λe=hT1q {\lambda _e} = {h \over {{T_1}}}q where h denotes the thickness of the specimen in the direction of heat flow analyzed.

2.2.

Numerical solution of the heat conduction problem

At the stage of analyzing the temperature field distribution within the specimen investigated, the finite element method (FEM) was applied. A key task in developing the finite element model is the discretization of the specimen domain. In this study, a scheme analogous to that proposed in [12] was used, in which the entire cross-sectional area of the composite matrix and fibers was discretized with 3-node triangular elements.

During the generation of the random fiber distribution, a simple algorithm was employed, in which a regular hexagonal fiber arrangement was randomly perturbed. A predefined number of fibers was introduced into the matrix domain so that they could uniformly fill the specimen’s cross section. The computational fiber radius was adjusted to ensure the required fiber volume fraction within the matrix. The fibers were introduced individually, based on their positions in the original hexagonal grid. From the node of the hexagonal grid with coordinates x_n the fiber was displaced in a random direction defined by a unit vector m and random fraction ξ ∈ (0,1] of the prescribed maximum displacement M to a new position x_t, such that: (2) xt=xn+ξ M m {{\rm{x}}_t} = {{\rm{x}}_n} + \xi \;M\;{\rm{m}}

The scheme described is illustrated in Figure 2.

Fig. 2.

Schematic of the fiber displacement to a random position

The admissibility of each generated fiber position was verified by checking for overlaps with previously generated fibers and with the external boundaries of the specimen. In the case of a conflict, a new fiber position was generated, but no more than a predefined number of times. If a fiber could not be generated in an admissible position, the maximum displacement M was reduced to hM, where η∈(0,1), and the generation process was repeated. The fiber generation process was terminated when all fibers were successfully assigned admissible positions.

To illustrate the quality of the algorithm proposed, several generated meshes are presented in Figures 3, 4, and 5. The meshes shown correspond to fiber volume fractions of ρ =0.20, 0.40 and 0.60.

Fig. 3.

Example of random fiber distribution for ρ = 0.20

Fig. 4.

Example of random fiber distribution for ρ = 0.40

Fig. 5.

Example of random fiber distribution for ρ = 0.60

The calculations were carried out for a fibrous composite consisting of an epoxy resin matrix reinforced with two types of carbon fibers. The thermal conductivities were assumed to be λ_m = 0.2 W/(m·K) for the resin and λ_f = 10 and 100 W/(m·K) for the carbon fibers, while the fiber volume fractions in the composite were taken as ρ = 0.2, 0.3, 0.4, 0.5 and 0.6.

Each specimen consisted of 36 fibers. The number of fibers adopted in the calculations was a consequence of both computational efficiency and the possible variability of fiber arrangements within the specimen. The base specimen was built with six rows, each consisting of six hexagonal cells. A representative (i.e., appropriate for the study conducted) number of fibers is difficult to determine; consequently, the simulation results obtained can be considered approximate.

For specimens with different random fiber arrangements and for each combination of the input data considered, 1000 thermal conductivity tests were performed. Subsequently, the probability density functions were fitted, and the modes and standard deviations were determined.

Figures 6–9 present four example normalized histograms and the corresponding probability density function plots for the following data sets: λ_f=10, ρ=0.3; λ_f=10, ρ=0.6; λ_f=100, ρ=0.3; λ_f=100, ρ=0.6.

Fig. 6.

Histogram and gamma distribution for λ_f =10, ρ=0.3

Fig. 7.

Histogram and gamma distribution for λ_f =10, ρ=0.6

Fig. 8.

Histogram and gamma distribution for λ_f =100, ρ=0.3

Fig. 9.

Histogram and gamma distribution for λ_f =100, ρ=0.6

Tables 1 and 2 present the calculated distribution parameters α and β and the shift σ, as well as the computed modes and standard deviations for all fiber volume fractions ρ of the composite matrix analyzed. Table 1 contains data for cases with fibers of thermal conductivity λ_f=10 W/(m·K), while Table 2 refers to fibers with λ_f=100 W/(m·K).

It is worth noting the relatively large standard deviations, particularly for higher fiber volume fractions in the composite matrix. Figure 10 visualizes the variation of the modes as a function of the fiber volume fraction in the composite matrix.

Fig. 10.

Modes as a function of the fiber volume fraction in the composite matrix

The small differences observed in the effective thermal conductivity values (shown in Figure 10) of the composite materials analyzed, despite the large differences in the thermal conductivity of the reinforcing fibers, are a typical characteristic of fibrous composite materials. For high ratios of λ_f/λ_m the effective thermal conductivity depends only on the fiber volume fraction in the composite matrix. This property was demonstrated, among others, in [21] for a regular fiber distribution within the composite matrix.

To provide a more comprehensive illustration of the properties of the composites analyzed, the range of variation of the effective thermal conductivity was determined for the entire population of values obtained from 1000 tests conducted for each set of composite parameters analyzed. In addition, the relative difference (RD) between the maximum and minimum thermal conductivity values calculated was determined. The corresponding minimum, maximum, and relative error values for the two fiber thermal conductivities analyzed are presented in Tables 3 and 4.

To determine which fiber arrangements result in the material exhibiting the maximum and minimum thermal conductivity observed, Figures 11 and 12 show their configurations obtained for ρ=0.4 and λ_f=100 W/(m·K).

Fig. 11.

Fiber arrangement corresponding to the maximum thermal conductivity observed

Fig. 12.

Fiber arrangement corresponding to the minimum thermal conductivity observed

In the case of the maximum thermal conductivity of the composite, the fibers are arranged in chains aligned with the direction of heat flow, practically connecting the thermally loaded edges of the sample (Figure 11). Considering that the thermal conductivity of the fibers is much higher than that of the composite matrix, such a fiber arrangement appears natural—the fibers effectively form pathways for heat flux.

In the case of the minimum thermal conductivity of the composite, the fibers are arranged in chains perpendicular to the direction of heat flow, maintaining certain distances between each other and between themselves and the thermally loaded edges. This fiber arrangement introduces insulating layers within the composite structure, formed by the matrix material.

The study was conducted on composites consisting of typical phases - an epoxy resin matrix reinforced with carbon fibers. Due to the wide range of thermal conductivity properties of carbon fibers, two types of fibers with significantly different thermal conductivities were selected, while the composite matrix remained the same in all cases.

The simulation studies conducted enabled the determination of the most probable values of effective thermal conductivity and the corresponding standard deviations for fibrous composite materials with parallel, yet randomly distributed fibers. Analyses were carried out for a relatively large number of material samples (1000), which enhanced the reliability of the statistical evaluation of the thermal conductivity properties of these materials.

Attention should be paid to the relatively large variations in the effective thermal conductivity, which result from differences in the spatial distribution of fibers within the composite matrix. This behavior is reflected in the relatively high standard deviations and large relative differences of the extreme values of the composite’s thermal conductivity, as shown in Tables 1–4. The significant potential differences in effective thermal conductivity are particularly pronounced at higher fiber volume fractions in the composite matrix. Depending on the type of reinforcing fibers used, the relative difference between the maximum and minimum values of the composite thermal conductivity reached nearly 23% for fibers with a thermal conductivity of λ_f=10 W/(m·K) and 34% for fibers with λ_f=100 W/(m·K), as shown in Tables 3 and 4.

Relatively large differences in thermal conductivity for different fiber distributions indicate the possibility of affecting them through the appropriate arrangement of fibers within the composite matrix.

Furthermore, the experimental results indicated a small dependence of the thermal conductivity of the composites studied on the increase in the thermal conductivity of the reinforcing fibers (for high λ_f/λ_m ratios), while a stronger dependence was observed on the fiber volume fraction in the composite matrix, as shown in Figure 10.