John Armitage

Dept. of Physics, Carleton University, Ottawa K1S 5B6, Canada

Concrete, once poured and set, is permanent. If you want to change its shape or properties, it has to be demolished and the work started again. This project arose out of a discussion on what is the best type of concrete to use for the floor of a building, that will facilitate its heating by solar energy. The floor in question lies on top of sub-floor also made of concrete. In between the floor and sub-floor, a set of heating wires of minimal thickness is laid. What properties should the sub-floor have? Should the floor be made of concrete with low or high thermal conductivity? There are arguments in favour of both approaches – low thermal conductivity concrete has a greater thermal mass and might stay warm longer, whereas the high conductivity concrete will transmit more heat to the interior. To examine the advantages or disadvantages of these different materials an electrical model was constructed.

The correspondence between the equations governing heat flow, and those describing the flow of electricity has been known for a long time. There are journal references stretching back over eighty years [1],[2],[3]. Unfortunately, many of those early references are now difficult to access, but a useful survey has been given by Gerald Parnis [4], and this is readily available on the web.

Figure 1. Volume element with heat flows in and out

On the electrical side, Fig 2, we express the electrical resistance: R, in terms of the resistivity r, and substitute it into Ohm’s Law. Both sides of the equation are differentiated by dx to get the spatial derivative of the current. Then, defining the capacitance: C, in terms of the capacitance per unit volume: c, we substitute that into the equation relating charge and voltage on a capacitor and differentiate it with respect to time. Recognizing that electrical current is time derivative of charge we change variables to current, then differentiate again by dx. We now have two equations for the spatial derivative of current, one from the flow of current through a resistor and one from the flow of charge into a capacitor. We eliminate di/dx between them and obtain:

Comparing equation 1 with equation 2, we can draw the analogy between the description of heat flow and that of electricity flow. In this analogy, temperature is the analog of voltage, and the time constant ‘rc’ is then identified with the inverse of the thermal diffusivity:

(3)

One caveat: only the heat conduction process has been considered. This is a linear process whereas convection and radiation are not. While it is possible to model those processes with non-linear components [5],[6], for this study, a linear approximation is adequate. The temperature changes dealt with here are small enough (±5 degrees C) that a linear approximation suffices. Also, the model is used to compare different materials under the same circumstances, rather than trying to determine an absolute measure of the final temperature. The differences between solutions are of the order of one degree C.

where N is the number of stages, and R_SI* is the metric (Système International) ‘R-value’ of the panel or sheet. Substituting typical values for the thermal constants we find resistor values in the milli-Ohm range and capacitances in the Mega Farad range. For example, take a 1m² (area A) block of concrete 1m in depth: - x. From the above formulae 1/k = 1/1.88 = 0.532 m ^oC/W and s = 1.98 x 10⁶ J/m^{3 o}C – see Table 3 for values. For 10 stages (N=10) a resistance of 0.0588 Ohms is required which is not practical as it approaches the resistance of the connecting wire. (24 gauge copper wire has a resistance of 0.025 ohms per foot.) Also, capacitors of 2 x 10⁶ Farads are not practical to build. In mathematical modelling [8],[9],[10], where connecting wires are not required, and size is immaterial, a straightforward conversion of the thermal quantities into electrical ones can be carried out. On the lab bench, reality intervenes.

The association of thermal quantities with electrical ones does not have to be on a ‘one to one’ basis. We are free to add in multiplicative factors m and n as:

For an N=10 stage filter, applied to a 1m² (A) block of concrete 0.1m (x) thick, using m = 1.13 x 10^-7, produces a resistor value of 50 k. Similarly for the same number of stages of the same block, a value of n = 4.21 x 10¹³ produces a capacitance of 0.47nF. There are readily available components close in value to these.

A consequence of this scaling is that the time recorded in the electrical system, t_e, is not the same as time recorded in the thermal system t_th. This can be seen from the above equation relating the electrical time constant: rc, to the thermal time constant of 1/ . Therefore, elapsed time in the two systems is related by:

Hence the time in the thermal system is reduced by a factor of just under a million. A time period of one hour (3,600s) in the thermal system is represented by a period of 0.75 milliseconds (ms) in the electrical system. This greatly facilitates the representation of the current flow on an electronic device: the oscilloscope.

Using this analogy, a mapping of thermal to electrical quantities is given in Table 1.

Table 1. Comparison of thermal and electrical units.

where q = electrical charge and i = electrical current, n, m are the scaling factors: n = 4.21 x 10¹³ and m = 1.13 x 10^-7 in this example, and N is the number of stages.

