Acceleration of the Measurement Time of Thermopiles Using Sigma-Delta Control

This work presents a double sliding mode control designed for accelerating the measurement of heat fluxes using thermopiles. The slow transient response generated in the thermopile, when it is placed in contact with the surface to be measured, is due to the changes in the temperature distributions that this operation triggers. It is shown that under some conditions the proposed controls keep the temperature distribution of the whole system constant and that changes in the heat flux at the thermopile are almost instantaneously compensated by the controls. One-dimensional simulations and experimental results using a commercial thermopile, showing the goodness of the proposed approach, are presented. A first rigorous analysis of the control using the Sliding Mode Control and Diffusive Representation theories is also made.


Introduction
Measurement of heat fluxes on surfaces using thermopiles is an extended practice for many applications [1][2][3][4][5][6][7]. When estimating heat fluxes on surfaces by contact, though, these measurements can be slow. For instance, as shown later in this paper, very long time responses are typical in commercial thermopiles and also even with those used in micro or nanoscale calorimetry applications [8,9]. The main reason for these slow response times is that the temperature distribution in the complete system can have large stabilization times depending on the associated thermal structure [6]. The measurement process can distort the temperature in the body to be measured [5] and reaching a steady state is a slow process.
Algorithms and circuits have been proposed to obtain the response of the sensors without having to wait for their thermal stabilization. For example, in Reference [10] prediction algorithms are proposed using the thermal circuit associated to the structure. Additionally, in Reference [11] hybrid methods combining spatial and temporal measurements are used to improve the time response and the accuracy of sensors. Alternatively, it is also possible to obtain quantitative results using non-invasive methods, such as optical measurements [12,13].
This paper focuses on improving the measurement time of heat fluxes using thermopiles placed in contact with surfaces. To this end, we propose applying a sliding mode control to avoid thermal transients by keeping an almost constant temperature distribution in the thermopile. This type of control is extensively used in many applications [14][15][16] mainly due to their robustness against external disturbances and model uncertainties, such as variations of the thermal circuit of the system. In our case the control is implemented as a double-loop circuit based on thermal sigma-delta modulation [17,18]. This modulation provides an implicit analog-to-digital conversion, that makes it suitable for many sensing applications, even in the case of heterogeneous systems [19,20]. Moreover, this approach also inherits other benefits of sigma-delta converters, such as quantization noise shaping [21].
The aim of the control is to avoid the thermal transients in the thermopile. As it will be explained later in more detail, this is done by keeping constant the state variables of the thermopile, in what we call constant state operation (CSO). This implies that any external variations, such as the heat fluxes generated by the environment on the thermopile, must be compensated by the actions of the control, which constitute the output of the sensor. To achieve the CSO condition, the control has two specific objectives: • To ensure that at the moment of placing the thermopile on the surface to be measured, the contact surfaces are at the same temperature.

•
To guarantee that the spatial temperature distribution in the thermopile remains constant in time once the surfaces are in contact.
The mathematical analysis included in the paper is necessary to understand what are the capabilities and limits of the controls being proposed. For many applications, it is no longer possible to dissociate the design of the sensor from the controls being implemented during its operation. The reason is that sensor performance without using the controls would not be within the expected operational limits. For example, our research group has participated in three instruments for three NASA missions to Mars: Rover Environmental Monitoring Station, REMS (Mars Science Laboratory, 2011), TWINS (InSight, 2018) and Mars Environmental Dynamics Analyzer, MEDA (for Mars2020). Constant temperature operation has been applied in all these sensors [22,23]. The mathematical thermal modeling and analysis of the control dynamics was presented in Reference [24]. Open loop operation (constant power operation in those cases) would provide unreasonable time responses.
Additionally, controls based on sigma-delta modulation have been also proposed for MOX gas sensors [21,25] and a first mathematical analysis has been presented in Reference [26]. This means that the proposed techniques are transversal and can be applied to a number of different sensor fields. In all cases it is necessary to analyze the control dynamics from a mathematical point of view, and sensor operation and even design cannot be dissociated from this analysis. This work applies these techniques for the first time to heat flux sensors. This paper is organized as follows. Section 2 explains the proposed control, which is analyzed in Section 3 using the Diffusive Representation and Sliding Mode Control theories. Section 4 presents and discusses 1D simulations. Finally the experimental results obtained with a commercial thermopile are presented and discussed in Section 5.

