# Thermal Modeling Using Two-Port Network Impedance Fractional-Order Approximations

Jean-François Duhé, Stéphane Victor, Pierre Melchior, Youssef Abdelmounen, François Roubertie

► To cite this version:

Jean-François Duhé, Stéphane Victor, Pierre Melchior, Youssef Abdelmounen, François Roubertie. Thermal Modeling Using Two-Port Network Impedance Fractional-Order Approximations. 17th ASME/IEEE MESA 2021 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference, Aug 2021, Virtual, United States. 10.1115/DETC2021-69968. hal-03548461

HAL Id: hal-03548461

https://hal.science/hal-03548461v1

Submitted on 10 Jan 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
THERMAL MODELING USING TWO-PORT NETWORK IMPEDANCE
FRACTIONAL-ORDER APPROXIMATIONS

Jean-François Duhé, Stéphane Victor,∗ Pierre Melchior

Univ. Bordeaux, CNRS, IMS-UMR 5218,
Bordeaux INP/enseirb- matmeca,
351 cours de la Libération,
33405 Talence Cedex France

Email: jean-francois.duhe@ims-bordeaux.fr,
stephane.victor@ims-bordeaux.fr
pierre.melchior@ims-bordeaux.fr

Youssef Abdelmounen, François Roubertie

IHU Liryc, Electrophysiology and Heart Modeling Institute,
Fondation Bordeaux Université 33000 Bordeaux Cedex France
INSERM, Centre de recherche
Cardio-Thoracique de Bordeaux,
33000 Bordeaux Cedex France

Email: dr.a.youssef@free.fr, rouberf@me.com

ABSTRACT

Sufficiently accurate thermal modeling is necessary for many applications such as heat dissipation, melting processes, building design or even bio-heat transfers in surgery. Circuit models help modeling heat transfer dynamics: this method is simple and is often used to model thermal phenomena. However, such models well approximates low and high frequency behavior but they are not accurate enough in the middle band of interest, thus lacking of precision in dynamical terms. A more complete and accurate description of conductive heat transfer can be obtained by using a two-port network. The resulting analytical expressions are complex and nonlinear in the frequency ω. This complexity in the frequency domain is difficult to handle when it comes to control applications and more specifically in real-time applications such as surgery. Consequently, an analysis of this thermal two-port network in the frequency domain directly leads to fractional-order systems. A frequency domain analysis of the series and shunt impedances will be presented and different approximations will be explored in order to obtain simple but sufficiently precise linear fractional transfer function models. The series impedances are approximated by using asymptotic and pole-zero approximations and the shunt impedance is approximated by using a capacitance approximation and two fractional model approximations.

NOMENCLATURE

ρ Density

k Thermal conductivity

c Specific heat capacity

T Temperature

Q̇ Heat flow

S<sub>w</sub> Surface

s Laplace variable

<footnote>∗ Address all correspondence to this author.</footnote>
# INTRODUCTION

Thermal modeling of systems may help to control situations in which temperature is critical, as it can be the case of microchip cooling or a bio-medical application. It is common to use finite element methods in order to solve the heat equation, but this models may be far more complex than required in order to meet a desired criteria in the application. It may be possible to obtain a simpler model that allows to attain specific objectives.

Simple RC circuits have been widely used to model the dynamics of heat transfer in different applications, such as microchip cooling [1], building simulation [2] and even to model human heat loss [3]. This type of models are also present in similar applications like lithium-ion battery models [4]. Most of this applications will suppose that thermal fluctuations are in the low-frequency range. The thermal two-port network [5] still offers a circuit model to analyze thermal dynamics, but doesn't suppose low frequency behavior.

If one considers a 1D plane wall with heat conduction, the thermal two-port network represents the wall as a T circuit with identical series impedance $Z_1(s) = Z_2(s)$ and a shunt impedance $Z_3(s)$. The Laplace-domain expressions of these impedance are complex and the transfer function models obtained are not a ratio of s polynomials.

For a case with the end of the wall isolated ($\phi_{out} = 0$), it may be proved that the input thermal impedance will behave as a half-order integrator for high frequency ([6]):

$$Z_{thermal}(s) = \frac{T(s)}{\phi(s)} = \frac{K}{\sqrt{s}} \qquad (1)$$

