Simplified thermoelectric generator (TEG) with heatsinks modeling and simulation using Matlab and Simulink based-on dimensional analysis

Abbreviations and Symbols: a: TEG p-n junction thermocouple area in m²; A₁: TEG hot side heatsink1 total surface area in m²; A₂: TEG cold side heatsink2 total surface area in m²; Ab: Heatsink base area in m²; ΔT: TEG(s) temperature difference (T_h – T_c) in ℃ or kelvin; DT_s: TEG dimensionless temperature difference; Eff*: TEG dimensionless conversion efficiency; h₁: TEG hot side heatsink1 convection coefficient; h₂: TEG cold side heatsink2 convection coefficient; HFD: TEG heat flux density in W/m²; HS: Heatsinks; I: TEGs output current in ampere through the TEG; k: TEG/TEC effective thermal conductivity in W/mK; K: TEG/TEC thermal conductance in (W/K); K₁: TEG hot side heatsink1 thermal conductance; K₂: TEG cold side heatsink2 thermal conductance; L: TEG/TEC p-n junction thermocouple length in meter; LHS: Left Hand Side; n: TEG/TEC manufacturer p-n thermocouples amount used in a TEG/TEC; N: TEG number of modules required; N_h: TEG dimensionless convection conductance; N_i: TEG-HS dimensionless output current; N_k: TEG dimensionless thermal conductance; N_v: TEG-HS dimensionless output voltage (V_os); ɳ: TEG thermal/electrical/conversion efficiency; Ƞ₁: TEG hot side heatsink1 fin efficiency; Ƞ₂: TEG cold side heatsink2 fin efficiency; ɳhA: TEG convection conductance; ɳ₁h₁A₁: TEG hot side heatsink1 convection conductance; ɳ₂h₂A₂: TEG cold side heatsink2 convection conductance; ρ: TEG/TEC electrical resistivity in Ωm; ρ_e: TEG/TEC effective electrical resistivity in Ωm; P_n: TEG normalized output power; P_o: TEG output power in watt, is the difference between Q_h and Q_c; P_os: TEG dimensionless output power, is the difference between Q_h*(Q_s1) and Q_c*(Q_s2); Q_c: TEG heat emitted on its cold side in watt; Q_h: TEG heat absorbed on its hot side in watt; Q_s1: TEG hot side dimensionless heat absorbed (Q_h*); Q_s2: TEG cold side dimensionless heat released (Q_c*); $\mathrm{r}$: TEG/TEC thermocouples p-n junction unit resistance in ohm; RHS: Right Hand Side; R_L: TEG electrical load resistance in Ω connected to the TEG output; R_t: TEG/TEC module internal resistance in ohm—also denoted as R_R in Figure 2; R_r: TEG dimensionless internal electrical resistance—also denoted as rr in Figure 2; S: TE device Seebeck coefficient in V/K; S_e: TEG/TEC effective Seebeck coefficient in V/K; $\overline{\mathrm{T}}$: TE device average temperature (T_h + T_c)/2 in K or ℃; T₁: TEG hot side temperature (T_h); T₂: TEG cold side temperature (T_c); TE: Thermoelectric; TEC: Thermoelectric cooler; TEG: Thermoelectric generator; T_i: TEG hetasinks fluid dimensionless temperatures (Tis); T_i1: TEG hot side heatsink1 fluid dimensionless temperature; T_i2: TEG cold side heatsink2 fluid dimensionless temperature; T_s1: TEG hot side dimensionless temperature; T_s2: TEG cold side dimensionless temperature; T_c: Temperature on TEG/TEC cold side in ℃; T_h: Temperature on TEG/TEC hot side in ℃; T_p: TEG/TEC modules quantity in parallel; T_s: TEG/TEC modules quantity in series; V_o: TEG module output voltage in volt; V_oc: TEG ideal or open circuit output voltage (consider it as EMF = electromotive force) in volt; V_os: TEG dimensionless output voltage (N_v); V_n: TEG normalized output voltage; Z: TE device figure of merit in per K; Z_e: TEG/TEC effective figure of merit in per K; $\mathrm{Z}\overline{\mathrm{T}}$: TE device average dimensionless figure of merit; ZT_i2: TEG dimensionless figure of merit at temperature T_i2; ZTA: TEG dimensionless figure of merit at temperature TA ($\overline{\mathrm{T}}$)

1.   Introduction

