1. Introduction
Thermoelectric material can transform thermal energy into electrical energy and vice versa. Its applications can be found in various fields such as temperature measurement, power generation, temperature control, etc. Conventional thermoelectric devices consist of both n-type and p-type thermoelectric legs. Within each leg, the heat flux and electrical current flow are parallel to each other [
1,
2]. Transverse thermoelectric devices make use of the transverse Seebeck effect so that the electrical current and heat flux flow perpendicular to each other [
3]. There are four main types of transverse thermoelectric materials, which are anisotropic single-crystal material, polycrystal material with engineered anisotropy, anisotropic organic thin-film thermoelectrics, and anisotropic thermoelectric composites. The single-crystal material anisotropy is the result of an unsymmetrical lattice structure [
4]; the material anisotropy in polycrystal is contributed to by both the lattice structure of the grain crystal and the grain geometry [
5]. The organic materials anisotropic property is caused by the preferred alignment of polymer chains [
6,
7,
8]. Studies have shown that the anisotropic thermoelectric composites can provide the best system efficiency among all of the candidates [
9].
The transverse thermoelectric composite can further be categorized into layered and fibrous thermoelectric composites. These two types of composites can also be seen as two-dimensional (2D) and one-dimensional (1D) inclusion composites. For layered composites with isotropic components, Babin et al. [
10] established a mathematical model for fast prediction of dimensionless transverse thermoelectric figure of merit (Z
transT), while other researchers performed finite element simulations and experimental tests and validated such a mathematical model [
11,
12,
13,
14,
15]. Similar to layered composites, fibrous composites can also provide appreciable Z
transT values. A mathematical model was established by Qian [
16] to study the thermoelectric performance of fibrous composites.
Until now, most studies on transverse thermoelectric composites assumed isotropic material properties in the component phase. However, many studies have already shown that thermoelectric materials can exhibit anisotropic properties, such as Bi
2Te
3 [
5,
17], SnSe [
18], and organic thermoelectric PEDOT:PSS (poly(3,4-ethylenedioxythiophene) polystyrene sulfonate) [
19,
20], etc. As a result, anisotropic material properties were applied into the component phase of layered transverse thermoelectric composites, and the results showed that the maximum Z
transT could be improved by introducing material anisotropy in a polycrystal [
21]. However, the effect of material anisotropy on fibrous composite still remains unknown.
In this study, anisotropic material properties were applied to the component phase of fibrous transverse thermoelectric composites. A mathematical model was established for predicting the effective material properties and thermoelectric cooling capacities of the composite. Finite element simulations were conducted to validate the mathematical model. A case study was conducted on a fibrous composite material containing Copper fibers and an anisotropic Bi2Te3 matrix, where the possible aggregate (anisotropic) properties for Bi2Te3 were calculated and applied into the component phase. The influence of the material anisotropy and the material property axis alignment in the composite were then discussed with respect to the maximum ZtransT and maximum cooling capacity of the composites. In the next step, the cooling performances of both the layered and the fibrous composites were compared. The thermoelectric properties of three representative composites, which were Bi2Te3/Copper, In4Se2.25/Copper, and SnSe/Copper, were applied into the component phase of the composites under different material property axis alignment configurations and extreme aggregate properties of the anisotropic component. The influence of material anisotropy on the enhancement of the maximum cooling capacity and maximum ZtransT of the transverse thermoelectric composites was demonstrated. The cooling capacities of the layered and fibrous composites were also compared. The derived mathematical model can serve as an efficient tool for selecting high-performance fibrous transverse thermoelectric composites, while the comparison between fibrous and layered thermoelectric composites can inspire and assist thermoelectric researchers in designing higher performance thermoelectric devices.
2. Mathematical Model for Fibrous Thermoelectric Composites with Anisotropic Components
In this section, analytical equations are derived for the effective properties of a transverse thermoelectric composite. In the mathematical model, unit cell structure is used to reduce the complexity of this problem. This unit cell model has been proven to be a convenient tool in studies on the effective properties of a material with periodical structures [
16,
22]. The schematic of a fibrous thermoelectric composite unit cell is shown in
Figure 1, where a fiber (F) with a square cross-section was placed at the corner of a cubic matrix (M). The square fiber was used in the unit cell model for mathematical simplicity. The contact between the matrix and fibrous material was assumed to be ideal so that electrical and thermal contact resistance were not considered.
When anisotropic material properties are introduced into the component phase (i.e., fibrous phase and matrix phase) of a composite material, the alignment of the material local property axis with respect to the composite material coordinate system plays a key role in the effective properties of the composite. Many 1D inclusions, such as carbon nanotubes [
23,
24], have different material properties along and perpendicular to the fiber axis direction. Therefore, in this study, we used the subscript “p” to represent the properties on the cross-sectional plane of the fiber and subscript “a” to represent properties along fiber axis (
Figure 1a). As for the matrix phase, we used the material local axis system “uvw” to describe the anisotropic material properties (
Figure 1b). There are four sets of coordinate systems in
Figure 1, which are the local material coordinate system of the fibrous and matrix phase, the coordinate system of composite C
1(F + M
1), and the coordinate system of unit cell C
2(F + M
1 + M
2). In
Figure 1, the coordinate axes u, x
1, and x
2 align parallel to each other; v, y
1, and y
2 align parallel to each other, and a, w, z
1, and z
2 align parallel to each other. The cross-sectional plane of fiber (p–p) is parallel with the u–v plane of the matrix material.
The effective properties of the unit cell can be derived in two steps as illustrated in
Figure 1c. In the first step, the square fiber (F) is combined with part of the matrix phase (M
1) to form a rectangular composite C
1, where the matrix block has the same width and height as the fibrous phase. In the second step, the composite C
1 is combined with the rest of the matrix to form the unit cell. The effective thermal and electrical conducting properties of the unit cell can be calculated using Kirchhoff’s law, and the effective Seebeck coefficient can be calculated using Thevenin’s theorem [
25]. The material properties of the composite C
1 can be derived according to Equation (1):
The symbols
ρ,
λ,
S, and
f stand for electrical resistivity, thermal conductivity, Seebeck coefficient, and fiber volume fraction. The subscripts p, a, u, v, w, x
i, y
i, and z
i stand for component material’s properties along each material’s local coordinate systems. The material properties of the unit cell can be derived according to Equation (2):
In
Figure 1, the fiber axis aligns parallel to the z
2-axis of the composite. When fibers are tilted aligned in the composite, the effective properties of the composite can be calculated through matrix transformation. For example, if the fibers in
Figure 1 are rotated by an angle of
θ around x
2-axis into the configuration in
Figure 2, the effective properties of the composite can be calculated as [
26]:
In a transverse thermoelectric material, heat flux and electrical current flow perpendicular to each other. Based on
Figure 2, we assume the electrical current flows in the y-direction and heat flux flows in the z-direction. The transverse thermoelectric figure of merit is defined as
[
10,
27], where
Szy is the transverse Seebeck coefficient,
ρyy is electrical resistivity in the y-direction,
λzz is thermal conductivity in the z-direction, T is the operating temperature. The term
Szy,
ρyy, and
λzz values can be calculated using Equation (1) to Equation (3).
Figure 2 can also represent a fibrous transverse thermoelectric device operating under cooling mode, where the top surface serves as a cooling surface with a temperature of T
c and the bottom surface serves as a heat sink with a temperature of T
h. When the electrical current flows in the y-direction, the transverse Peltier effect will trigger heat flux along the z-direction. The maximum cooling capacity of the system is related to the Z
transT value, according to previous studies [
10,
28]. In this study, the T
h was set to be the measuring temperature of material properties, which are shown in
Table 1.
In the next step, experimentally measured material properties are implemented into the derived equations. Previous studies [
9] have shown that high-performance transverse thermoelectric composites usually consist of one semiconducting thermoelectric phase and one highly conductive phase, where the thermoelectric phase provides high Seebeck coefficient and low thermal conductivity, the conducting phase provides low electrical resistivity. Therefore, we chose anisotropic Bi
2Te
3 polycrystal as the matrix phase and Copper as the fibrous phase for the composite. The properties of these materials are shown in
Table 1. According to
Table 1, the Copper phase was isotropic and Bi
2Te
3 was anisotropic with planar (u–v plane) isotropy.
When anisotropic material properties are implemented in the component phase of a composite, the alignment of component material’s local property coordinate system, with respect to the composite coordinate system, will affect the effective properties of the composite. In
Figure 1, the u–v plane of the matrix material aligns perpendicular to the fiber axis. If the material property coordinate in
Figure 1b is rotated around the u-axis 90 degrees, the anisotropy plane of Bi
2Te
3 will be parallel to the fiber axis direction. For the convenience of analysis, we shall refer to the two configurations mentioned above as configuration I and configuration II. In configuration I, the material property axes u, v, and w are parallel to the composite axes x
2, y
2, and z
2, respectively. In configuration II, the material property axes u, v, and w are parallel to the composite axes x
2, z
2, and y
2, respectively. The corresponding Z
transT and maximum cooling capacities (∆T
max) based on these two configurations were calculated using the derived mathematical model with respect to fiber rotation angle,
θ, and fiber volume fraction,
f. The results are shown in
Figure 3.
Based on
Figure 3, the Z
transT and ∆T
max values exhibited similar trends for both configurations. The maximum values appearred at rotation angles between 75 and 85 degrees. The change in fiber volume fraction under a constant rotation angle had a minor influence on the Z
transT and ∆T
max values. The maximum Z
transT and ∆T
max values in configuration I were 0.29 and 34 K, while the values of maximum Z
transT and ∆T
max in configuration II were 0.22 and 27 K. These results indicate a 26% difference in the peak thermoelectric performance among the two configurations. Hence, the alignment of the material’s local property coordinates with respect to the composite coordinate system had significant influence on the composite’s cooling performance.
Many thermoelectric materials are polycrystals. These polycrystals may exhibit certain degrees of texture, and thus, anisotropic properties as a result of the fabrication processes. There exists theoretical models that can correlate the aggregate properties of polycrystals with its single-crystal material property. Among these theoretical models, the Voigt model and Reuss model can provide upper and lower bounds for the aggregate properties of polycrystals, respectively [
32,
33,
34]. In this study, we take experimentally measured values from a single-crystal Bi
2Te
3 (
Table 1) and calculate possible aggregate material properties using Reuss and Voigt models. The Bi
2Te
3 polycrystal has isotropic Seebeck coefficients and exhibits planar isotropy (u–v plane) in electrical resistivity and thermal conductivity based on existing experimental results [
5,
18,
35,
36]. The anisotropy ratio terms
and
are used to relate the material properties in different axis directions, such that
,. By importing these hypothetical aggregate properties into the Bi
2Te
3 phase, the maximum cooling capacities were calculated for the Bi
2Te
3 matrix/Cu fiber composite. The results are shown in
Figure 4. According to the definition of
and
, the bottom-left corner of
Figure 4a,b represents the maximum cooling capacity of a polycrystal Bi
2Te
3 with isotropic material properties. The results in
Figure 4 show that by introducing material anisotropy in a fibrous transverse thermoelectric composite, the cooling capacity of the composite can be improved by 7% (
Figure 4b, Voigt model) and 21% (
Figure 4a, Reuss model).
3. Finite Element Simulation for Fibrous Thermoelectric Composite with Anisotropic Components
In order to validate the effectiveness of the mathematical model, finite element simulations were carried out using COMSOL (COMSOL Multiphysics, COMSOL Inc., Burlington, MA, USA). A unit cell model based on
Figure 1c was constructed using polycrystal Bi
2Te
3 as the matrix phase and Copper as the fiber phase assuming the isotropy plane of Bi
2Te
3 was perpendicular to the fiber axis. For all finite element models involved in this study, mesh convergence studies were carried out and the convergence error was limited to below 3%.
According to the unit cell model shown in
Figure 1c, square fiber was used for mathematical simplicity. In the finite element simulations, both cylindrical fiber and square fiber were investigated so as to explore the effect of fiber shapes. The effective properties of the unit cells are shown in
Figure 5.
According to
Figure 5, the finite element simulation results and mathematical model agreed well with each other. The conductivities of the unit cell increased as the volume fraction of the fibers increased, since the fiber was more conductive than the matrix material. In
Figure 5f, the effective Seebeck coefficient had negative values at a low fiber volume fraction, but became positive as the fiber volume fraction increased. This was because the Bi
2Te
3 phase had a negative Seebeck coefficient and Copper had a positive Seebeck coefficient. The effective Seebeck coefficient was dominated by the matrix phase at low fiber volume fractions, and by the fibrous phase at high fiber volume fractions.
Some minor deviations were observed between the mathematical model and the finite element simulation model. This was because the mathematical model assumes a one-dimensional flux flow, i.e., the electrical and heat flux only flow along the applied field, whereas fluxes in other directions are not considered. These fluxes are caused by the material inhomogeneity and the shape feature of fibers. When the fiber volume fraction is very small or very large, homogeneity of the composite will inevitably be compromised. Moreover, when the direction of the applied potential is not perpendicular to the interface between the matrix and fibrous phase, the secondary dimensional flux will occur. Both situations above are the cause of deviation between the mathematical model and finite element simulation results.
Finite element models with tilted aligned fibers in a matrix were also constructed to study the cooling capacity of the composite. The properties of polycrystal Bi
2Te
3 and Copper were applied. The fiber rotation angle was fixed at 80 degrees. According to
Figure 2, the composite can be seen as fiber arrays periodically aligned in the x-direction. Therefore, to simplify the finite element model, only one fiber array was constructed inside the matrix block, and periodic boundary conditions were applied on the surfaces parallel to the y-z plane. The geometry of the model was 2 mm × 500 mm × 10 mm, while the fiber diameter was adjusted to suit different fiber volume fractions. The bottom surface of the device was fixed at 300 K to serve as a heat sink, while the top surface was subjected to natural convective heat transfer. Two 1 mm thick Copper blocks were added on both ends of the y-direction of the composite to serve as electrodes. Twenty cases were studied with different fiber volume fractions, fiber shapes, and with respect to the two configurations mentioned above. Within each case, the input electrical current density was adjusted until the maximum temperature difference between the top and bottom surface was reached.
Figure 6 provides temperature distributions on the y–z plane for both cylindrical fiber model and square fiber model at 30% fiber volume fraction. Enlarged views of the top surface are shown in the insets. For both composites in
Figure 6, a temperature gradient caused by the transverse Peltier effect can be clearly observed in the vertical direction. In the horizontal direction, the temperature was uniformly distributed except at the ends of the electrodes. In the insets of
Figure 6, the abrupt temperature change at the interfaces of the fiber and matrix phases was the result of the Peltier effect.
The surface temperature was averaged on the top surface of each finite element model. The maximum temperature difference in each case was calculated by subtracting the heat sink temperature by the averaged cooling surface temperature. The results are presented in
Table 2. In general, the mathematical model results agreed well with the finite element simulation results with at most a 10% difference. These small discrepancies were mainly caused by the secondary dimensional flux that was previously discussed. The only exception in
Table 2 was the case for the cylindrical fiber with f = 0.1 under configuration I, where the ∆T
max calculated using the mathematical model was more than 10% higher than the finite element simulation result. This was because under small or large fiber volume fraction cases, the homogeneity of the composite was very poor, which conflicts with the assumption in the mathematical model. The differences between the two configurations also agreed with the results shown in
Figure 3, which emphasized the significant effect of the component material’s local property axes alignment in the composite coordinate system. During the finite element analysis, it was found that the size of the fibers has a slight influence on the TE performance. In anisotropic thermoelectric composites, smaller fibers provide the composite with better homogeneity. Although, under the same volume fraction of fibers the cooling capacity of the composite may not vary much with the size of the fibers, large fiber geometry is likely to cause higher temperature fluctuations in between the boundaries of the fibrous phase and the matrix phase.
4. Cooling Capacity Comparison for 1-Dimensional and 2-Demensional Inclusion Thermoelectric Composites with Anisotropic Components
The cooling performance of fibrous transverse thermoelectric composites were investigated with both a mathematical approach and a simulation approach. While the fibrous transverse thermoelectric composite can be treated as a 1D inclusion composite, the layered transverse thermoelectric composite can be treated as a 2D inclusion composite. The cooling capacity of the layered transverse thermoelectric composites with anisotropic components have been investigated thoroughly in previous studies [
21]. Hence, it is informative to compare the cooling performances of transverse thermoelectric composites between 1D and 2D inclusion composites. In this section, three types of anisotropic thermoelectric crystals were chosen as the semiconducting phase of the composite, while Copper was chosen as the conducting phase of the composite. The material properties are listed in
Table 1. The thermoelectric properties of Bi
2Te
3 single crystal were measured at 300 K, the thermoelectric properties of In
4Se
2.25 single crystal were measured at 600 K, and the thermoelectric properties of SnSe single crystal were measured at 700 K. The Copper served as the fibrous phase in the fibrous composite thermoelectric materials; its thermoelectric properties at 600 K and 700 K are also given in
Table 1. It should be noticed that, in transverse thermoelectric composites, the transverse ZT matrix is anisotropic. The study made use of engineered anisotropy in composite materials and optimized the transverse ZT with respect to one specific element in the anisotropic ZT matrix. The anisotropic ZT values for other elements in the matrix are not discussed.
The semiconducting materials used in this comparison, i.e., Bi
2Te
3, In
4Se
2.25, and SnSe, all belong to planar isotropic material. Therefore, for both 1D and 2D inclusion composites consisting of semiconducting crystal and isotropic Copper, two representative configurations were used. Assuming the local material property axis as u
iv
iw
i (where the subscript i differentiates the type of materials), and the isotropy plane of material lies in the u
iv
i plane. For fibrous composite materials, the definition of the two configurations were explained in previous sections of this study. For layered composite materials, in configuration III, the u
iv
i plane aligns parallel to the layered plane of materials; in configuration IV, the u
iv
i plane aligns perpendicular to the planar surface of the layers. Besides the four configurations, the conducting properties of isotropic polycrystal Bi
2Te
3 and SnSe, listed in
Table 1, can be calculated using Reuss and Voigt models. Since Bi
2Te
3 and SnSe single crystals have isotropic Seebeck properties, their isotropic polycrystal Seebeck properties are the same as anisotropic single crystals. The conducting properties of the isotropic polycrystals are the average of the results calculated by the Reuss and Voigt models, which is referred as ‘Average RV’ in
Table 3, respectively.
The results in
Table 3 present both maximum Z
transT as well as the maximum cooling capacity (∆T
max) for composites with different anisotropic components and configurations using the derived mathematical model. It can be seen that, for composites which consist of the same combination of component materials but under different configurations, the alignment of component material’s property axis in the composite has a major influence on Z
transT as well as the ∆T
max values. For Bi
2Te
3/Copper composite, the difference in thermoelectric performances caused by the alignment of the anisotropic component’s property axis was 50%. For the SnSe/Copper composite, this difference exceeded 400%. This difference was caused by the variance of properties among different material property axes in the anisotropic component phase, where larger variance usually leads to bigger variance in the Z
transT and ∆T
max values of a composite.
When the cooling performance of the transverse thermoelectric was compared under the same component combination but with different inclusion types, it was found that for the composite with the same component volume fraction, the fibrous and layered composites yielded very similar ∆T
max as well as Z
transT values. However, it should be noted that the peak performances of the layered and fibrous composites under the same volume fraction were not identical to each other. The ∆T
max values were rounded up to the nearest integer, and the maximum Z
transT values were rounded up to the second decimal in
Table 3. Although according to
Table 3, the Z
transT and ∆T
max values in between different configurations have same values, it should be noted that none of the results in
Table 3 are exactly identical to one another.
By comparing the results between the composite with anisotropic component material and isotropic component material, it was found that the anisotropic component in single-crystal form in both the layered transverse thermoelectric composite and the fibrous transverse thermoelectric composite can provide a better cooling performance compared to their isotropic polycrystal counterpart. For the Bi2Te3/Copper composite at 300 K, the composite with the anisotropic Bi2Te3 single crystal improved the maximum ZtransT and ∆Tmax by as much as 15% compared to the composite with the isotropic Bi2Te3 polycrystal. For the SnSe/Copper composite at 700 K, the improvement in the maximum ZtransT and ∆Tmax for composite with SnSe single crystal component was as much as 51% compared to composite with isotropic SnSe polycrystal components. For the In4Se2.25/Copper composite, the single-crystal form of In4Se2.25 had anisotropic Seebeck properties, but there are no existing studies regarding the effective properties of polycrystals with anisotropic Seebeck properties. Therefore, the comparison on the cooling performances for the In4Se2.25/Copper composite between composites with isotropic and anisotropic In4Se2.25 phase properties was not presented.