High-frequency heat transfer may then be approximated by the means of constant-phase elements (see [7] and [8]). This type of element has been used to model irregular systems, such as porous films [9] or intestine tissue [10] and proves the growing interest in fractional calculus for modeling.

Even though fractional-order thermal models are present in literature, in some cases the basic notion of the half-order derivative of equation (1) is used as a starting point to establish a general mathematical structure, such as in [11] in which a commensurate-order fractional transfer function is proposed (with order 0.5). Other models are only loosely based on the fact that thermal phenomena implies fractional-order derivatives and general fractional transfer function black-box models are used [12]. Fractional-order transfer functions that even lack the classic polynomial structure may also be found, as is the case with the lung model by means of Havriliak-Negami function [13]. The thermal two-port network has also been presented and fractional models have been established to approximate the results obtained with the exact impedance expressions [14]. However, the fractional transfer function is still a black-box model and the flexibility of a two-port representation is not fully exploited as an equivalent impedance model for a specific case is used.

The approximations presented try to keep the flexibility of the thermal two-port network by separately analyzing each impedance of the general thermal circuit (which will be explained in section 2) and combining the use of circuit elements and mathematical corrective terms. Two approximations will be presented for the series impedance $Z_1(s)$ and one for the shunt impedance $Z_3(s)$, which keep (at least partially) circuit elements in their expressions.

This paper is organized as follows. The presentation of the thermal two-port network and its frequency-domain analysis will be presented in section 2. Section 3 includes the approximation propositions and simulations for $Z_1(s)$ and $Z_3(s)$. Finally, the conclusions will be given in section 4.

# 1 The analytical thermal two-port network

## 1.1 Heat equation on a unidirectional medium

The case considered will be of a simple plane wall with heat conduction in its longitudinal direction $x$. This plane wall will have a thermal conductivity $k$, density $\rho$ and heat capacity $c$. The length is $L$ and the cross section is $S_w$, as shown in figure 1.

Starting by the heat equation:

$$\rho c \frac{\partial T}{\partial t} = k \nabla^2 T(x, y, z, t), \qquad (2)$$

by considering conduction only goes in the $x$ direction and taking the Laplace transform of the heat equation:

$$sT(x, s) = \frac{k}{\rho c} \frac{\partial^2 T(x, s)}{\partial x^2}. \qquad (3)$$
<figure><img src="image_2.png" /><figcaption>FIGURE 1. 1D thermal system</figcaption></figure>

By considering a heat input at $x = 0$ and an output at $x = L$:

$$ 
\begin{aligned}
\dot{Q}_{in}(s) &= -kS_w \left. \frac{\partial T(x,s)}{\partial x} \right|_{x=0} \\
\dot{Q}_{out}(s) &= -kS_w \left. \frac{\partial T(x,s)}{\partial x} \right|_{x=L}
\end{aligned}
\qquad (4) $$

put under matrix form, one gets:

$$ \begin{bmatrix} T_{in}(s) \\ \dot{Q}_{in}(s) \end{bmatrix} = M \begin{bmatrix} T_{out}(s) \\ \dot{Q}_{out}(s) \end{bmatrix} \qquad (5) $$

where:

$$ M = \begin{bmatrix} \cosh(\delta L) & \frac{1}{kS_w\delta}\sinh(\delta L) \\ kS_w\delta\sinh(\delta L) & \cosh(\delta L) \end{bmatrix} \qquad (6) $$

with $\delta = \sqrt{\frac{s}{a}}$ and $a = \frac{k}{\rho c}$.

It should be noted that $\delta$ involves the presence of fractional calculus, as its time-domain representation is a half-order derivative. This matrix equation can be represented by a T circuit model, such as represented in figure 2 (see Appendix A).

<figure><img src="image_3.png" /><figcaption>FIGURE 2. Thermal two-port network</figcaption></figure>

The expressions for the impedances are given by:

$$ Z_1(s) = Z_2(s) = \frac{1}{kS_w\delta} [\coth(\delta L) - \cosh(\delta L)] \qquad (7) $$

$$ Z_3(s) = \frac{1}{kS_w\delta}\cosh(\delta L) \qquad (8) $$
where coth and csch respectively are the hyperbolic cotangent and cosecant functions.

## 1.2 High and low-frequency limits

If the frequencies are supposed to be low, the impedance of the T circuit may then be reduced to a resistance and a capacitance, namely:

$$ \lim_{\omega \to 0} Z_1(j\omega) = \frac{L}{2kS_w} = R \qquad (9) $$

