5. Numerical examples

Ð

280 Bringing Thermoelectricity into Reality

Ð

of domain Ω<sup>e</sup> of boundary Γ<sup>e</sup> become: R<sup>V</sup> <sup>b</sup> <sup>¼</sup> <sup>Ð</sup>

> R<sup>T</sup> <sup>b</sup> <sup>¼</sup> <sup>Ð</sup>

Ω<sup>e</sup> BΤ

<sup>Ω</sup><sup>e</sup> <sup>B</sup><sup>Τ</sup>

(a) The time interval is divided into small time increments Δt.

<sup>k</sup> ¼ �∂R<sup>a</sup> ∂gb � � � � k dgk <sup>b</sup> ) �

as backward or forward finite differences, see [19].

Raj

capacity matrix C are derived at each iteration as:

ab ¼ � <sup>∂</sup>R<sup>V</sup>

ab ¼ � <sup>∂</sup>R<sup>V</sup>

ab ¼ � <sup>∂</sup>R<sup>T</sup>

ab ¼ � <sup>∂</sup>R<sup>T</sup>

ab ¼ � <sup>∂</sup>R<sup>T</sup>

a ∂V~ <sup>b</sup>

a ∂T~<sup>b</sup>

a ∂V~ <sup>b</sup>

a ∂T~<sup>b</sup>

a ∂~\_ Tb ¼ � ð Ω<sup>e</sup>

¼ � ð Ω<sup>e</sup> BΤ a ∂j ∂V~ <sup>b</sup>

¼ � ð Ω<sup>e</sup> BΤ a ∂j ∂T~<sup>b</sup>

¼ � ð Ω<sup>e</sup> BΤ a ∂q ∂V~ <sup>b</sup>

¼ � ð Ω<sup>e</sup> BΤ a ∂q ∂T~<sup>b</sup>

KVV

KVT

KTV

KTT

CTT

spatial derivatives:

involves three steps:

<sup>Ω</sup>∇δV � jdΩ � ∮ <sup>Γ</sup>δVj dΓ ¼ 0,

v) In this step, similar discretization to those of (10) are applied to the test functions and their

Introducing these expressions in the weak form of (11), the FE residuals for each finite element

N <sup>b</sup> j dΓe,

<sup>b</sup> <sup>þ</sup> <sup>N</sup> <sup>b</sup> <sup>j</sup> � <sup>∇</sup><sup>V</sup> <sup>þ</sup> <sup>r</sup><sup>m</sup> <sup>c</sup> <sup>T</sup>\_ � � � � <sup>d</sup>Ω<sup>e</sup> � <sup>∮</sup> <sup>Γ</sup><sup>e</sup>

vi) Finally, the solution is calculated by solving a set of nonlinear transient equations that

(b) The time derivatives are replaced by discrete forms using time integration techniques such

(c) Since the thermoelectric problem is nonlinear due to both the temperature-dependency of the material properties and the Joule heating, the resulting nonlinear algebraic problem is linearized at each time increment by using a numerical algorithm. For instance, the Newton-

> ∂R<sup>a</sup> ∂gb � � � � k

where a, b refer the local numbering of two generic FE nodes, the parameters c<sup>1</sup> and c<sup>2</sup> depend on the time integration scheme, see [19], and the consistent tangent matrices K and the nonzero

dΩe,

dΩe,

þ N <sup>a</sup>

þ N <sup>a</sup>

N <sup>a</sup> r<sup>m</sup> c N <sup>b</sup> dΩe,

where the derivatives are reported in [24] and are not repeated here for the sake of brevity.

� �

� �

<sup>∂</sup> <sup>j</sup> � <sup>∇</sup><sup>V</sup> � � ∂V~ <sup>b</sup>

<sup>∂</sup> <sup>j</sup> � <sup>∇</sup><sup>V</sup> � � ∂V~ <sup>b</sup>

dΩe,

(15)

dΩe,

δV ≈ N <sup>a</sup>δV~ a, δT ≈ N <sup>a</sup>δT~a,

<sup>b</sup> j dΩ<sup>e</sup> � ∮ <sup>Γ</sup><sup>e</sup>

Raphson algorithm use k iterations for the linearization of the residuals:

<sup>Ω</sup> <sup>∇</sup>δ<sup>T</sup> � <sup>q</sup> <sup>þ</sup> <sup>δ</sup><sup>T</sup> <sup>j</sup> � <sup>∇</sup><sup>V</sup> <sup>þ</sup> <sup>r</sup><sup>m</sup> <sup>c</sup> <sup>T</sup>\_ � � � � <sup>d</sup><sup>Ω</sup> � <sup>∮</sup> <sup>Γ</sup>δ<sup>T</sup> q d<sup>Γ</sup> <sup>¼</sup> <sup>0</sup>: (11)