Proposed Controls
The double control proposed is based on sigma-delta modulation. Sigma-delta modulators are discrete-time systems in which the control variables are changed periodically, with period T S . In our case, two signals are monitored at sampling times nT S : the voltage generated by the thermopile, proportional to the temperature difference between its top and bottom surfaces: T g (nT S ) = T top (nT S ) − T bottom (nT S ); and the average temperature, defined as T a (nT S ) = 1/2(T top (nT S ) + T bottom (nT S )), see Figure 1. Note that T a is simply the average temperature of the top and bottom surfaces and therefore it must not generally coincide with the average of the temperature distribution in the thermopile.
Additionally, when some current flows between the two terminals of the thermopile, two separate effects are generated [27,28]: • Peltier effect: current flowing in one direction generates heat injection at one of the thermopile surfaces and heat extraction at the other. The heat changes sign when the current is reversed. This is the main effect in a thermopile and the heat injected/extracted at the interfaces is proportional to the value of the applied current.
• Joule effect: this is normal Joule heating and therefore it is proportional to the square of the applied current. This is a volumetric phenomenon, which happens on the bulk of the thermopile. Joule effect obviously does not produce any cooling, only heating. According to all this, the objective of the control is to enforce the conditions T g (nT S ) = T th g and T a (nT S ) = T th a , where T th g and T th a are previously chosen target values, by applying adequate Peltier and Joule excitations. To this effect, depending on the instantaneous value of T g (nT S ) the T a (nT S ), the control applies to the thermopile one of the current waveforms shown in Figure 2 during the next sampling period, t ∈ [nT S , (n + 1)T S ).
In the next section it is shown that, if the sampling period T S is chosen to be much smaller that the relevant time constants in the structure, the effective excitations that these current waveforms produce are: • a Peltier excitation, which is proportional to the average value of the current density, J and • a Joule excitation, which is proportional to the average value of the power dissipated, that is, proportional to J 2 .
Then, by selecting between the two left or two right, waveforms in Figure 2 it is possible to increase or decrease the average current density J in the thermopile. Since J can be positive and negative (left or right excitations) this allows to heat or cool down the top surface (and cool or heat the bottom). Additionally, by choosing between the upper or lower waveforms, it is possible to increase or decrease the average value of the square of the current density, J 2 , which is equivalent to increase or decrease Joule heating in the thermopile.
The application of the four waveforms, governed by the double sigma-delta control loop, will allow to control the desired temperatures, T a and T g . Table 1 summarizes some of the works in the literature on the measurement of heat flux and comments on the main differences between the methods proposed in them with regard to this work. To the best knowledge of the authors, there are no works on control or feedback methods in which the thermopile itself is used as sensing and actuation. The works in References [3,29] proposed active methods for the cancellation of the heat flux across the thermopile but using external heaters. The time responses in these works are very small in some cases but most of them involve no contact (radiation heating) and/or micromechanized sensors.