$$ \lim_{\omega \to 0} Z_3(j\omega) = \frac{1}{\rho c L S_w j\omega} = \frac{1}{j\omega C_t} \qquad (10) $$

where $C_t = \rho c L S_w$.

If one wants to analyze the high-frequency limit, one gets the following for the gains and arguments of the circuit elements:

$$ \lim_{\omega \to \infty} |Z_1(j\omega)| = 0, \lim_{\omega \to \infty} \arg|Z_1(j\omega)| = -45^\circ \qquad (11) $$

$$ \lim_{\omega \to \infty} |Z_3(j\omega)| = 0, \lim_{\omega \to \infty} \arg|Z_3(j\omega)| = -\infty. \qquad (12) $$

Both elements tend to block high-frequencies (both high-frequency gains are zero), but the series impedance has a constant $-45^\circ$ phase at high frequency. This allows us to deduce a high frequency model for this impedance:

$$ Z_{1-HF}(s) = \frac{1}{C_s s^{0.5}} \qquad (13) $$

The shunt capacitance seems coherent with its capacitance approximation in terms of its gain, but its high-frequency argument suggests that there is an additional low-pass filtering added as frequency increases:

$$ Z_{3-HF}(s) = \frac{1}{C_t s} H_{\text{filter}}(s). \qquad (14) $$

*Remark:* It is not possible to obtain a finite-order filter which fully corrects $Z_3(s)$ as its high-frequency argument is negative infinity.

## 2 Impedance approximations

Even though it is impossible to obtain a model valid for any frequency by using transfer function models for the quadrupole, a truncated frequency-model may be a fairly useful approximation with a wider validity compared to RC circuits. A model going from low-frequency (or even a static case) to a high but finite upper frequency will be developed. An academic example will be used to illustrate the main ideas of sections 3 and 4, and the parameters are presented in Table 1.

The proposed approximations will require the determination of different parameters in order to be as close as possible to the real impedance behavior. A common quadratic error criteria will be used by taking the error of the gain in dB:

$$ J(\theta) = \frac{\sum_{i=1}^{N} \lambda^{2(N-i)} [|Z(j\omega_i)|_{dB} - |Z_{app}(j\omega_i, \theta)|_{dB}]^2 [\log \omega_{i+1} - \log \omega_i]}{\log \omega_N - \log \omega_1} \qquad (15) $$

where $\omega$ is a frequency vector that will go from 0.1 rad/s to 100 rad/s. This criteria measures the error in the gain curve given in dB. $\lambda$ will be a factor which may add more weight to higher frequencies.

The approximations for the series impedance will be optimized by using the Flower Pollination Algorithm [15], which is a novel meta-heuristic optimization method. The parameters to be tuned in order to apply FPA are pretty low, which makes it faster to tune and test when compared to other algorithms. The main parameters to be tuned are the population size $n$ (normally between 10-25) and a switching probability $p$ which states the possibility of doing one of the two principal operations of the algorithm in an iteration: global
<table><caption>TABLE 1. Academic example simulation parameters</caption><thead><tr><td>Parameter</td><td>Value</td></tr></thead><tbody><tr><td>a</td><td>1 m<sup>2</sup> &middot; s<sup>-1</sup></td></tr><tr><td>k</td><td>1 W &middot; m<sup>-1</sup> &middot; K<sup>-1</sup></td></tr><tr><td>L</td><td>1 m</td></tr><tr><td>S<sub>w</sub></td><td>1 m<sup>2</sup></td></tr></tbody></table>

<table><tr><td>Algorithm 1: Flower Pollination Algorithm</td></tr><tr><td>Objective function f(x), x = (x<sub>1</sub>, x<sub>2</sub>, ..., x<sub>d</sub>);</td></tr><tr><td>Initialize a population of n flowers with random solutions;</td></tr><tr><td>Find the best solution g<sub>best</sub> in the initial population;</td></tr><tr><td>Define a switch probability p;</td></tr><tr><td>while t &lt; MaxIterations do</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;for i = 1 : n do</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;if rand &lt; p then</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Calculate L (Levy flight distribution);</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;x<sub>i</sub><sup>t+1</sup> = x<sub>i</sub><sup>t</sup> + L(x<sub>i</sub><sup>t</sup> - g<sub>best</sub>);</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;else</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;x<sub>i</sub><sup>t+1</sup> = x<sub>i</sub><sup>t</sup> + &epsilon;(x<sub>j</sub><sup>t</sup> - x<sub>k</sub><sup>t</sup>);</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;end</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;end</td></tr></table>