<sup>∇</sup>δ<sup>V</sup> <sup>≈</sup> <sup>∇</sup><sup>N</sup> <sup>a</sup>δV<sup>~</sup> <sup>a</sup> <sup>¼</sup> <sup>B</sup>aδV<sup>~</sup> a, <sup>∇</sup>δ<sup>T</sup> <sup>≈</sup> <sup>∇</sup><sup>N</sup> <sup>a</sup>δT~<sup>a</sup> <sup>¼</sup> <sup>B</sup>aδT~a: (12)

<sup>N</sup> <sup>b</sup> <sup>q</sup> <sup>d</sup>Γe, (13)

¼ c<sup>1</sup> Kab þ c<sup>2</sup> Cab, (14)

As commented, the present numerical tool can be used as a "virtual laboratory" that allows to numerically design and optimize devices made out of thermoelectric materials. In this connection, the aim of this section is to show several features of the FE code. In particular, a commercial Peltier cell for cooling applications is simulated and several variables such as heat extracted, voltage and temperature distributions are analyzed.

A CP1.4-127-045 thermoelectric cooling module composed of 127 thermocouples electrically connected in series, as observed in Figure 4 (left), and manufactured by LairdTech [32] is modeled. The technical specifications of this device are: maximum intensity I ¼ 8:7 (A), maximum extracted heat Qc <sup>¼</sup> <sup>82</sup>:01 (W) with voltage drop <sup>V</sup> = 15.33 (V) and under Th <sup>¼</sup> Tc <sup>¼</sup> <sup>50</sup><sup>∘</sup> C, where Tc and Th refer to the temperature at cold and hot ends, respectively.

Numerically, only half of the thermocouple needs to be modeled due to its symmetry and, consequently, the computational time is reduced; the FE mesh composed of 12,670 eight noded elements is shown in Figure 4 (right). On the one hand, the Neumann boundary conditions, namely, the electric current is enforced by the specific two-dimensional FE developed in [21]. This element can be also used to take into account convection and radiation phenomena among thermocouples; however, both phenomena are ignored in the present chapter. On the other hand, the Dirichlet boundary conditions are applied by setting the voltage reference

Figure 4. Left: Sketch of the thermoelectric cooling module. Right: Finite element mesh and boundary conditions applied to half of the thermocouple.

V ¼ 0 and by prescribing temperatures at cold and hot ends, see Figure 4 (right). Finally, the coefficients of the temperature-dependent material properties of (7) are also reported in [21].

As observed, analytical and numerical solutions for the voltage disagree due to the electrical loss at the corners of the thermoelements. This drop loss is captured by the FE but not by the one-dimensional analytical solution, as was concluded in [21] by comparing analytical, numerical and experimental solutions. Therefore, for a proper design of thermoelectric cells the numerical tool is more appropriated than the simple analytical solutions. For the extracted heat, despite the fact that the slope is slightly different, both analytical and numerical solutions

Computational Thermoelectricity Applied to Cooling Devices

http://dx.doi.org/10.5772/intechopen.75473

283

Regarding the drop loss at the corners of the thermoelement, Figure 7 shows contour plots of the electric currents along horizontal (left) and vertical (right) axes. As observed, the fluxes produce vortices-like phenomena and, consequently, part of the prescribed intensity does not

Figure 6. Voltage drop (left) and extracted heat (right) vs. temperature at the cold end for <sup>I</sup> <sup>¼</sup> <sup>1</sup>:7 (A), Th <sup>¼</sup> <sup>50</sup><sup>∘</sup> C.

Figure 7. Contour plots of electric current along horizontal (left) and vertical (right) axes for <sup>I</sup> <sup>¼</sup> <sup>1</sup>:7 (A), Th <sup>¼</sup> Tc <sup>¼</sup> <sup>50</sup><sup>∘</sup> C.

agree.

go to the thermoelement.

#### 5.1. Convergence tests

In order to study the accuracy of the FE simulation and since the problem is nonlinear, convergence tests must be carried out. In particular, this section presents two tests: h and residual convergences. The former refers to the improvement of the solution with decreasing the size of the finite elements—this convergence is called h-refinement by the computational mechanics community—and the latter guarantees a proper performance of the Newton-Raphson algorithm to solve the nonlinearities.

For this purpose, the most detrimental operating condition is simulated: the boundary conditions are Tc <sup>¼</sup> Th <sup>¼</sup> <sup>50</sup><sup>∘</sup> C and I ¼ 8:7 (A). Figure 5 shows the h-refinement (left) and residual convergence (right) for the potential drop V (top) and the extracted heat at the cold end Qc (bottom).