4.1 Implementation of the electrical circuit.

Using a square wave input, the signal after stage 10, is shown in fig 5. The square wave input (blue) is 11.37 V amplitude, with a duration of 30 milliseconds. In the output signal at stage 10 (red), some of the high frequencies have effectively been filtered out by the network resulting in a measurable rise time from 0% to 100% of 2.5 milliseconds before the output reaches a saturation value of 5.77 V. It then falls back to zero in a similar time before beginning the cycle again. The 30 ms pulse in Fig 5 corresponds to about 40 hours in the thermal system.

In previous publications [2],[3], the time constants were chosen such that the output did not have time

to fall back to zero before the next pulse arrived. Consequently, they implemented circuits to discharge the capacitors after each pulse. As can be seen from the waveforms on figure 5, that is not necessary here.

For comparison with those publications [2],[3], the output of a 20-stage filter is also shown, Fig 6. This shows an 11.3V input pulse of 30 ms duration producing an output pulse that rises during 5 ms to 3.95 V. In particular this result is in approximate agreement those shown by Robertson and Gross [3], based on measurements made of the temperature difference across a firewall. The comparison is only approximate because they used an input pulse that was not flat topped – but had a constant rise to it.

4 . 2 . Modelling of floor section

The dimensions of the floor and sub-floor are 1m² in area, and 10cm thick. In the model, no thickness is assigned to the heating coils – they are just a 2-D layer between the floor and sub-floor. The equivalent electrical circuit is shown in Fig 8. In this diagram, RF and CF designate the resistors and capacitors representing the floor, while R and C designate those representing the subfloor. This allows different materials, with different material properties to be modelled for floor and sub-floor. The sub-floor network is connected to ground (zero Volts) – representing a constant temperature of 0 degrees C. The upper floor surface is connected to a load resistor RA representing the thermal resistance encountered through conduction from the floor surface to the surrounding air. The air itself is held at a contant temperature of 0 ^oC.

The thermal properties of concrete can be modified by the inclusion of additives [13],[14]. The properties used in this comparison are taken from ref. [15] and re-printed here:

Table 2. Physical and thermal properties of concrete. (Ref. [15])

The admixtures with glass and steel fiber are each 0.75% by weight. As is expected, the presence of steel increases the conductivity whereas glass fiber decreases it. What is not perhaps intuitive is that the heat capacity is greater for the steel fiber filled concrete than for the glass fiber filled.

Based on the formulae given above for converting thermal properties to electrical properties, we take three cases:

1. Floor and sub floor made from unmodified concrete (BM)

2. Floor made from glass fiber reinforced concrete (GF) and sub floor made from unmodified concrete (BM).

3. Floor made from steel fiber reinforced concrete (SF) and sub-floor made from unmodified concrete (BM).

4. Floor made from steel fiber reinforced concrete (SF); sub-floor made from glass fiber reinforced concrete (GF)

These four cases are each modelled with two 10-stage filters, one for the floor and one for the sub floor.

Using the thermal values given in Table 2, and using the scaling factors given earlier (n = 4.2 x 10¹³ and m = 1.13 x 10^-7) and taking the case of a 10-stage filter, N=10, we arrive at the following resistor and capacitor values for each of the three scenarios mentioned above:

Table 3. Resistor and capacitor values used to model the four cases.

A useful set of data is provided by Rowlands et.al.[16]. They provide measured and modelled data from the output of a photovoltaic (PV) panel set up in Ottawa, Canada, and orientated at several tilt angles to the incoming radiation. Their range of angles is fairly limited, from 30 to 44 degrees, but they are in a location comparatively close to the test house discussed in this work. Unfortunately, they only provide yearly average PV output values. The usefulness of this data set is that it can be compared with calculations from the Pacific North-West Lab. (PNW) [17] which provides both yearly and daily solar radiation values, but with a specified air attenuation factor of 0.9. This is a very low attenuation factor and measurements made in Arizona [18] indicate that a factor of 0.78 is more applicable there. Arizona has a very dry atmosphere whereas Ottawa has much more humid summers. Therefore, one would expect the air attenuation coefficient to be much larger. Comparison of the data from Rowlands data to the calculations from PNW allows the air attenuation factor to be extracted for which a value of 0.47 was found. This is quite low for air attenuation alone, but it is likely that this includes the effects of reflectance and absorptability on the front surface of the PV cell [19]. Applying this factor to the daily radiation figure from PNW gives us: 0.62kWh/m² and 0.43kWh/m² per day on June 22^nd for tilt angles of 30 degrees and 60 degrees. A representative figure to use for June is therefore 0.5 kWh/m² day (corresponding to a tilt angle of about 45 degrees). Alternately this can be expressed as 1.8MJ per m² per day.