pollination or local pollination. The most critical parameter is the switching probability, but its creator suggests a value of $p = 0.8$. The main steps for classic FPA are shown in Algorithm 1.

*Remark:* The choice of FPA as an optimization tool was mainly due to its simple structure and low parameters to be tuned. However, it should be noted that other optimization algorithms could also be used to approximate the parameters, such as genetic algorithm (GA) or particle swarm optimization (PSO).

## 2.1 Approximation of impedances $Z_1$ and $Z_2$

As it was seen in the previous section, the series impedance behaves like a resistance for low frequencies and like a constant phase element at high frequency. This suggests a simple model of this element as the parallel combination of a resistance and a fractance:

$$Z_{asymp}(s) = \frac{R}{1 + RC\sqrt{s}} \qquad (16)$$

This first proposition will be named “asymptotic approximation”, as it is obtained by only considering both extremes of the frequency domain. A simple expression can be obtained for the bandwidth:

$$\omega_b \approx \frac{0.2679}{(RC)^2} \qquad (17)$$

Gain comparison between the exact expression and the asymptotic approximation is shown in figures 3. This plots show a larger frequency band to clearly show the asymptotic behaviour of this impedance.
<figure><img src="image_4.png" /><figcaption>FIGURE 3. Gain frequency responses for Z₁ and Z₁-asymp</figcaption></figure>

Both curves tend to get closer on both low and high-frequency extremes, but the mid-band error is non-negligible and thermal phenomena occurring in this range may be inaccurately modelled with this approximation. It may also be noted that asymptotic approximation bandwidth is almost a decade lower than exact impedance's one, which is another drawback to this approach.

In order to better fit the mid-band frequencies, an additional correction term is proposed as the asymptotic approximation can already correctly handle both low and high frequency limits. As the correction term should not make any phase or gain contribution outside the mid-band range, a pole-zero cell structure is proposed:

$$Z_{1-approx}(s) = \frac{R}{1 + RC_s \sqrt{s}} \prod_{i}^{N_{cells}} \frac{1 + \frac{s}{z_i}}{1 + \frac{s}{p_i}} \quad (18)$$

The parameter vector is given by:

$$\theta = [p \quad z] \quad (19)$$

where p and z are pole and zero vectors, respectively.

The total number of correction cells $N_{cells}$ required to get an optimal approximation may vary depending on the parameters. It should be noted that increasing $N_{cells}$ may reduce the error, but the number of added parameters should not be neglected (each cell adds 2 parameters to be determined). Simulations were run for different values of $N_{cells}$ and error criterion $J(\theta)$ was calculated for each model. $\lambda = 1$ to give equal weights for all frequencies as all asymptotic slope are constant in low and high frequency. The plot is shown in figure 4.

The Pham Information Criterion as proposed in [16] was calculated in order to choose an adequate number of cells. The criterion is defined as:

$$PIC = SSE + k \left[ \frac{n-1}{n-k} \right] \quad (20)$$

where SSE is the sum of squared errors, k the number of parameters and n is the number of observations. The criterion is plotted in figure 5.

The optimal choice according to PIC criterion is found for $N_{cells} = 3$, which is the chosen cell number for further simulations. It should be noted that the expression for this impedance may be fairly complex, as the parameter dimension for the correction term
<figure><img src="image_5.png" /><figcaption>FIGURE 4. J criterion versus the number of cell approximations $N_{cells}$</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>FIGURE 5. PIC versus the number of cell approximations $N_{cells}$</figcaption></figure>

is $dim(\theta) = 6$ for this case. Figure 6 shows the gain curves for the exact expression and the pole-zero cell approximation. Figure 7 compares the errors obtained with the asymptotic and the one obtained with the zero-pole correction, where the error is computed as:

$$ \varepsilon(j\omega) = \left| Z_{1-exact}(j\omega) \right|_{dB} - \left| Z_{1-approx}(j\omega) \right|_{dB}. \qquad (21) $$