Indicated in ^[1], South Africa has pledged its commitments to become a carbon reduced and green economic system, by considering a mixture of energy sources to secure energy sustainability. Upon ratifying the Paris agreement on climate change, renewable energy and energy efficiency frameworks as well as the regulatory policies have been developed and continuously fine-tuned taking into cognizance the current energy dynamics. In this regard, there is necessity and increasing demand for renewable energy to supplement and stabilize the unstable national grid, as well as for private use. These circumstances, warrant our research for an alternative energy rooted in thermoelectricity, with focus on basic residential energy sources and energy efficient loads. Basically, TEGs convert heat to DC power, whereas their dual thermoelectric coolers (TECs), reversibly converts DC power to cold and heat depending on the applied voltage polarity. If TEGs and TECs are properly designed, both can be relatively energy efficient and helpful for essential domestic energy use such as low-voltage DC power, lighting, heating and cooling. Theoretical frameworks and applications of both TEG and TEC were comprehensively presented in ^[2,3,4,5]; however, they fell short of some practical aspects—which this article seeks to advance as articulated next.

In ^{[4,5,6,7,8,9]}, thermoelectricity or thermoelectric (TE) devices (TEGs and TECs), require heatsinks or heat-exchangers to physically and reliably function and maintain a working temperature difference on the TEG and TEC hot and cold sides and most significantly to discharge the internal heat resulting from Joule (Ohmic) heating cause by the current flowing through the TEG/TEC. This Ohmic heating is irreversible and will become excessive if not properly managed and consequently, will cause drop in the output/cooling powers and inefficiency (change in entropy-heat flow will change direction due to increase in the cold side junction temperature). This can as well damage the TEG and TEC as a result of over-heating exceeding the hot and cold sides maximum temperature limit—which will melt the thermoelectric devices p-n junction thermocouples solder joints. The inclusion of heatsinks as a solution to alleviate this over-heating practical limitations in thermoelectric devices, unfortunately as well add/increase the TEGs and TECs hot and cold sides thermal resistances, which then oppose the ideal heat flow, hence making the TE devices inefficient.

To circumvent this thermal resistance problem resulting from adding heatsinks/heat-exchangers, various approaches have been investigated and some are stated briefly as follows. Simplified in ^[6] is a dimensional analysis by converting thermal resistance to convection conductance when using TEC with heatsinks—which was researched originally in ^[7]; however, the presentation lacks key analytical formulas. In ^[8], a TEG with heatsink for waste energy harvesting was analysed and optimised as a Thevenin equivalent circuit and investigated in ^[9] is a dimensionless model of a TEC functioning at real heat transfer state. Examined in ^[10], a multi-physics simulations were conducted to examine the thermal and electrical performances of a TEG module sand-witched between hot and cold blocks, whereas in ^[11], a parametric thermal analysis of TEG performance was done using dimensional analysis and presented in ^[12] is a numerical simulation and performance analysis of a TEC based-on the lattice Boltzmann method. Studied in ^[13] is a multi-parameters analysis and optimization of a typical TEC using dimensional analysis. Proposed in ^[14] is a thermal resistance matching for TEC systems whereas examined in ^[15], is an optimal design of a multi-couple TEG, meanwhile ^[16] employed a graphical approach to design TEC systems. Demonstrated in ^[17] is an electron-transparent TEC using nano-particles and condensation thermometry and investigated in ^[18] is an optimal heat-exchanger in different automobile exhaust temperatures for TEG system using dimensional analysis. From the reviewed studies, the dimensional analysis technique, especially the unique and practical approach presented in ^[7] to improve TE devices efficiency when used with heatsinks, is of interest. This paper with keen focus on TEG as exemplified in Figure 1, further develops the dimensional analysis, by now employing a simplified Matlab implementation to derive novel analytical formulas to directly calculate a TEG hot and cold sides dimensionless temperature.

Our study was structured around the following points. We begin with the TEG with heatsinks applicable maths developed from ^[6,7] with now emphasis on the new additional analytical formulas. We derived the TEG hot and cold sides dimensionless temperature, as well as the further optimal simplifications contributed relative to the shortcomings of the formulas/approach used in ^[7]. In a next step, we depict the TEG with heatsinks simulator easy user interface modeled with Matlab/Simulink. The presentation of our results is followed by discussions as well as validations and lastly we conclude our study with a summary.

2.   Dimensional analysis mathematics

The relevant mathematics for a TEG with heatsinks is expressed herein and we then introduce the new dimensionless temperature difference (DTs)—which is employed to simplify the apt approach given by ^[7]. Dimensional analysis is a technique that enables parameters with the same unit to be normalized within a minimum and maximum value, thus making them dimensionless and easier to work with, without worrying about their dimensions or unit of measurement. Nonetheless, novel in ^[7] is the conversion of TEG heatsinks thermal resistance to their convection conductance. This optimization technique gets rid of the TEG heatsinks thermal resistance which are harder to work with, in favour of the TEG heatsinks fluid convection conductance which are related to the TEG heatsinks fluid temperature and physically simpler to work with. This approach is pragmatic and the mathematical analysis is articulated in what follows.