A waveform which approximates the time variation of energy from a PV cell during the course of one day was used. The AWG does have a rectified sin wave generator which would represent 12 hours of daylight quite well, but unfortunately each ‘lobe’ is followed immediately by another. That would be equivalent to the PV cell receiving the same radiation over the course of a night as during the day. Instead, a Gaussian shape was used which did provide a time gap between successive pulses. However, the AWG does not provide separate control over the width and the period. Once the period is fixed the width is also. Setting the period to be the electrical equivalent of one day, or 18ms, sets the Full Width at Half Maximum (FWHM) to be 5ms or 6hrs 40min.

we can take the area under the voltage-time waveform and when that is divided by the resistance, it gives the charge. Then from Table 1, we can set the electrical charge q, according to the thermal energy required:

So the 1.8MJ per day mentioned in the previous section, corresponds to an electrical charge of 0.043C.

We use:

to determine the peak voltage (V_p). Since the FWHM has been set to 5ms, by the frequency, the height (V_p) is then adjusted so that the area under the curve is 0.043C which is equivalent to a thermal energy of 1.8MJ. This sets V_p to be about 10V. The input pulse is shown in Fig 8 below.

Figure 8. Input pulse (blue) injected in the center of two 10-stage RCfilters.

The red pulse shows the waveform obtained at the floor surface. These were recorded for case 1 with the BM concrete used in both sub floor and floor.

where v_i is the voltage recorded at time t_i. The sum was taken over the digitizations that occurred during the ‘daylight’ hours. From this the delay between input and output waveforms was calculated.Then 10 digitizations on either side of t_p were averaged to calculate the peak output voltage.

Table 4. Output voltages and delays for the four Cases

The resistors and capacitors representing the floor were then changed to the values shown in Table 3 for Case 2. The measurements were repeated, and the output voltage and delay recorded. This procedure was repeated using the resistor and capacitor values shown for Case 3 and 4. The results are tabulated in Table 4.

Results were obtained for the temperature (voltage) obtained on the top surface of the main floor, together with the delay between the peak of the applied thermal input and the peak output. Case 1 using the base material for both floor and sub-floor produces an output voltage in the middle of the range, whereas the two configurations with steel fibers in the main floor (Cases 2 and 4) produce the largest output. Case 3 with glass fibers in the main floor reduces the output voltage. The differences are not great: + 1.3 % and -6.1% respectively from the value obtained for the base material. Note that these differences will be magnified if larger pieces of concrete or larger PV cells are used. The results quantify what might be expected from the numerical values for the thermal conductivity.

Since there is little difference between Case 2 and Case 4, where steel fibers are used for the floor in each, the results of this study would imply that it is not worth the extra cost of using glass fibers in the sub floor (Case 4). The most beneficial effect at the lowest cost would be obtained with Case 2.

The minimum delay of the output peak is found in case 3. Using glass fibers in the main floor reduced the delay of the output peak to 0.51ms. Glass fibers also have the minimum heat capacity. Steel fibres slightly reduce the delay in case 2 where the base material is used for the sub floor, but only show a marginal increase in delay (0.02ms) over the BM Case 1. The difference in the delays does not seem to be a large effect. All mixes delay the output by just inder an hour.

In order to give absolute results, a better model of the floor air interface is required that includes both convection and radiation. This can be implemented with non-linear devices like the varister [2],[3].

The method can also be used to model a room – with four circuits in parallel, to represent the walls, floor and ceiling. [20]. This opens up the possibility of studying the effect of windows for example – which may act as a source of radiant heat and also act as a sink for convective heat. Another advantage is that the analogy facilitates the use of more than one source, for example radiant heat coming in through a window as well as underfloor heating. This can be achieved numerically, but with multiple sources and multiple sinks, the calculations become more complex. Indeed, Ljung and Gland [21] indicate severe problems when the number of parameters exceeds 10. Apart from the intricacies of wiring up circuits, this method can handle the complex layouts required. Electrical modelling provided a quantitative answer to thermal problems and is far less costly and time consuming than pouring several different concrete floors.

ACKNOWLEDGEMENTS

The author acknowledges many useful conversations with Federico Fernandez of Algonquin College, Ottawa. I would also like to thank my wife, Beverley, for her patience.