As it can be seen, the proposed approximation, namely the pole-zero cell approximation, fits better the exact impedance than the asymptotic one: the pole-zero approximation has a maximum error of 0.4dB, whereas the asymptotic approximation presents a maximum error of 5dB. Moreover, over the whole chosen frequency band, the asymptotic approximation error is above 1dB, whereas it never even reaches 0.5dB for the pole-zero cell approximation.
<figure><img src="image_7.png" /><figcaption>FIGURE 6. Gain diagrams of Z₁ and approximation Z₁-approx with N<sub>cells</sub> = 3</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>FIGURE 7. Gain errors for asymptotic Z₁-asymptotic and pole-zero approximations Z₁-approx</figcaption></figure>

## 2.2 Shunt impedance Z₃

The shunt impedance behaves as a pure capacitance for low frequencies, but its high frequency behavior is more complex than that of the series impedance, namely with $C_t = kLS_w/a = 1$:

$$Z_{3-cap}(s) = \frac{1}{C_t s}. \qquad (22)$$

The gain frequency response for the shunt impedance is shown in figure 8 and the corresponding gain error to the exact Z₃ impedance in figure 9 in red (—).

The low frequency gain slope is -20dB/dec, which corresponds to a capacitance model. It can be seen that this slope is not constant and increases with frequency, which limits the validity of the capacitance model. Therefore, it is suggested to modify the capacitance model to include the slope increment.

As the final slope in for the frequency range may not be a multiple of -20dB/dec, a correction proposition would be to add a
<figure><img src="image_9.png" /><figcaption>FIGURE 8. Gain diagrams for Z3 and the proposed approximations</figcaption></figure>

fractional low-pass filter to the capacitance model:

$$Z_{3-frac1}(s) = \frac{1}{C_t s} \left[ \frac{1}{1 + (\tau s)^\phi} \right]. \qquad (23)$$

$\phi$ may have any positive value and is not constrained to be an integer and $\tau$ indicates the breaking point for the addition of this extra term. The parameter vector for this case is:

$$\theta_{frac-slope} = [\tau \quad \phi]. \qquad (24)$$

Fractional stability conditions should be checked *a priori* for this type of model. The proposed filter corresponds to the following general structure:

$$G(s) = \frac{1}{s^{\alpha} + b}. \qquad (25)$$

This system may be seen as a commensurate order system (with commensurate order $\alpha$) and therefore, Matignon's stability theorem may be applied [17]:

$$|\arg(\lambda_i)| > \alpha \frac{\pi}{2} \qquad (26)$$

where $\lambda_i = s^\alpha$ and $\lambda_i$ are the roots of the characteristic polynomial in $\lambda$. Figures 8 and 9 respectively show the gain curves and errors (in purple ... ) obtained by approximating the shunt impedance to a pure capacitance and by adding the fractional slope. For the general structure presented in equation (25):

$$\lambda_i = -b \Rightarrow |\arg(-b)| = \pi. \qquad (27)$$

Therefore, the stability condition reads $\alpha < 2$ for $b > 0$.
<figure><img src="image_10.png" /><figcaption>FIGURE 9. Gain errors for capacitance $Z_{3-cap}$, one cell fractional $Z_{3-frac1}$ and two cell fractional $Z_{3-frac2}$ approximations</figcaption></figure>

It may be interesting to estimate the low and high frequency slopes on figure 8 to have an initial guess on this value:

$$ \left. \frac{d|Z_3(j\omega)|}{d\omega} \right|_{\omega=0.1} = -20.00 \, \text{dB/dec} \qquad (28) $$

and

$$ \left. \frac{d|Z_3(j\omega)|}{d\omega} \right|_{\omega=100} = -78.28 \, \text{dB/dec.} \qquad (29) $$

The slope difference is of $-68.28 \, \text{dB/dec}$, which may allow us to calculate an initial guess for fractional order $\phi \approx -68.28/20 = 2.90$ which will give a result close to the exact final slope. Unfortunately, such a fractional order value violates the stability condition of commensurate systems. By limiting the $\phi$ to an upper bound of 2, the FPA algorithm has estimated $\phi = 1.55$. Note that the weighting coefficient $\lambda = 0.995$ to give more weight for high frequencies as the asymptotic slope increases gradually in high frequency. Normally, one could have expected the optimization to saturate around the upper boundary of this parameter ($\phi_{max} = 2$) which may have been able to further reduce the errors around higher frequencies as it can add a bigger slope. This would have been the case for a weighted criteria with a higher weight on the high-frequency spectrum. However, we may see this result as a trade-off: higher errors around the upper frequency and overall curve fitting. Error curve shows that this model still manages to significantly reduce errors produced by the capacitance approximation. The gain curves and errors are shown in figures 8 and 9 in green (-).