TEG with heatsinks general heat flow equations:

Herein, the TEG heat flow equations, the heatsinks thermal resistance, their corresponding thermal conductance and convection conductance relationships developed from ^[6] are constituted as follows:

With Q₁ and Q₂ being respectively the heat flow rates on the TEG heatsinks hot and cold sides, K₁ and K₂ are respectively the TEG hot and cold sides heatsinks thermal conductance, T_i1 and T_i2 are respectively the temperatures of the heatsinks fluid on the TEG hot and cold sides and lastly T₁ and T₂ are respectively the TEG hot and cold sides p-n junction temperatures. It is worthy of note that thermal resistance is the reciprocal of thermal conductance K—which corresponds to the convection conductance (ȠhA) and thus, Eqs 1 and 2 in respect of the convection conductance, can be re-written as:

With Ƞ₁ being the fin efficiency of heatsink1 (on the TEG hot-side), h₁ is the convection coefficient of heatsink1 and A₁ is the total surface area of heatsink1. Similarly, Ƞ₂ is the fin efficiency of heatsink2 (on the TEG cold-side), h₂ is the convection coefficient of heatsink2 and lastly A₂ is the total surface area of heatsink2.

The TEG standard ideal heat flow equations are defined as:

where Q_h is the heat absorbed on the TEG hot side, n is the number of p-n junction thermocouples used in the TEG, S being the Seebeck coefficient, I being the output current from the TEG, K is the thermal conductance (computed as ak/L, with a being the area of the TEG p-n junction thermocouple, k is the TEG thermocouple p-n junction thermal conductivity and L is the TEG thermocouple p-n junction length). ΔT = T₁ − T₂ is the temperature difference between the TEG hot and cold sides, R = nr is the TEG module resistance, with r being the resistance of the TEG thermocouple p-n junction and Q_c is the heat released on the TEG cold side. The three terms on the right side of Eqs 5 and 6 are the Seebeck, Fourier and Ohmic terms respectively; with S, K and R considered as temperature invariant. Now, taking into account the energy balance of the TEG with heatsinks system, Eqs 1, 3 and 5 are respectively equivalent now to Eqs 2, 4 and 6—which boils down to (with parameters T₁, T₂ and I being the unknowns):

With recognition to ^[7], TEG with heatsinks (HS) is optimised by defining the dimensionless parameters with regards to fluid 2 (water or air on the TEG cold side of heatsink2) and because the optimization is rendered dimensionless relative to fluid 2, fluid 2 temperature (T_i2) and convection conductance must be given initially. The following dimensionless maths is developed further from ^[6,7] as:

TEG-HS dimensionless thermal conductance (N_k):

This is the ratio of the thermal conductance K and the convection conductance ȠhA in fluid 2, deduced as:

TEG-HS dimensionless convection (N_h):

This is the ratio of fluid 1 and fluid 2 convection conductances, denoted as:

TEG-HS dimensionless current (N_i):

TEG dimensionless temperatures (T_s1, T_s2, T_i and DT_s):

TEG dimensionless temperatures are presented as:

TEG-HS dimensionless heat absorbed (Q_s1):

TEG-HS dimensionless heat released (Q_s2):

TEG dimensionless output power (P_out* or P_os):

where P_out is the TEG output power (power delivered to the electrical load)

TEG dimensionless output voltage (V_os or N_v):

TEG dimensionless conversion efficiency (E_ff^*):

NB: TEG conversion efficiency by default is dimensionless. The mentioned of dimensionless in front conversion efficiency is just for emphasis on the dimensionless analysis technique used.

TEG-HS dimensionless heat absorbed (Q_s1) in terms of T_s1:

TEG-HS dimensionless heat released (Q_s2) in terms of T_s2:

TEG-HS dimensionless internal electrical resistance (R_r):

where R_L is the electrical load connected to the TEG output and R_t is the TEG module internal electrical resistance. NB: R_r is denoted as rr in the TEG-HS simulator in Figure 2.

As per ^[7], using the dimensionless Eqs 9–14, Eqs 3–6 reduce to the following two expressions (26) and (27) having five unknowns that must be found for T_s1 and T_s2 in terms of five independent dimensionless parameters of ZT_i2 (dimensionless figure of merit at T_i2), T_i, N_h, N_k and R_r.