As observed in the left figure, a mesh composed of 12,670 elements gives a perfect convergence for both electrical and thermal variables. The potential drop converges faster than the extracted heat; due to the fact that the voltage is directly calculated while the heat depends on spatial derivatives. Regarding the residual convergence and as observed in the right figure, both electrical and thermal variables require three iterations of the Newton-Raphson algorithm (the residual norm becomes 10�16) to minimize the residual due to the relative nonlinearity of the problem. Notice that each iteration requires 0.46 (s) in a 1.6 (GHz) Intel Core i5 processor.

#### 5.2. Comparisons of extracted heat and voltage drop

The aim of this section is to compare analytical and FE solutions for the commercial cell under <sup>I</sup> <sup>¼</sup> <sup>1</sup>:7 (A) and Th <sup>¼</sup> <sup>50</sup><sup>∘</sup> C. The analytical solutions are obtained from [33] and are typically used by designers. For this aim, Figure 6 shows the voltage drop (left) and the extracted heat (right) versus the temperature at the cold end.

Figure 5. Left: mesh convergence, voltage drop (top) and heat extracted (bottom) vs. number of elements. Right: residual convergence, residual norm of voltage (top) and heat extracted (bottom) vs. number of iterations.

As observed, analytical and numerical solutions for the voltage disagree due to the electrical loss at the corners of the thermoelements. This drop loss is captured by the FE but not by the one-dimensional analytical solution, as was concluded in [21] by comparing analytical, numerical and experimental solutions. Therefore, for a proper design of thermoelectric cells the numerical tool is more appropriated than the simple analytical solutions. For the extracted heat, despite the fact that the slope is slightly different, both analytical and numerical solutions agree.

V ¼ 0 and by prescribing temperatures at cold and hot ends, see Figure 4 (right). Finally, the coefficients of the temperature-dependent material properties of (7) are also reported in [21].

In order to study the accuracy of the FE simulation and since the problem is nonlinear, convergence tests must be carried out. In particular, this section presents two tests: h and residual convergences. The former refers to the improvement of the solution with decreasing the size of the finite elements—this convergence is called h-refinement by the computational mechanics community—and the latter guarantees a proper performance of the Newton-

For this purpose, the most detrimental operating condition is simulated: the boundary conditions are Tc <sup>¼</sup> Th <sup>¼</sup> <sup>50</sup><sup>∘</sup> C and <sup>I</sup> <sup>¼</sup> <sup>8</sup>:7 (A). Figure 5 shows the h-refinement (left) and residual convergence (right) for the potential drop V (top) and the extracted heat at the cold end Qc

As observed in the left figure, a mesh composed of 12,670 elements gives a perfect convergence for both electrical and thermal variables. The potential drop converges faster than the extracted heat; due to the fact that the voltage is directly calculated while the heat depends on spatial derivatives. Regarding the residual convergence and as observed in the right figure, both electrical and thermal variables require three iterations of the Newton-Raphson algorithm (the residual norm becomes 10�16) to minimize the residual due to the relative nonlinearity of the problem. Notice that each iteration requires 0.46 (s) in a 1.6 (GHz) Intel Core i5 processor.

The aim of this section is to compare analytical and FE solutions for the commercial cell under

used by designers. For this aim, Figure 6 shows the voltage drop (left) and the extracted heat

Figure 5. Left: mesh convergence, voltage drop (top) and heat extracted (bottom) vs. number of elements. Right: residual

convergence, residual norm of voltage (top) and heat extracted (bottom) vs. number of iterations.

C. The analytical solutions are obtained from [33] and are typically

5.1. Convergence tests

282 Bringing Thermoelectricity into Reality

(bottom).

<sup>I</sup> <sup>¼</sup> <sup>1</sup>:7 (A) and Th <sup>¼</sup> <sup>50</sup><sup>∘</sup>

Raphson algorithm to solve the nonlinearities.

5.2. Comparisons of extracted heat and voltage drop

(right) versus the temperature at the cold end.

Regarding the drop loss at the corners of the thermoelement, Figure 7 shows contour plots of the electric currents along horizontal (left) and vertical (right) axes. As observed, the fluxes produce vortices-like phenomena and, consequently, part of the prescribed intensity does not go to the thermoelement.

Figure 6. Voltage drop (left) and extracted heat (right) vs. temperature at the cold end for <sup>I</sup> <sup>¼</sup> <sup>1</sup>:7 (A), Th <sup>¼</sup> <sup>50</sup><sup>∘</sup> C.

Figure 7. Contour plots of electric current along horizontal (left) and vertical (right) axes for <sup>I</sup> <sup>¼</sup> <sup>1</sup>:7 (A), Th <sup>¼</sup> Tc <sup>¼</sup> <sup>50</sup><sup>∘</sup> C.