In order to further improve high-frequency approximation, a double fractional-slope model is proposed:

$$ Z_{3-frac2}(s) = \frac{1}{C_t s} \left[ \frac{1}{1 + (\tau_1 s)^\nu} \right] \left[ \frac{1}{1 + (\tau_2 s)^\nu} \right]. \qquad (30) $$

The high frequency slope added by the double fractional slope model is $-40 \, \text{dB/dec}$, which means expected $\nu$ should be around 1.45, thus satisfying the stability condition (26). The gain and error curves are shown in figures 8 and 9.

By limiting $\nu$ to an upper bound of 2 for stability reasons, the FPA algorithm has estimated $\nu = 1.15$ which is rather surprising and may suggest that this particular scenario may also be well simulated by using a series of integer poles without increasing the number of parameters (as only 3 poles would be required). Again, the same trade-off is obtained as higher errors around the upper frequency are obtained for overall curve fitting.

In this example (see table 1), the slope difference suggests that a 3-pole correction could have been used. However, fractional order models allow to reduce parameter vectors and still correctly handle high frequency. For example, a scenario which would require 6 poles could be handled by 3 fractional slopes.
<table><caption>TABLE 2. Models used for the time simulation</caption><thead><tr><th>Model</th><th>Z<sub>1</sub></th><th>Z<sub>3</sub></th></tr></thead><tbody><tr><td>A</td><td>Z<sub>1-asymp</sub></td><td>Z<sub>3-cap</sub></td></tr><tr><td>B</td><td>Z<sub>1-approx</sub></td><td>Z<sub>3-frac1</sub></td></tr><tr><td>C</td><td>Z<sub>1-approx</sub></td><td>Z<sub>3-frac2</sub></td></tr></tbody></table>

<figure><img src="image_11.png" /><figcaption>FIGURE 10. Time responses of the thermal model approximations</figcaption></figure>

### 3 Time simulations of the proposed approximations

By using the previous approximations proposed for $Z_1$ and $Z_3$ (see Section 2), a time simulation with a simple scenario is proposed. For a thermal wall insulated at $x = L$ ($\phi_{out} = 0$), a sine temperature fluctuation input is applied:

$$T_{in}(t) = T_0 \sin(\omega t) \qquad (31)$$

with $T_0 = 1K$ and $\omega = 6rad/s$. The transfer function which relates input and output temperatures is given by:

$$\frac{T_{out}(s)}{T_{in}(s)} = \frac{Z_{3}(s)}{Z_{1}(s) + Z_{3}(s)}. \qquad (32)$$

The combination for $Z_1$ and $Z_3$ approximations are shown in table 2 and the simulation results are shown in figure 10.

Model A is made of the two asymptotic impedances: $Z_{1-asymp}$ for $Z_1$ and $Z_{3-cap}$ for $Z_3$. It presents the highest magnitude as it does not considers the additional gain loss for $Z_3$. Model A shows the worst time response as the corresponding approximations are the worst: indeed, both approximation errors in gain and phase are directly present in the time simulation.

Model B is made of the two approximations impedances: $Z_{1-approx}$ for $Z_1$ and $Z_{3-frac1}$ for $Z_3$. A slight improvement is present for its magnitude, but the resonance peak present in $Z_{3-frac1}$ is close to the sine input frequency, which mitigates its gain adjustment. Model B is a compromise in the time response results as the corresponding approximations are better fit: indeed, the errors are more present in gain and less in the phase which are directly present in the time simulation.

Model C is made of the two approximations impedances: $Z_{1-approx}$ for $Z_1$ and $Z_{3-frac2}$ for $Z_3$. The magnitude is the lowest for model C, as there is no resonance peak for the $Z_{3-frac2}$ model and $Z_{1-approx}$ is not overestimated as in model A. At the beginning, model
A answers faster as compared to the other models, which have an initially slower transient (from 0 to 0.3s approximately). This suggests that models B and C more successfully try to reproduce typical delay phenomena related to heat accumulation than model A. Model C behaves as the true system as its corresponding approximations have the best fit: indeed, the errors in gain and phase are the smallest, which gives promising results.

## 4 Conclusion and final remarks