As noticeable, Eqs 26 and 27 are cumbersome and unsolved further in ^[7] in terms of T_s1 and T_s2; therefore, there are no exact analytical equations to directly compute T_s1 and T_s2. Eqs 26 and 27 are awkward closed-form expressions of T_s1 and T_s2, which can only be solved by using numerical analysis; as a result, a numerical method using iterations, tables, graphs and approximations were employed by ^[7], as well as using a computer programme (NEDO) was further recommended. Consequently, because Eqs 26 and 27 could not be simplified further by ^[7], T_s1 and T_s2 could only be expressed as in Eqs 28 and 29 as functions of ZT_i2, T_i, N_h, N_k and R_r for solving numerically.

This entire numerical process is tedious to determine T_s1 and T_s2; thus, this limitation motivated our study to seek for a direct analytical and better solution(s) which are asserted next.

We now introduce Eq 15—the dimensionless temperature difference DT_s, which is used to replace T_s1 − T_s2 in Eqs 26 and 27 to give Eqs 30 and 31 as:

As instantly apparent, Eqs 30 and 31 can respectively and directly be solved for T_s1 and T_s2 independently in terms of the five dimensionless parameters ZT_i2, T_i, N_h, N_k and R_r. Now, making T_s1 and T_s2 subjects of their respective Eqs 30 and 31, we derived new analytical equations for T_s1 and T_s2 as given in Eqs 32 and 33 as:

Furthermore, by adding Eqs 26 and 27 or Eqs 30 and 31; we can easily get rid of the nonlinear quadratic terms as shown in Eqs 34 and 35 and by re-arranging and solving for T_s1 and T_s2, we can further derive a simpler T_s1 and T_s2 novel analytical formula as shown in Eqs 36 and 37.

Equation (34) was re-arranged to give Eq 35, such that the left hand side (LHS) of Eq (35) constitutes T_s1 without T_s2 and the right hand side (RHS) of Eq 35 constitutes T_s2 without T_s1. Finally, the LHS and RHS were both solved independently by equating each to zero, to easily obtain T_s1 in terms of ZT_i2, T_i, N_h, N_k, R_r and DT_s (NB: now six unknown parameters) and similarly to obtain T_s2 in terms of ZT_i2, N_k, R_r and DT_s (NB: now four unknown parameters) to respectively give Eqs 36 and 37 expressed as:

As seen, T_i and N_h are absent in Eqs 33 and 37, meaning these two dimensionless parameters, don't directly affects T_s2 or affects T_s2 indirectly through T_s1 and DT_s; thus, making these equations insightful.

Finally, we furthermore obtained the simplest optimal analytical relationship to calculate T_s1 and T_s2 by utilizing a set of optimal numerical values given in ^[7] as: ZT_i2 = 1, T_i = 2.6, N_h = 1, N_k = 0.3, R_r = 1.7 and DT_s = 0.8 and by substituting these values in Eq 35 and solving further, we obtain Eqs 38 and 39 as:

It should be noted that, the ZT_i2 = 1, T_i = 2.6, N_h = 1, at N_k = 0.3, R_r = 1.7 and now DT_s = 0.8 optimal values set are not arbitrarily chosen but specifically deduced optimal values, that are inter-linked to each other and obtained under the same specified operating conditions. Furthermore, parameters ZT_i2, T_i and N_h are initially provided or can be calculated from a chosen TEG-HS operating parameters and finally the optimal combinations of N_k, R_r and now DT_s, can all be established from the above specified dimensionless equations. NB: other optimal values set can be used in Eq 35 to get new Eqs 38 and 39 also.

Now, to validate our new T_s1 and T_s2 formulas, we substitute this the same chosen set of optimal values (given in ^[7]) of ZT_i2 = 1, T_i = 2.6, N_h = 1, N_k = 0.3, R_r = 1.7 and DT_s = 0.8 to finally determine, compare and verify the T_s1 and T_s2 values by using our newly derived simple, simpler and simplest analytical formulas; respectively using Eqs 32 and 33, Eqs 36 and 37 and Eqs 38 and 39 as follows:

Using the now simple formulas (Eqs 32 and 33), T_s1 and T_s2 are optimally directly computed as:

Using the now simpler formulas (Eqs 36 and 37), T_s1 and T_s2 are optimally directly deduced as:

As seen, both these new simple and simpler optimal solutions are approximately the same. We can then also easily verify both with Eqs 38 and 39 to get the simplest solution, by substituting T_s2 to calculate T_s1 or vice versa. We've now introduced three new analytical sets of formulas to easily compute T_s1 and T_s2 without going through the cumbersome numerical process. With T_s1 and T_s2 easily deduced, the other dimensionless parameters and ultimately all their physical dimensional parameters, can then be worked out using Eqs 1–25.