Analysis of the Control
The first step towards understanding the proposed controls is to analyze their dynamics when they are applied to the thermopile before placing it on the surface. To this purpose, let us consider the case described in Figure 3, in which the thermopile is suspended on air, with negligible radiation losses. In this case all heat is lost through convection at the surfaces of the thermopile, modeled as two thermal resistances R cv . For the sake of simplicity, let us use T ∞ , the temperature of the air, as temperature reference. In this case the heat equation is: where T(x, t), x ∈ [0, L], is the temperature and C p , ρ p and k p are the thermal capacitance, the density and the thermal conductivity of the thermopile, respectively. We define: as the term related to the Joule effect, being p the electrical resistivity of the thermopile and J(t) the current density applied. The boundary conditions at the thermopile surfaces x = 0 and x = L are: where the left terms correspond to the heat losses in the surfaces of the thermopile and g 0 S (t) = α p J(t)T(t, 0 + ) and g L S (t) = α p J(t)T(t, L − ) are the terms related to the Peltier effect, which is proportional to the current applied J(t) [28]. There, α p is the Seebeck coefficient of the thermopile. This formulation based on the Seebeck coefficient is similar to that used in other works [36,37].
Additionally, in our case the temperature variation in the thermopile is very small when compared to its absolute temperature and therefore we consider g 0 In the Laplace domain, Equations (1)-(4) become: Now, defining the average temperature, T a (t) := 1/2(T(t, L) + T(t, 0)) and the gradient temperature, T g (t) := T(t, L) − T(t, 0), together with their respective Laplace transforms, T a (s) and T g (s), it is possible to find that: This result implies that the Joule input changes only the average temperature, whereas the Peltier input produces only variations on the gradient temperature. In both cases, the poles are simple, real and non-positive. Therefore, the impulse responses of (10) are: η a n e −ω a n t t > 0 where −ω a n and −ω g n are the poles of H a S (s) and H g J (s). Besides, {η a n } and {η [38]. These diffusive symbols are the pole residues of (10). For n 1, the pole frequencies can be approximated as: and the values of η a n are positive and decreasing as 1/n 2 . Additionally, η g n are also positive but tending to a constant value 8r/Lk p for large n values.
The exact values of the poles of (10) and their corresponding residues have been computed to obtain the diffusive symbols η a n and η g n shown in Figure 4. The thermopile parameters listed in Table 2 and R cv = 25 mΩ have been used in the calculations. The table also includes the thermal parameters of a solid used later in the simulations. As expected, all the symbols are real and positive, with values of η a n decreasing as 1/n 2 with frequency. Also as expected, η g n rapidly tends to 8r/Lk p = 6.236×10 −4 K m 2 Hz/W. Table 2. Parameter values used in the theory and simulations.

Diffusive Representation of the Operators
As it has been shown, there are two impulse responses, h a J (t) and h g S (t), with their respective sets of poles {−ω a n } and {−ω g n }, n > 0. In order to analyze the dynamics of the system under the infinite sampling approximation, it is convenient to use the Diffusive Representation theory [38]. Taking into account that the impulse responses present negative simple real poles, we can define a function ψ a n (t), associated to each pole −ω a n , in such a way that: ψ a n = −ω a n ψ a n (t) + g J (t) T a (t) = ∑ ∞ n>0 η a n ψ a n (t) (14) In the same way, for the gradient temperature, T g , we may define functions ψ g n (t): with g J (t) and g S (t) defined in Equations (2) and (5). It must be noted that ψ a n and ψ g n are the amplitude of the excited modes of the heat equation and therefore they can be seen as state variables. This is in fact the standard approach in Diffusive Representation theory. Since the modes that both g S and g J excite are different, the state variables of both systems are independent. Joule actuation activates only symmetric temperature distributions, whereas Peltier actuation activates only antisymmetric temperature distributions.
On the other hand, the actuation signals, g J (t) and g S (t), are linked, because both depend on the applied current, J(t). For example, it is not possible to have g J (t) = 0 with g S (t) = 0.
Simple Diffusive Systems-that is, systems with numerable negative single real poles-can be approximated, with arbitrary accuracy, using a finite number of poles [38]. Therefore, from here on we will consider a finite set of poles {−ω a n } and {−ω g n } for 0 < n ≤ N. For convenience, we will assume the same cardinality in both cases.

InfinitE Sampling Frequency Approximation
As it has been explained in Section 2, the control applies the waveforms in Figure 2 depending on the instantaneous values of T a (t) and T g (t), which are the average and the difference, respectively, of the temperatures at the top and bottom surfaces of the thermopile. Since the control is discrete in time, control decisions are taken at multiples of the sampling period, T S and hence, T a and T g are discrete variables.
Following the same approach as in Reference [19], we will now analyze the control dynamics under the infinite sampling approximation. Under this approximation we consider that T S → 0 and that therefore the system only responds to the average value of the excitation waveforms, which have an infinitesimal duration. We may then define: Depending on the value of these variables, the control will apply the waveforms in Figure 2. In this case, T a (t) and T g (t) are continuous in time because of the infinite sampling approximation. Now we can define the average Joule excitation applied by the control, u a (t) = p J 2 (t), as: (17) and the average Peltier excitation, u g (t) = α p J(t)T ∞ , as: Then, under the infinite sampling approximation, systems in (14) and (15) become: ψ a n = −ω a n ψ a n (t) + u a (t) T a (t) = ∑ N n>0 η a n ψ a n (t) In general we will consider that G on + > G off + > 0 and G on -< G off -< 0. In the case of the average excitation, we will consider that G on J > G off J .

Sliding Mode Analysis of the Average Control
Now we will analyze the dynamics of the system described in (19) and (20), with the average actuations defined in (17) and (18), which are provided by the waveforms generated by the control. As it has been mentioned before, the state variables of the system are ψ a n (t) and ψ g n (t), n = 1, · · · , N. Therefore we can generate two vectors, Ψ a = (ψ a 1 , · · · , ψ a N ) T and Ψ g = (ψ g 1 , · · · , ψ g N ) T , representing the state of the complete system at any given time.
Since the objective of the controls is to set values to T a and T g we may speak of two control surfaces: S a (T th a ) := {Ψ a ∈ R N : σ a (Ψ a ) = ∑ n η a n ψ a n − T th a = 0} and S g (T th It is desired then that the control places Ψ a and Ψ g on the respective surfaces, which is equivalent to controlling T a and T g .

Attractive Set on S a
The control surface for the average temperature is S a (T th a ) := {Ψ a ∈ R N : σ a (Ψ a ) = ∑ n η a n ψ a n − T th a = 0}, with Ψ a = (ψ a 1 , · · · , ψ a N ) T . We will generally assume that G off J < G on J . The first result is that the set: is forward invariant. This follows from having single negative poles with input u a (t) ∈ {G off J , G on J }. On the other hand, we have that: where Γ a = ∑ N n η a n . Now, the intersection of the following subset of the state space: n η a n ω a n ψ a n < G on J with the control surface S a , that is, Ω a ∩ S a , is attractive. This is due to the fact that if Ψ a (t) ∈ Ω a and σ(t) < 0 we will haveσ(t) > 0, whereas if σ(t) > 0 then we will haveσ(t) < 0. Therefore it is σσ < 0 in a neighbourhood of Ω a ∩ S a , excluding the control surface, which means that the set is attractive.
, there is a maximal interval I a ⊂ R such that, for all T th a ∈ I a , the set S a (T th a ) ∩ Ω ∞ a is an attractive non empty sliding set.
we have that S a (T th a ) ∩ Ω ∞ a is a non empty attractive sliding set. The interval is maximal since for any T th a outside the interval (24) we have that S a (T th a ) ∩ Ω ∞ a = ∅.

Equivalent Control
We now may calculate the equivalent control u eq a (t): σ a (t) = N ∑ n η a n (ω a n ψ a n (t) − u eq a (t)) = 0 and therefore: u eq a (t) = Γ −1 a N ∑ n η a n ω a n ψ a n (t) (26) The global dynamics on the control surface is then: ψ a n (t) + ω a n ψ a n (t) = Proposition 2. If η a ∈ R N + \ {0} and T th a ∈ I a the system presents a unique stable equilibrium point at S a (T th a ) ∩ Ω ∞ a . (27) can be rewritten asψ = Mψ with:

Proof. The system in
where diag(ω 1 , · · · , ω N ) is the diagonal matrix with components ω n and 1 N = (1, · · · , 1) T ∈ R N . M is a proper Metzler matrix because its off diagonal elements are nonnegative and at least there is a nonzero entry. This means that for some > 0, matrix A = M + I is a positive matrix and we can then apply the Perron-Frobenius Theorem: the eigenvalue, λ 0 , of largest absolute value of a positive square matrix, A, is both simple and positive and belongs to a positive eigenvector, f ∈ R N + \ {0}. All other eigenvalues are smaller in absolute value [39].
Matrix M has the same eigenvalue spectrum of A displaced to the left . This implies that M has a real simple eigenvalue λ 0 = λ 0 − such that for any other eigenvalue, λ, it is Re(λ) < λ 0 [39].
It is easy to check that for x = (ω −1 1 , · · · , ω −1 N ) T , we have Mx = 0. Now, we will follow a procedure similar to that of [40]. If there is an eigenvalue, λ, with Re(λ) > 0 then λ 0 ≥ Re(λ) > 0 and therefore there is a left eigenvector, f ∈ R N + \ {0}, such that: and since x, f ∈ R N + \ {0} we have that λ 0 = 0. The left eigenvector of the zero eigenvalue is obviously (η 1 , · · · , η N ) and, therefore, the trajectory followed by the system is on the control surface S a (T th a ). This trajectory will tend to an asymptotic stable equilibrium on S a (T th a ), driven by all the remaining eigenvalues of the spectrum of M, since all of them have a negative real part.

Sliding Mode Analysis of the Gradient Control
The control surface for the gradient temperature is S g (T th We will generally assume that G on + > G off + > 0 and G on -< G off -< 0. The gradient temperature actuation (18) can be rewritten as: where λ(t) reflects the state of the average temperature control: Remark 1. The gradient control uses switching between G on + |G on or G off + |G off depending on the state of the average control. This means that the effective actuation of the gradient control depends on the instantaneous value of the actuation determined by the average control. This is what the terms in (31) reflect.
Let us now define the sets: where Γ g = ∑ N n η g n .

Proposition 3.
If η g ∈ R N + \ {0}, there is an interval I g ⊂ R such that, for all T th g ∈ I g the set S g (T th g ) ∩ Ω r g is an attractive non empty sliding set.
then, following the same approach as in Proposition 1, S g (T th g ) ∩ Ω r g is a non empty attractive sliding set.

Equivalent Control
Under the conditions of Proposition 3, we now may calculate the equivalent control u The dynamics on the control surface is then: Proposition 4. If η g ∈ R N + \ {0} and there is an equilibrium point of the system in S g (T th g ) ∩ Ω r g , it will be unique and asymptotically stable.
Proof. It can be proven using the same approach as in Proposition 2.
This proposition shows that, regardless of the state of the average control, it is possible to guarantee under some conditions the stability of the system in a set, Ω r g , smaller than the one that can be reached under more favourable conditions of the average control.

Discussion
The results in Section 3.4 show that, regardless of the state of the average control, it is possible to guarantee under some conditions the stability of the system in a set, Ω r g , smaller than the one that can be reached under more favourable conditions of the average control. The dynamics of the thermopile while it is in contact with the surface is more complex and subject of future work. In order to understand what is to be expected in the experimental measurements, we will proceed with simulations.
In this Section we have analyzed the dynamics of the control before the thermopile is placed on the surface. The fact that the structure is symmetrical (convection on the two sides of the thermopile) implies that the Joule excitation generates symmetrical modes of the heat equation, whereas the Peltier excitation generated antisymmetric modes. Under this approach it has been possible to obtain rigorous results with regard to the stability of the average and gradient temperature controls.
The dynamics of the thermopile while it is in contact with the surface will be more complex and subject of future work. In order to understand what is to be expected in the experimental measurements, we will proceed with simulations in the next Section.

Simulation Results
To evaluate the proposed controls, discrete-time simulations with the 1D structure shown in Figure 5 have been performed. The solid corresponds to the segment x ∈ [0, L 1 ]. The heat fluxQ is injected at x = 0, whereas x = L 1 represents the surface on which the thermopile will be placed. We consider that this last surface dissipates heat initially through convection, R cv , before placing the thermopile. Once the thermopile is on the surface, convection takes place through the other side of the thermopile, x = L 2 , see Figure 5. The heat equation including the Joule effect is: where C 1 the thermal capacitance, ρ 1 the density and k 1 the thermal conductivity of the solid. Three different boundary conditions are applied. The first one is: whereQ(t) is the heat flux being injected from the maximum depth in the material conforming the surface. In general we will assume thatQ(t) is constant or presents very slow time-variations. The second boundary condition is at the interface between the solid and the thermopile: Finally, the third boundary condition is at the interface between the thermopile and the convection thermal resistance: where h is the convection coefficient and T ∞ the temperature of the air at infinity. Table 2 summarizes the values of all parameters used in the simulations. The thermal parameters of the solid are those of the polypropylene used later in the experiments. The parameters of the thermopile are those of the commercial sensor used. In the simulations, the Crank Nicolson finite difference method was used for solving the heat equations [41]. The solid and the thermopile segments were discretized in 100 points each.

Results and Discussion
The first simulation focuses on the evolution of the temperatures in the structure without applying any control. The initial temperature distribution is the dotted line of Figure 6. An initial difference of 0.05 K between the temperature of the surface and that of the thermopile has been introduced. The solid line of Figure 6 is the temperature distribution when, after approximately 120 min, a stable regime is reached. As it can be observed, there is no appreciable temperature discontinuity at the solid-thermopile interface. The temperature slope within the thermopile is noticeably smaller, since k p > k 1 . Figure 6. Temperature distribution in the 1D simulation domain at two different moments when no control is applied. The heat flux injected at the left of the structure is 0.1 W/m 2 . The initial temperature difference between the thermopile and the surface is 0.05 K. Steady-state is reached when the temperature distribution at the solid-thermopile interface is almost continuous.
One can conclude that the time necessary for measuring the heat flux on a surface can be generally quite large, due to the slow transients in the temperature distributions both in the surface (and underlying materials) and in the thermopile. Now, with regard to the control and as discussed in Section 2, the objective of the controls is to set T(x, t) = T const , with x ∈ {L 1 , L 2 }. Additionally, the temperature at the surface of the solid and at the contact surface of the thermopile should be the same, T(L − 1 ) = T(L + 1 ) = T const at the moment of placing the thermopile in contact with the solid. To this effect, two consecutive stages can be contemplated: • Setup: Before placing the thermopile. Controls are active. • Operation: once the temperature of the thermopile is the same as that of the surface, the thermopile is placed on the surface.
In the control simulations, the target average temperature T th a is set slightly above the temperature of the solid surface, therefore, near saturation. On the other hand, the target of the gradient control, T th g is set to zero. The top graph in Figure 7 shows the time evolution of T g , being t = 0 the moment at which the thermopile contacts the surface of the solid. It can be seen that the gradient temperature fluctuates around zero with no appreciable initial transient. As it is also shown in Figure 7, the average temperature increases until approximately t = 17 min, where T th a is reached. Additionally, the bottom graph if Figure 7 shows the evolution of the gradient control flux. The time response is almost instantaneous (≈6 s). The main function of the average temperature control is to set the thermopile at the same temperature as that of the surface to be measured. The gradient temperature control, on the other hand, is used to enforce an almost null temperature gradient in the thermopile. The simulations show that the controls allow keeping constant the temperature distributions in the whole structure: solid-thermopile. When the thermopile contacts the surface, there is a sudden change in the heat injection at the bottom of the thermopile, which is almost immediately canceled by the control.

Effect of Contact Resistance and Side Heat Losses
Let us now investigate the effect of two non-ideal behaviours typically present while measuring heat fluxes using thermopiles in contact with the surface to measure:

•
Contact resistance: depending on the roughness of the surfaces in contact, or the applied force, a thermal contact resistance R c may appear between the thermopile and the surface of the solid. This effect is taken into account by modifying the boundary condition (39) as: Any heat flux across the interface generates a temperature discontinuity drop.
• Side heat losses: this is taken into account in 1D models by adding a −h L (T − T ∞ ) term, accounting for the side losses to the heat equation on the thermopile (37): In the simulations we have used R c = 0.1 K/(W/m 2 ) and h L = 500 W/(m 3 K). The result of simulating the structure in open loop configuration can be observed in Figure 8. The contact resistance generates a temperature drop at the interface of approximately 0.01 K. On the other hand, the asymptotic value reached by the heat flux, as estimated in the thermopile, is 0.097 W/m 2 , which is slightly below the desired value, 0.1 W/m 2 , due to the side heat losses in the thermopile.  Figure 9 shows the time evolution of the average temperature in the thermopile and the estimated heat flux from a simulation using the control. The main difference in the final temperature distribution is that now there is no temperature gradient in the thermopile. The dynamical response is no longer as fast as it was with no contact resistance but it is still much faster than in the open loop case.

Experimental Results
In order to test the proposed control strategy an experimental setup including a polypropylene cylinder of 15 cm diameter and 10 cm high as the solid whose surface heat flux is to be measured has been used. Two class-A accuracy PT500 resistors continuously monitor the temperatures at the top and bottom surfaces of the cylinder. The cylinder is placed on top of a Printed Circuit Board (PCB) in which a copper resistor has been drawn, which works as a constant heat source.
The thermopile used is a HFP01 commercial heat flux sensor from Hukseflux, which has a sensitivity S = 60 µV/(W/m 2 ), a thermal resistivity R th = 0.0071 K/(W/m 2 ), 8 cm 2 of sensing area and typical response times around 180 s. A third class-A PT500 allows to monitor the temperature at the bottom-contact surface of the thermopile, as seen in Figure 10. In the experiments, the thermopile has been placed directly on top of the cylinder, with no other contact forces than its own weight. To implement the controls, only the temperature of the solid surface and the temperature at one of the surfaces of the thermopile are needed. The sensor at the bottom of the cylinder is used only to have an independent measurement to estimate the heat flux along the cylinder. Finally, 1D simulations have been made to estimate the side heat losses in the cylinder. Assuming a convection heat coefficient of 2 W/(m 2 K) (value obtained from Reference [42] in the case of natural convection between a person and the surrounding room air), the simulations show that the difference between the heat fluxes estimated by the PT500s measurements and the one that the thermopile should be having is in the order of 10%. For the heat fluxes measured, this difference falls in the tolerance range of the PT500 resistors used. It must be noted that more than about half the power generated at the PCB is lost to the supporting structure. The thermopile is placed in a rest position and then, at approximately t = 3000 s, it is placed on the cylinder. As it can be observed, there is a slow transient in the estimated heat flux that is stabilized after approximately 30 min. Both the heat flux estimated by the thermopile and that estimated by the PT500s placed on the cylinder are very similar and around 2.5 W/m 2 .

Closed Loop Experiment
In the second experiment the thermopile was under closed loop actuation with the proposed control. The current value used is J 0 = 0.42 mA/cm 2 and the sampling frequency is 1/T S = 1.5 KHz. Each sampling period is divided in 256 time intervals, allowing the generation of the waveforms shown in Figure 2. In this case P + 1 = 215/256T S , P + 2 = 235/256T S and P − 1 = 40/256T S , P − 2 = 45/256T S . The thermopile was placed on top of the cylinder at approximately t = 907 s and removed at t = 1430 s. The target values of the controls are T th g = 0 and T th a = T top for the whole experiment. This means that the target temperature for the top and bottom surfaces of the thermopile are equal to the temperature at the surface of the cylinder. The heat flux is calculated, as in the case of the open loop experiment, directly from the sampled output voltage of the thermopile. Figure 12a shows the instantaneous heat flux in the thermopile near the moment the thermopile is placed on top of the cylinder, t = 907 s. The heat flux is successfully kept near zero, showing the typical chattering signal coming from a sliding control based on sigma-delta modulation. These quick variations are due to the small but fast changes in the heat flux of the thermopile generated by the application of the different current waveforms, see Figure 2. Figure 12b shows the temperatures measured by the PT500s placed at the top surface of the cylinder and the bottom surface of the thermopile. There is a small change in the waveforms at approximately t = 907 s.  Figure 13 shows the evolution of the control heat flux, calculated with α p = 8.5 mV/K and A = 8 × 10 −4 m 2 . In the thermopile, we define the control heat flux as negative when the current applied by the control removes heat from its bottom surface and injects heat at its top surface. Then, when at t = 907 s the thermopile is placed on top of the cylinder, it suddenly receives the heat flux coming from the cylinder. This heat flux is compensated by the control, by cooling the bottom surface of the thermopile, that is, generating a negative heat flux. Besides, before placing and after removing the thermopile the control had to compensate an average heat flux from the environment of −0.25 W/m 2 but between t = 907 s and t = 1460 s the estimated heat flux, approximately 1.5 W/m 2 , is very close to the heat flux in the cylinder calculated from the PT500s at its top and bottom surfaces.

Conclusions
The objective of the double control proposed is to avoid slow thermal transients by achieving the CSO condition in the thermopile. This implies keeping a constant temperature distribution in the whole structure once the thermopile is placed in contact with the surface to measure. The gradient control ensures that the difference between the temperatures at the bottom-contact and top surfaces of the thermopile is constant. The average control is set to ensure that the bottom surface is at the same temperature as the surface to measure. An analysis of the dynamics of the control using the Sliding Mode Control and Diffusive Representation theories has been performed. Finally, simulations and experimental results, in which response times of a few seconds have been obtained, confirmed the goodness of the proposed approach.