The frequency behavior of thermal impedance was analyzed and the possibility of higher frequency fluctuations of heat and temperature were considered. It was then shown that the classic RC circuit model is limited to low-frequencies and different approximations were proposed to take into account higher frequencies for the circuit models. The asymptotic approximation succeeds in approximating both extremes in frequency domain for the series impedance, but its error in the mid-band frequency is not be negligible. A series of poles and zeros were added to reduce the mid-band frequency error and a more accurate approximation of the impedance was obtained.

The shunt impedance has a more difficult frequency behavior as it was only possible to approximate it by a known element for the lower frequencies. The pure capacitance model is not accurate enough in high frequency if there are heat or temperature variations. The fractional slope helped in the reduction of the error around the higher frequency range and kept the number of additional parameters low. A better approximation was obtained by using a double fractional slope model, which only introduced a single additional parameter.

Research perspectives may focus on analyzing this two port-networks for more complex geometries, such as radial transfer in cylinders and spheres as well as trying to extend this models for multi-dimensional scenarios.

## REFERENCES

[1] Künzi, R., 2015. "Thermal design of power electronic circuits". arXiv: Accelerator Physics, 3, pp. 311-327.

[2] Parnis, G., 2012. Building thermal modelling using electric circuit simulation.

[3] Jiang, G., Qu, T., Shang, Z., and Zhang, X., 2004. "A circuit simulating method for heat transfer mechanism in human body". In The 26th Annual International Conference of the IEEE Engineering in Medicine and Biology Society, Vol. 2, pp. 5274-5276.

[4] Zhang, L., Peng, H., Ning, Z., Mu, Z., and Sun, C., 2017. "Comparative research on RC equivalent circuit models for lithium-ion batteries of electric vehicles". Applied Sciences, 7(10).

[5] Maillet, D., André, S., Batsale, J., Degiovanni, A., and Moyne, C., 2000. *Thermal Quadrupoles: Solving the Heat Equation through Integral Transforms*. Loyola Symposium Series. John Wiley & Sons.