3.   Model simulator

The relevant TEG with heatsinks mathematics based on dimensional analysis was articulated in Section 2 and the formulas for T_s1 and T_s2 were developed and expressed analytically in respect of the now six independent parameters ZT_i2, T_i, N_k, N_h, R_r and DT_s. Furthermore, T_s1 and T_s2 were optimally simplified and validated in terms of T_s1 and T_s2 using Eqs 38 and 39. As apparent, the entire numerical as well as analytical procedures could be very tedious and as well subject to errors if manually done. In view of this, a Matlab and Simulink TEG model with heatsinks on both the TEG hot and cold sides, was then implemented using the formulas asserted in Section 2. Enclosed in the simulated TEG with heatsinks model, are the followings: i) The TEG parameters analysed in Section 2, which include; S the Seebeck coefficient, 273.15 the absolute temperature in kelvin, T_i1 fluid temperature of heatsink1, T_i2 fluid temperature of heatsink2, if required the TEG configuration in series T_s and or in parallel T_p, r the TEG thermocouple p-n junction resistance, k the TEG thermocouple p-n junction thermal conductivity, Ῥ the TEG thermocouple p-n junction electrical resistivity, a the TEG thermocouple p-n junction area, L the TEG thermocouple p-n junction length and P_req the TEG total output power required. ii) The heatsinks parameters, which include; Ab the heatsink base area, n₁ fins efficiency of heatsink1, h₁ convection coefficient of heatsink1 fluid, A₁ total area of heatsink1 with fins sides, n₂ fins efficiency of heatsink2, h₂ convection coefficient of heatsink2 fluid and A₂ total area of heatsink2 with fins sides. All these parameters can easily be input in the model simulator and the TEG optimal parameters value computed as depicted in Figure 2.

Additionally, the TEG dimensionless parameters such as; N_i the dimensionless output current, N_k the dimensionless thermal conductance, Eff* (n_th) the conversion efficiency, P_os (W_s1) the dimensionless output power, T_s1 the TEG hot-side dimensionless temperature, T_s2 the TEG cold-side dimensionless temperature and N_v the dimensionless output voltage, can all be either manually entered (as depicted in Figure 2 on the left side of the TEG with heatsinks model) or computed automatically. They're then used to finally calculate the TEG practical parameters values in SI units (as shown in Figure 2 on the right side of the TEG with heatsinks model). Some of the TEG physical parameters of interest that are computed include; Z the TEG figure of merit, n the TEG amount of p-n junction thermocouples, Q_h (Q₁) the TEG hot side heat absorbed, Q_c (Q₂) the TEG cold side heat released, P_o the TEG output power, T₁ (T_h) the TEG hot side temperature, T₂ (T_c) the TEG cold side temperature, R_R (R_t) the TEG internal resistance, V_o the TEG output voltage, V_oc is the TEG ideal output voltage, I the TEG output current, N the total amount of TEG modules required, DT or ΔT the TEG temperature difference and QPD is TEG with heatsinks heat flux density. The computed heatsinks parameters are: nhA the dimensionless convection conductance of heatsinks 1 or 2 fluid, T_i1 and T_i2 are fluids 1 and 2 hot and cold temperatures respectively, T_i the dimensionless fluid temperature, N_h the dimensionless convection conductance, ZT_i2 the dimensionless figure of merit at T_i2, ZTA the TEG dimensionless figure of merit at TA. The TEG with heatsinks normalised, maximum and effective parameters are also computed.

In terms of the aesthetics, the TEG-HS output and inputs parameters were sorted together based on their commonalities, as well as labeled and colour coded accordingly to make the TEG with heatsinks model simple to comprehend and as well user friendly as shown in Figure 2.

Finally, in addition to the TEG with heatsinks numeric model computations, miscellaneous characteristics curves of some of the TEG crucial parameters of interest such as Q_s1, Q_s2, P_os, Eff*, T_s1 and T_s2; were plotted against N_k, R_r and DT_s to graphically calculate Q_s1, Q_s2, P_os, Eff*, T_s1 and T_s2 different optimal values for ZT_i2 = 1, N_h = 1, T_i = 2.6, N_k = 0.1 − 0.4, R_r = 0.5 − 2 and DT_S = 0.1 − 1. These graphical results are in details, variously demonstrated next in Figures 3 to 6 and with the highlights summarised in Table 1.

4.   Simulation results, discussions and validation