Another important aspect to be considered for the proper operation of Peltier cells is the maximum temperature inside the thermoelements since an overheating produces thermal stresses that could break the semiconductors from a mechanical point of view. For that matter, Figure 8 shows the maximum and minimum temperature inside the thermoelements for two operating conditions: <sup>I</sup> <sup>¼</sup> <sup>1</sup>:7 and <sup>I</sup> <sup>¼</sup> 3 (A), at Th <sup>¼</sup> <sup>50</sup><sup>∘</sup> C.

Figure 9 shows the voltage drop (left) and the extracted heat (right) versus the thermoelement length for the detrimental case <sup>I</sup> <sup>¼</sup> <sup>8</sup>:7 (A) and Tc <sup>¼</sup> Th <sup>¼</sup> <sup>50</sup><sup>∘</sup> C. As observed, decreasing the size of thermocouple implies to reduce the heat extracted and increase the potential drop. Obviously, both results are detrimental since the coefficient of performance of the cell becomes worse. This fact is due to the Joule heating that is a volumetric effect and, consequently, decreasing the size and maintaining constant the prescribed intensity results in overheating of the thermoelements. Again, the present modeling tool could be used for calibrating the miniaturization of devices by setting the applied intensity in order to achieve an efficient functioning.

Computational Thermoelectricity Applied to Cooling Devices

http://dx.doi.org/10.5772/intechopen.75473

285

This chapter has presented a thermodynamically consistent formulation to obtain the governing equations of thermoelectricity, from a theoretical point of view. Then, a nonlinear numerical formulation within the finite element method is developed. The nonlinearities emerge from the presence of Joule heating and the temperature-dependence of the material properties and they are numerically solved by the Newton-Raphson algorithm. Notice that this material nonlinearity directly allows the inclusion of the Thomson effect. Furthermore, the formulation is dynamic and monolithic; the former feature is solved by backward finite differences and the latter is carried out by defining a coupled assembled matrix that increases the accuracy of the formulation. Finally, several examples are reported to show the capabilities of the numerical formulation that can be used as a "virtual laboratory." In particular, hrefinements and residual convergence tests are conducted to validate the codification. Then, comparisons between analytical and numerical solutions for cooling thermoelectric cells are reported in order to highlight the advantages of the simulations against the simple onedimensional analytical solutions. In conclusion, the use of the present numerical tool could be applied for a proper design and optimization of thermoelectric devices. For instance, this tool

was used to optimize the shape of a pulsed thermoelectric in [25].

Department of Mechanical Engineering and Construction, Universitat Jaume I, Spain

[1] Dresselhaus MS, Chen G, Tang MY, Yang RG, Lee HL, Wang DZ, Ren ZF, Fleurial JP, Gogna P. New directions for low-dimensional thermoelectric materials. Advanced Mate-

Roberto Palma\*, Emma Moliner and Josep Forner-Escrig

\*Address all correspondence to: rpalma@uji.es

rials. 2007;19:1043-1053

6. Concluding remarks

Author details

References

As observed, for <sup>I</sup> <sup>¼</sup> 3 (A) the maximum temperature becomes approximately 80<sup>∘</sup> C due to the increasing of the Joule heating that depends on the prescribed intensity. To sum up, the present numerical tool can be used for a proper design and optimization of cooling devices from both thermoelectric and thermomechanic interactions.

#### 5.3. Miniaturization of thermoelements

The last example deals with the minimization of thermoelements. Currently, there exists a trend towards minimization of thermoelectric devices for two reasons: decreasing power consumption and reducing the size of the devices.

Figure 8. Minimum and maximum temperature inside thermoelements vs. temperature at the cold end for I ¼ 1:7 (left) and <sup>I</sup> <sup>¼</sup> 3 (A) (right), Th <sup>¼</sup> <sup>50</sup><sup>∘</sup> C.

Figure 9. Voltage drop (left) and extracted heat (right) vs. the thermoelement length.

Figure 9 shows the voltage drop (left) and the extracted heat (right) versus the thermoelement length for the detrimental case <sup>I</sup> <sup>¼</sup> <sup>8</sup>:7 (A) and Tc <sup>¼</sup> Th <sup>¼</sup> <sup>50</sup><sup>∘</sup> C. As observed, decreasing the size of thermocouple implies to reduce the heat extracted and increase the potential drop. Obviously, both results are detrimental since the coefficient of performance of the cell becomes worse. This fact is due to the Joule heating that is a volumetric effect and, consequently, decreasing the size and maintaining constant the prescribed intensity results in overheating of the thermoelements. Again, the present modeling tool could be used for calibrating the miniaturization of devices by setting the applied intensity in order to achieve an efficient functioning.