[6] Malti, R., Sabatier, J., and Akçay, H., 2009. "Thermal modeling and identification of an aluminium rod using fractional calculus". In 15th IFAC Symposium on System Identification (SYSID'2009), pp. 958-963.

[7] Krishna, B., 2011. "Studies on fractional order differentiators and integrators: A survey". Signal Processing, 91(3), pp. 386-426.

[8] Nakagawa, M., and Sorimachi, K., 1992. "Basic characteristics of a fractance device". *IEICE Transactions on Fundamentals of Electronics, Communications and Computer Sciences*, **75**, pp. 1814-1819.

[9] Das, S., Sivaramakrishna, M., Das, S., and Biswas, K. and Goswami, B., 2009. "Characterization of a fractional order element realized by dipping a capacitive type probe in polarizable medium". In Symposium on Fractional Signals and Systems.

[10] Elwakil, A., 2010. "Fractional-order circuits and systems: An emerging interdisciplinary research area". *IEEE Circuits and Systems Magazine*, **10**(4), pp. 40-50.

[11] Maachou, A., Malti, R., Melchior, P., Battaglia, J., and Hay, B., 2012. "Thermal system identification using fractional models for high temperature levels around different operating points". *Nonlinear Dynamics*, **70**(2), pp. 941-950.

[12] Długosz, M., and Skruch, P., 2015. "The application of fractional-order models for thermal process modelling inside buildings". *Journal of Building Physics*, **39**(5), 2021/04/28, pp. 440-451.

[13] Victor, S., Melchior, P., Pellet, M., and Oustaloup, A., 2020. "Lung thermal transfer system identification with fractional models". *IEEE Transactions on Control Systems Technology*, **28**(1), pp. 172-182.

[14] Gabano, J.-D., and Poinot, T., 2011. "Fractional modelling and identification of thermal systems". *Signal Processing*, **91**(3), pp. 531-541. Advances in Fractional Signals and Systems.

[15] Yang, X., 2012. "Flower pollination algorithm for global optimization". In Unconventional Computation and Natural Computation, J. Durand-Lose and N. Jonoska, eds., Springer Berlin Heidelberg, pp. 240-249.

[16] Pham, H., 2019. "A new criterion for model selection". Mathematics, **7**(12).

[17] Matignon, D., 1998. "Stability properties for generalized fractional differential systems". *ESAIM proceedings - Systèmes Différentiels Fractionnaires - Modèles, Méthodes et Applications*, **5**.
## A Appendix A: Impedance expressions from the change of a matrix formulation to a T circuit

For a T network as shown in figure 2, the input-output relationship may be expressed by means of a transmission matrix:

$$ \begin{bmatrix} T_{in}(s) \\ \dot{Q}_{in}(s) \end{bmatrix} = \begin{bmatrix} A & B \\ C & D \end{bmatrix} \begin{bmatrix} T_{out}(s) \\ \dot{Q}_{out}(s) \end{bmatrix}. \qquad (33) $$

The transmission matrix can be expressed in function of impedance $Z_1$, $Z_2$ and $Z_3$:

$$ \begin{bmatrix} A & B \\ C & D \end{bmatrix} = \begin{bmatrix} 1 + \frac{Z_1}{Z_3} & Z_1 + Z_2 + \frac{Z_1 Z_2}{Z_3} \\ \frac{1}{Z_3} & 1 + \frac{Z_2}{Z_3} \end{bmatrix} \qquad (34) $$

Therefore, the following impedances can be expressed as:

$$ Z_1 = \frac{A-1}{C}, \quad Z_2 = \frac{D-1}{C}, \quad Z_3 = \frac{1}{C} \qquad (35) $$

which leads to expressions (7) and (8). Note that $Z_1 = Z_2$ in this geometry as $A = D$, but it is not a general assumption.

## Appendix B: Proof of fractional first-order system bandwidth

The transmittance of the asymptotic approximation for the series impedance is:

$$ Z_{asym}(j\omega) = \frac{R}{1 + RC_s\sqrt{j\omega}} \qquad (36) $$

from which one gets the gain:

$$ |Z_{asym}(j\omega)| = \frac{R}{\sqrt{(RC_s)^2\omega + RC_s\sqrt{2\omega} + 1}} \qquad (37) $$

It can be seen that the highest gain in dB for this system is:

$$ |Z_{asym}(j\omega)|_{dB-max} = 20\log_{10}(R) \qquad (38) $$

The -3 dB frequency is given by:

$$ |Z_{asym}(j\omega_b)|_{-3dB} = |Z_{asym}(j\omega)|_{dB-max} - 3dB. \qquad (39) $$

By replacing equation (37) in (39), one gets the polynomial:

$$ (RC_s)^2\omega_b + \sqrt{2RC_s\sqrt{\omega_b}} + (1 - 10^{0.3}) = 0. \qquad (40) $$

Solving the latter gives the positive and real solution for $\omega_b$:

$$ \omega_b \approx \frac{0.2679}{(RC_s)^2} \qquad (41) $$

which is the expression for bandwidth.
<table><caption>TABLE 3. Estimated models</caption><tr><td>Z<sub>1-asymp</sub></td><td><sup>0.5</sup>&frasl;<sub>1 + 0.5&radic;s</sub></td></tr><tr><td>Z<sub>1-approx</sub></td><td><sup>0.5</sup>&frasl;<sub>1 + 0.5&radic;s</sub> <sup>1 + <sup>s</sup>&frasl;<sub>0.2647</sub></sup>&frasl;<sub>1 + <sup>s</sup>&frasl;<sub>0.3787</sub></sub> <sup>1 + <sup>s</sup>&frasl;<sub>3.446</sub></sup>&frasl;<sub>1 + <sup>s</sup>&frasl;<sub>5.73</sub></sub> <sup>1 + <sup>s</sup>&frasl;<sub>23.72</sub></sup>&frasl;<sub>1 + <sup>s</sup>&frasl;<sub>11.21</sub></sub></td></tr><tr><td>Z<sub>3-cap</sub></td><td><sup>1</sup>&frasl;<sub>s</sub></td></tr><tr><td>Z<sub>3-frac1</sub></td><td><sup>1</sup>&frasl;<sub>s</sub> <sup>1</sup>&frasl;<sub>1 + (0.1097s)<sup>1.55</sup></sub></td></tr><tr><td>Z<sub>3-frac2</sub></td><td><sup>1</sup>&frasl;<sub>s</sub> <sup>1</sup>&frasl;<sub>1 + (0.1185s)<sup>1.15</sup></sub> <sup>1</sup>&frasl;<sub>1 + (0.0258s)<sup>1.15</sup></sub></td></tr></table>