In Section 3, the TEG with heatsinks simulated model was presented with highlights on its various parameters and the numeric results, as exemplified in Figure 2 and summarised in Table 2. These results as well as the numerous graphical results displayed extensively next, are based on the new simple and simpler analytical T_s1 and T_s2 formulas derived in Section 2. These are comparatively engaged and displayed next and various comparisons of Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* with different combinations of ZT_i2, N_h T_i, N_k, R_r and DT_S using the simple Eqs 32 and 33 as well as using the simpler Eqs 36 and 37 are shown and organised as follows. i) all the figures on the left are based on using the simple Eqs 32 and 33 T_s1 and T_s2 values, whereas all the figures on the right are based on using the simpler Eqs 36 and 37 T_s1 and T_s2 values. ii) Figure 3 graphs depict Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* results using ZT_i2 = 1, N_h = 1, T_i = 2.6, N_k = 0.1 − 0.4, R_r = 0 − 6 and DT_S = 0.8. iii) Figure 4 graphs depict Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* results using ZT_i2 = 1, N_h = 1, T_i = 2.6, N_k = 0.05 − 0.5, R_r = 0.5 − 2 and DT_S = 0.8. iv) Figure 5 graphs depict in 3D Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* results with ZT_i2 = 1, N_h = 1, T_i = 2.6, N_k = 0.3, R_r = 0.2 − 3.2 and DT_S = 0.1 − 1. v) Figure 6 graphs depict in 3D Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* results using ZT_i2 = 1, N_h = 1, T_i = 2.6, N_k = 0.05 − 0.5, R_r = 1.7 and DT_S = 0.1 − 1. Table 1 summarises the comparison of the simple Eqs 32 and 33 results and the simpler Eqs 36 and 37 results. Figures 2–6 results are also provided as supplementary picture files to enable every aspect of our findings to be viewed clearly in high quality without ambiguities.

Figure 3 graphs depict Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* results using ZT_i2 = 1, N_h = 1, T_i = 2.6, N_k = 0.1 − 0.4, R_r = 0 − 6 and DT_S = 0.8. It also compares the newly introduced simple vs simpler equations results.

Figure 4 graphs depict Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* results using ZT_i2 = 1, N_h = 1, T_i = 2.6, N_k = 0.05 − 0.5, R_r = 0.5 − 2 and DT_S = 0.8. It also compares the newly introduced simple vs simpler formulas results.

Figure 5 graphs depict in 3D Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* results with ZT_i2 = 1, N_h = 1, T_i = 2.6, N_k = 0.3, R_r = 0.2 − 3.2 and DT_S = 0.1 − 1. The newly introduced simple vs simpler equations results are also correlated.

Figure 6 graphs depict in 3D Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* results using ZT_i2 = 1, N_h = 1, T_i = 2.6, N_k = 0.05 − 0.5, R_r = 1.7 and DT_S = 0.1 − 1. The newly introduced simple vs simpler formulas results are compared.

In the first part of Section 4, the study results were variously and extensively displayed, revealing what combinations of R_r vs N_k vs DT_s values that will give optimal Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff^* outcomes; by employing our newly introduced analytical Eqs 32 and 33 which are termed the "simple equations" and comparing with the further simplified versions, Eqs 36 and 37 which as well are termed the "simpler equations". It should be noted that each set of figures results, validate each other set of figures results. That is, Figure 3 results validates Figure 4 results which in turn validates Figure 5 results and consequently which validates Figure 6 results. Furthermore, the figures on the left results (based on Eqs 32 and 33 exclusively) validate the figures on the right results (based on Eqs 36 and 37 exclusively). In the second part of Section 4, the highlights of our presented findings are discussed next and lastly the results are validated where applicable with reference to the original results presented in ^[7].