REFERENCES

[1] Paschkis and Baker,‘A method for determining Unsteady-State Heat Transfer by Means of an Electrical Analogy’, Trans. ASME, 64, New York, NY (1942).

[2] Lawson and McGuire, ‘The solution of transient heat flow problems by analogous electrical networks’, Proc. Inst. Mech. Engrs. Part C: Mech. Eng. Science, Vol A1, p167 (1952).

[3] Robertson and Gross, ‘An electrical analog method for transient heat flow analysis’, Institute of Research of the National Bureau of Standards, Vol.61 No.2, p105 (1958).

[4] G. Parnis, ‘Building Thermal Modelling using Electrical Circuit Simulation’, thesis submitted in 2012 to the School of Photovoltaic and Renewable Energy Engineering, University of New South Wales, Sydney, Australia. http://unsworks.unsw.edu.au/bitstreams/68b8c779-35c5-499d-bb5f-08dea3b5a023/download

[5] Khane, Vaibhav, ‘Analogy based modeling of natural convection’ (2009). Masters Theses. 4723.

[6] Ram´ırez-Laboreo, Sagu´es_and Llorente,’Thermal modelling, analysis and control using an electrical analogy’, Conference Paper · June 2014 DOI: 10.1109/MED.2014.6961423 https://webdiis.unizar.es/~ramirlab/papers/2014_MED.pdf

[7] Stephenson and Mitalas, ‘Lumping errors of analog circuits for heat flow through a homogeneous slab’, NRC Publications (1961) https://nrc-publications.canada.ca/eng/view/object/?id=b8123eb2-d6ff-4c15-92db-ca3057df2f92

[8] Merkaj, Dhamo and Kalluci, ‘Thermal Model of a House using Electric Circuits Analogy’, Proc. 10th International Conf. on Smart Cities and green ICT systems (SMARTGREENS2021) ISBN: 978-989-758-512-8.

[9] Bastida et. al. ‘Thermal dynamic modelling and temperature controller design for a house’, 10th Conf on Applied Energy (ICAE2018) Hong Kong, China. Elsevier. Energy Procedia 158(2018) 2800 – 2805.

[10] Tate, Cheneler and Taylor, "Simplified models for heating system optimisation using the thermal–electrical analogy," 2019 25th International Conference on Automation and Computing (ICAC), 2019, pp. 1-6, doi: 10.23919/IConAC.2019.8895035.

[11] Picoscope 2204A 2 channel 10MHz 100MSPS, with built-in Arbitrary Waveform Generator. Pico-Technologies.

[12] TL082: Wide Bandwidth Dual JFET-input operational amplifier. https://www.ti.com/lit/gpn/tl082-n

[13] Malek, Jackowski, Lasica and Kadel, ‘Influence of Polypropylene, Glass and Steel Fiber on the Thermal Properties of Concrete Materials’, Materials 2021, 14, 1888. https://doi.org/10.3390/ma14081888

[14] Nagy, Nehme, Szagri, ‘Thermal Properties and Modeling of Fiber Reinforced Concretes’, Energy Procedia, Volume 78, 2015, pp2742-2747, ISSN 1876-6102, https://doi.org/10.1016/j.egypro.2015.11.616.

[15] Pavlík, Poděbradská, Toman, and Černý ‘Thermal Properties of Carbon- and Glass Fiber Reinforced Cement Composites in High Temperature Range in a Comparison with Mortar and Concrete’, ResearchGate (2002)

[16] Rowlands, Kemery, Beausoleil-Morrison, “Optimal solar-PV tilt angle and azimuth.” Energy Policy 39(2011) p1397 Elsevier.

[17] Buffo, Fritschen and Murphy, ‘Solar radiation on various slopes and latitudes’, USDA Forest Service Research Paper PNW-142 (1972).

[18] Idso, ‘Atmospheric Attenuation of Solar Radiation’, J. Atmos. Sci. 26 p1088, (1969).

[19] Poruba, et.al. ‘Optical absorption and light scattering in microcrystalline silicon thin films and solar cells’, Journal of Applied Physics 88, issue 1, 148-160, 2000

[20] Kreider, Curtiss and Rabl,’Heating and Cooling for Buildings’, Ch8 p388, CRC Press (2010) ISBN 978-1-4398-1151-1 (hardcover : alk. paper).

[21] Ljung and Gland,’On the global identifyability for arbitrary model parameterizations’, Automatica, 30 No2, pp 265-276, 1994

BIOGRAPH Y OF AUTHOR