Beginning with Figure 3, parameters Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* are plotted against R_r with N_k = 0.1 − 0.4 and DT_s = 0.8 as well as ZT_i2 = 1, N_h = 1, T_i = 2.6—which were fixed throughout our entire study and therefore have no relevance in our results/discussions. NB: N_k = 0.1 − 0.4 instead of just N_k = 0.3 was used to see the effects of different N_k values on Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* at also different R_r values. As evident, Q_s1, Q_s2 and T_s2 are inversely proportional to R_r with N_k = 0.3 and DT_s = 0.8. T_s1 is directly proportional to R_r whereas P_os and Eff* are initially directly proportional to R_r, until P_os and Eff* reach their respective maximum value and then becomes inversely proportional to R_r and also with the rates of P_os and Eff* incline and decline, being faster using the "simpler" Eqs 36 and 37. With the "simple" equations, our optimal values for P_os is 0.045518 and for Eff* is 0.112407 and both occur at approximately R_r = 2 and not at exactly R_r = 1.7 as reported in ^[7]; however, at R_r = 1.7, our P_os and Eff* values of 0.0451302 and 0.10731 respectively, are exactly the same as those reported in ^[7]. With the "simpler" equations, our optimal values for P_os is 0.0893268 and for Eff* is 0.170213 and both occur at different R_r values—at R_r = 0.5 for P_os and at R_r = 1.1 for Eff*; however, at R_r = 1.7, both the P_os and Eff* values are respectively 0.0716775 and 0.16567, as well as respectively 0.0669887 and 0.161491 at R_r = 2. The rest parameters Q_s1, Q_s2, T_s1 and T_s2 values are closely the same using either the "simple" or "simpler" analytical equations. Further observable, the results dynamics have the following interesting aspects. i) Q_s1, Q_s2 and T_s2 values increase with increasing N_k values; however, they decrease with increasing R_r values as noticed and explained already. ii) T_s1 values increase with decreasing N_k values; as well as increase with decreasing R_r values as noticed and already explained. iii) Eff* values also increases with decreasing N_k values with the increase being more significant at lower values of R_r, especially at R_r values below 1.0. iv) P_os values are interesting—as the dynamics are irregular with different N_k values, with P_os values sharply increasing initially with decreasing N_k values and at very low R_r values and once after their respective optimal points, P_os values again increase but with now increasing N_k values and at high R_r values.

In Figure 4, parameters Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* are plotted against N_k with R_r = 0.5 − 2 and DT_s = 0.8 and using ZT_i2 = 1, N_h = 1, T_i = 2.6—which were fixed throughout the study. R_r = 0.5 − 2 instead of only R_r = 1.7 was used to see the effects of different R_r values on Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* at also over different N_k values. As can be seen, Q_s1, Q_s2 and T_s2 values are directly proportional to N_k. However, T_s1 and Eff* values are inversely proportional to N_k values. P_os initially increases proportionally to N_k till it reaches optimal point at N_k = 0.25 in Figure 4c (simple equation) and N_k = 0.3 in Figure 4i (simpler equation); after which P_os values becomes inversely proportional to N_k. Further noticeable, Q_s1, Q_s2 and T_s2 values are inversely proportional to R_r values, though directly proportional to N_k values as indicated earlier. Nevertheless, T_s1 is directly proportional to R_r values but indirectly proportional to N_k. Initially, with increasing N_k, P_os increases more with less R_r at mostly low N_k values till maximum P_os is attained and thereafter, P_os decreases more with less R_r at high N_k values. The optimal P_os of 0.0478788 using our simple equation occurs at R_r = 1.1 with N_k = 0.2, whereas at R_r = 1.7, the optimal P_os = 0.0474483 with N_k = 0.25—slightly contrary to ^[7] which reported an optimal P_os of 0.045 at R_r = 1.7 with N_k = 0.3—which we also have exactly the same result of P_os = 0.0451302 at R_r = 1.7 with N_k = 0.3; however, this result was not the exact optimal outcome in our case. Using our simpler equation, we have at R_r = 0.5, two P_os optimal values of 0.0966733 and 0.0977376 at respectively N_k = 0.2 and N_k = 0.25. Initially Eff* with increasing N_k, decreases more with more R_r at mostly low N_k values till optimal Eff* is attained and thereafter, Eff* decreases more with less R_r at high N_k values. This Eff* dynamics is mostly and clearly noticeable in Figure 4l using our simpler equation.

Figure 5 is a 3D plot of Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* parameters against R_r = 0.2 − 3.2 but with now DT_s = 0 − 1 and N_k fixed at 0.3. The results summarily reveal that P_os and Eff* increase proportionally with increasing DT_s, reach maximum and decrease, especially at lower R_r values but at higher R_r, P_os and Eff* increase linearly with DT_s. Also, Q_s1, Q_s2, T_s1 and T_s parameters exhibit comparable dynamics to those of Figures 3 and 4 results.

Figure 6 is another 3D plot of Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* parameters against N_k = 0.05 − 0.5 but now with DT_s = 0 − 1 and R_r fixed at 1.7. As apparent, Q_s1, Q_s2 and T_s2 have similar dynamics as already noted and discussed earlier in the previous figures. The only difference now is, these dynamics happen across different values of DT_s. The same dynamics applies to T_s1 and Eff*. P_os generally increases proportionally to DT_s at especially lower N_k values but at higher N_k values, P_os decreases with DT_s as revealed.

It's worth stating that, these different plots simply show the different dynamics involved and generally, it's up to designers to now choose optimal values depending on the design constraints and design objectives.

To validate our findings in Figures 2–6 and Tables 1 and 2 where applicable, Figure 7 and Table 3 results adapted from ^[7] were used. In summary, validating our findings using ^[7] results affirms i) exactly the accuracy of our simple equations results and ii) fairly the approximation of our simpler equations results. It should be noted that, just reference ^[7] was used to validate our findings because a) reference ^[7] is the original research now developed and simplified with direct analytical formulas to easily calculate T_s1 and T_s2 and therefore it's only prudent to validate the new developments with the original study and b) reference ^[7] is the only study we know with the aforementioned method and having the stated datasets that we can use.

5.   Conclusions

Sustainable energy is transcending to become the future of energy to complement the national grid and for private use. As a result, we've researched thermoelectricity (TEG) as an alternative energy source for household applications requiring low DC power and lighting. However, practical TEG use requires heatsinks to work efficiently and reliably but regrettably adding heatsinks add thermal resistances, which consequently degrades the TEG efficiency as well as its output power. We investigated various techniques, especially dimensional analysis and shortlisted the approach by ^[7]—which converts thermal resistance to convection conductance, making it easier for practical use. However, the formulations in ^[7] could not find the exact analytical formulas to directly calculate T_s1 and T_s2; as a result, ^[7] used numerical analysis (a bit cumbersome) which we developed further by introducing DT_s to simplify and derive new simple accurate analytical formulas that can be applied to compute T_s1 and T_s2 directly. Furthermore, these simple formulas were further simplified to obtain simpler analytical equations for T_s1 and T_s2. Finally, our simplest T_s1 and T_s2 formulas and their optimal relationships were established. These formulas were all verified with chosen optimal values and all gave approximate results that correlated each other. We further used these newly introduced formulas to model, numerically simulate and plotted various characteristics curves of Q_s1, Q_s2, P_os, T_s1, T_s2 and Eff* against R_r vs N_k vs DT_s using Matlab and Simulink. These results were articulated variously in both 2D and 3D plots and finally the equivalent outcomes were comparatively discussed and validated using each other and as well with results asserted in ^[7]. In ^[7], P_os of 0.045 is optimal with N_k = 0.3 at R_r = 1.7, contrary to our study which gave optimal P_os of 0.0474483 at N_k = 0.25 using our simple equations; however, using our simpler equations, it gave optimal P_os of 0.0716775 at N_k = 0.3. Notwithstanding, our simple equations gave P_os of 0.0451302 at N_k = 0.3 (though not optimal in our case) which corresponds exactly to optimal P_os value of 0.045 with N_k = 0.3 as in ^[7]. In addition, the magnitude of P_os and Eff* are more (almost doubled) using our simpler Eqs 36 and 37, compared to using our simple Eqs 32 and 33. Furthermore, the optimal values of P_os and Eff* using both our "simple" and "simpler" equations, respectively occurred at different values of R_r. In sum, Eff* is best at very low N_k values. The highlights of our study are marked in the various figures and also summarised in Tables 1 and 2. It's worth mentioning that, DT_s = 0.8 for this set of optimal values, is rightly only equal to T_s1 − T_s2 = ~0.8 only at R_r = 1.7 and N_k = 0.3 as summarised in Table 1 and not at any other values when used in both our simple and simpler equations. Finally, in as much as our newly introduced simpler analytical formulas didn't give exact P_os and Eff* results correlating those in i) ^[7] and ii) our newly introduced analytical simple equations (that gave accurate results), we reckon the simpler equations can be used as first approximation to quickly find T_s1 and T_s2 as well as Q_s1 and Q_s2 values with fair accuracy. In sum, T_s1 and T_s2 are the most vital physical parameters, as they're those that can be easily practically manipulated to achieve the other values. However, first ascertaining the electrical load or R_r value to ensure maximum power transfer is most paramount to achieve optimal result, as a TEG as well TEC, are non-linear devices whose dynamics must be well understood before embarking on a practical design. The next phase of our research is to validate and refine our study practically and then design a household TEG alternative DC power source.

Acknowledgments

Thanks to the Cape Peninsula University of Technology (CPUT) and HySA Systems at the University of the Western Cape (UWC), Cape Town, South Africa for funding the study.

Conflicts of interest

The authors declare no conflict of interest.

Supplementary

Figures 1–6 are provided in link below in picture format for better quality viewing/review.

https://drive.google.com/drive/folders/1VHN2mHJ_wtYaztunL_Sn21XJmppc0o8z?usp=sharing