An innovative method to determine optimum insulation thickness based on non-uniform adaptive moving grid

TL;DR

This paper introduces a novel non-uniform adaptive method for accurately determining optimal building insulation thickness, significantly reducing computational time while maintaining high accuracy in energy simulations.

Contribution

It presents an innovative adaptive technique that improves simulation efficiency and accuracy in optimizing insulation thickness for building envelopes.

Findings

Achieves 25% reduction in simulation run time.

Provides high accuracy with fewer spatial nodes.

Efficiently simulates heat conduction in building envelopes.

Abstract

It is well known that thermal insulation is a leading strategy for reducing energy consumption associated to heating or cooling processes in buildings. Nevertheless, building insulation can generate high expenditures so that the selection of an optimum insulation thickness requires a detailed energy simulation as well as an economic analysis. In this way, the present study proposes an innovative non-uniform adaptive method to determine the optimal insulation thickness of external walls. First, the method is compared with a reference solution to properly understand the features of the method, which can provide high accuracy with less spatial nodes. Then, the adaptive method is used to simulate the transient heat conduction through the building envelope of buildings located in Brazil, where there is a large potential of energy reduction. Simulations have been efficiently carried out for different wall and roof configurations, showing that the innovative method efficiently provides a gain of 25% on the computer run time.

Tables

Table 1. Table 1. Computational time required for the numerical schemes to perform the nonlinear case for ε ∞ ≅ 𝒪 ​ ( 10 −  4 ) subscript 𝜀 𝒪 superscript 10 4 \varepsilon_{\,\infty}\cong\mathcal{O}\,(10^{\,-\,4})\, .

Simulation time	Crank-Nicolson	QUNT	IMEX
$2$ days	$1.43𝗌$	$1.12𝗌$	$1.08𝗌$
$1$ month	$21.36𝗌$	$16.2𝗌$	$15.2𝗌$
$1$ year	$259.32𝗌$	$205.1𝗌$	$185.9𝗌$
Spatial nodes	$1001$	$51$	$501$
CPU time ($\%$)	$100\%$	$75\%$	$72\%$

Table 2. Table 2. Climate and geographic characteristic of the Brazilian cities.

City	$T_{\min }$ $\left[{}_{}{}^{\circ }𝖢\right]$	$T_{\text{mean}}$ $\left[{}_{}{}^{\circ }𝖢\right]$	$T_{\max }$ $\left[{}_{}{}^{\circ }𝖢\right]$	Altitude $\left[𝗆\right]$	Location	Clim. zone
Curitiba	$-2$	$16.3$	$30.9$	$935$	${25.5}^{\circ }$ S ${49.3}^{\circ }$ W	$3$
Rio de Janeiro	$13$	$23.5$	$38.2$	$2$	${22.9}^{\circ }$ S ${43.2}^{\circ }$ W	$1$
São Paulo	$7.5$	$18.8$	$32.8$	$760$	${23.5}^{\circ }$ S ${46.6}^{\circ }$ W	$2$
Salvador	$14.2$	$25.3$	$33.5$	$8$	${12.9}^{\circ }$ S ${38.5}^{\circ }$ W	$1$

Table 3. Table 3. Thermophysical properties of the materials [ 27 ] .

Material	${\rho }_{\mathrm{\hspace{0.17em}0}}\left[{𝗆}^{\mathrm{𝟥}}/\mathrm{𝗄𝗀}\right]$	$c_p\left[𝖩/\mathrm{𝗄𝗀}/𝖪\right]$	$k\left[𝖶/𝗆/𝖪\right]$
Brick	$1800$	$840$	$0.69$
Concrete	$2200$	$840$	$2$
Extrude Polystyrene	$25$	$1470$	$0.0275$

Table 4. Table 4. Economic parameters.

Insulation price	Heating/cooling system efficiency	Energy price
$p_i$	$\eta$	$p_e$
$100\$/{𝗆}^{\mathrm{𝟥}}$	$0.8$	$0.218\$/\left(\mathrm{𝗄𝖶𝗁}\right)$

Table 5. Table 5. Optimum insulation thickness according to wall orientation for inside (or outside) insulation.

Façade orientation	Curitiba	Rio de Janeiro	Salvador	São Paulo
South	$5\mathrm{𝖼𝗆}$	$4\mathrm{𝖼𝗆}$	$5\mathrm{𝖼𝗆}$	$3\mathrm{𝖼𝗆}$
North	$4\mathrm{𝖼𝗆}$	$5\mathrm{𝖼𝗆}$	$6\mathrm{𝖼𝗆}$	$2\mathrm{𝖼𝗆}$
West	$4\mathrm{𝖼𝗆}$	$5\mathrm{𝖼𝗆}$	$6\mathrm{𝖼𝗆}$	$2\mathrm{𝖼𝗆}$
East	$4\mathrm{𝖼𝗆}$	$4\mathrm{𝖼𝗆}$	$6\mathrm{𝖼𝗆}$	$3\mathrm{𝖼𝗆}$

Table 6. Table 6. Optimum insulation thickness according to roof insulation position.

Insulation position	Curitiba	Rio de Janeiro	Salvador	São Paulo
Outside	$3\mathrm{𝖼𝗆}$	$6\mathrm{𝖼𝗆}$	$7\mathrm{𝖼𝗆}$	$2\mathrm{𝖼𝗆}$
Inside	$4\mathrm{𝖼𝗆}$	$6\mathrm{𝖼𝗆}$	$7\mathrm{𝖼𝗆}$	$3\mathrm{𝖼𝗆}$

Nomenclature
	Latin letters
$c_p$	material specific heat capacity	$\left[𝖩/\left(\mathrm{𝗄𝗀}\cdot 𝖪\right)\right]$
$c$	energy storage coefficient	$\left[𝖩/\left({𝗆}^{\mathrm{𝟥}}\cdot 𝖪\right)\right]$
$E$	heat transmission load	$\left[\mathrm{𝖬𝖩}/{𝗆}^{\mathrm{𝟤}}\right]$
$h$	convective heat transfer coefficient	$\left[𝖶/\left({𝗆}^{\mathrm{𝟤}}\cdot 𝖪\right)\right]$
$k$	thermal conductivity	$\left[𝖶/\left(𝗆\cdot 𝖪\right)\right]$
$l$	wall length	$\left[𝗆\right]$
$q$	heat flux	$\left[𝖶/{𝗆}^{\mathrm{𝟤}}\right]$
$T$	temperature	$\left[𝖪\right]$
$t$	time	$\left[𝗌\right]$
$x$	thickness coordinate direction	$\left[𝗆\right]$
	Greek letters
$\alpha$	solar absorptivity of outside surface wall	$[-]$
$ϵ$	emissivity	$[-]$
$\rho$	material density	$\left[\mathrm{𝗄𝗀}/{𝗆}^{\mathrm{𝟥}}\right]$
$\sigma$	Stefan-Boltzmann constant	$\left[𝖶/\left({𝗆}^{\mathrm{𝟤}}\cdot {𝖪}^{\mathrm{𝟦}}\right)\right]$
	Dimensionless parameters
$\mathrm{Bi}$	Biot number	$[-]$
$c^{\star }$	storage coefficient	$[-]$
$\mathrm{Fo}$	Fourier number	$[-]$
$k^{\star }$	diffusion coefficient	$[-]$
$q_{\infty }^{\star }$	short-wave radiation source term	$[-]$
$R_{\text{lw}}^{\star }$	long-wave radiation coefficient	$[-]$
$u$	temperature field	$[-]$

Equations

$$\displaystyle c\ \dfrac{\partial T}{\partial t}\ =\ \dfrac{\partial}{\partial x}\Biggl{(}\,k\ \dfrac{\partial T}{\partial x}\,\Biggr{)}\,,$$

$$T(x=0,t)=T_L(t)$$

$$\displaystyle-\,k\,\dfrac{\partial T}{\partial x}\biggr{|}_{x\,=\,l}$$

$$\displaystyle k\,\dfrac{\partial T}{\partial x}\biggr{|}_{x\,=\,0}$$

$$T(x,t=0)$$

$$\displaystyle q_{\,j}\,(t)\mathop{\stackrel{{\scriptstyle\ \mathrm{def}}}{{:=}}\ }-\,k\,\dfrac{\partial T}{\partial x}\,(x\,,t)\biggr{|}_{\,x\,=\,x_{\,j}}\,.$$

$$E={\int }_{t_0}^{t_1}q(\tau )\mathrm{d}\tau ,$$

$$\displaystyle u\$$

$$\displaystyle c^{\,\star}\$$

$$\displaystyle\mathrm{Fo}\$$

$$c^{\star }\frac{\partial u}{\partial t^{\star }}$$

$$u(x^{\star }=0,t^{\star })$$

$$u(x^{\star }=1,t^{\star })$$

$$\displaystyle k^{\,\star}\ \dfrac{\partial u}{\partial x^{\,\star}}\biggr{|}_{x^{\,\star}\,=\,0}$$

$$\displaystyle-\ k^{\,\star}\ \dfrac{\partial u}{\partial x^{\,\star}}\biggr{|}_{x^{\,\star}\,=\,1}$$

$$u(x^{\star },t^{\star }=0)$$

$$q_j^{\star }(t^{\star })$$

$$c(u)\frac{\partial u}{\partial t}$$

$$\displaystyle c_{\,j}^{\,n}\,\dfrac{u_{\,j}^{\,n+1}\ -\ u_{\,j}^{\,n}}{\Delta t}\ =\ \dfrac{1}{\Delta x}\,\Biggl{[}k_{\,j+{\frac{1}{2}}}^{\,n}\,\dfrac{u_{\,j+1}^{\,n+1}\ -\ u_{\,j}^{\,n+1}}{\Delta x}\ -\ k_{\,j-{\frac{1}{2}}}^{\,n}\,\dfrac{u_{\,j}^{\,n+1}\ -\ u_{\,j-1}^{\,n+1}}{\Delta x}\Biggr{]}\,.$$

$$x(q,t):\mathcal{Q}\times \text{mathds}{R}^{+}⟼\mathcal{I}.$$

$$x(0,t)=0$$

$$x(q_j,t^n)=x_j^n.$$

$$\displaystyle J\,(q\,,t)\,c\,(v)\,\dfrac{\partial v}{\partial t}\ =\ \dfrac{\partial}{\partial q}\Bigg{[}\dfrac{k(v)}{J\,(q\,,t)}\,\dfrac{\partial v}{\partial q}\Bigg{]}\ +\ c\,(v)\,\dfrac{\partial v}{\partial q}\,\dfrac{\partial x}{\partial t}\,,$$

$$\displaystyle J_{\,j}^{\,n}\,c_{\,j}^{\,n}\,\dfrac{v_{\,j}^{\,n+1}\ -\ v_{\,j}^{\,n}}{\Delta t}\ =\$$

$$\displaystyle\ +\ c_{\,j}^{\,n}\,\Bigg{(}\dfrac{v_{\,j+1}^{\,n+1}\ -\ v_{\,j-1}^{\,n+1}}{2\,h}\Bigg{)}\cdot\Bigg{(}\dfrac{x_{\,j}^{\,n+1}\ -\ x_{\,j}^{\,n}}{\Delta t}\Bigg{)}\,.$$

$$w(x,t)=1+{\alpha }_1|u{|}^{{\beta }_1}+{\alpha }_2|\frac{\partial u}{\partial x}{|}^{{\beta }_2},$$

$$w_j^n=1+{\alpha }_1|u_j^n{|}^{{\beta }_1}+{\alpha }_2|\frac{u_{j+1}^n-u_j^n}{x_{j+1}^n-x_j^n}{|}^{{\beta }_2}.$$

$$\displaystyle\dfrac{\partial}{\partial q}\biggl{(}w\,(x\,,t)\,\dfrac{\partial x}{\partial q}\biggr{)}\ =\ \beta\,\dfrac{\partial x}{\partial t}\,,$$

$$\displaystyle\dfrac{1}{h}\ \Biggl{(}w^{\,n}_{j+{\frac{1}{2}}}\dfrac{x^{n+1}_{j+1}\ -\ x^{n+1}_{j}}{h}\ -\ w^{\,n}_{j-{\frac{1}{2}}}\dfrac{x^{n+1}_{j+1}\ -\ x^{n+1}_{j}}{h}\Biggr{)}\ =\ \beta\ \dfrac{x^{n+1}_{j}\ -\ x^{n}_{j}}{\Delta t}\,,$$

$$\displaystyle\dfrac{\partial}{\partial q}\biggl{(}w\,(x\,,0)\,\dfrac{\partial x}{\partial q}\biggr{)}\ =\ 0\,,$$

$$\displaystyle\dfrac{1}{h}\ \Biggl{(}w^{\,0}_{\,j+{\frac{1}{2}}}\dfrac{x^{\,0}_{\,j+1}\ -\ x^{\,0}_{\,j}}{h}\ -\ w^{\,0}_{\,j-{\frac{1}{2}}}\dfrac{x^{\,0}_{\,j+1}\ -\ x^{\,0}_{\,j}}{h}\,\Biggr{)}\ =\ 0\,,$$

$$\displaystyle\bar{w}_{\,j\,+\,{\frac{1}{2}}}\ =\ w_{\,j\,+\,{\frac{1}{2}}}\ -\ \sigma\,\bar{w}_{\,j\,+\,{\frac{1}{2}}}\ +\ \dfrac{\sigma}{2}\,\Bigl{(}\bar{w}_{\,j\,-\,{\frac{1}{2}}}\ +\ \bar{w}_{\,j\,+\,\frac{3}{2}}\Bigr{)}\,,$$

$${\bar{w}}_{\frac{1}{2}}=w_{\frac{1}{2}},{\bar{w}}_{N-\frac{1}{2}}=w_{N-\frac{1}{2}}.$$

Full text

Key words and phrases: Optimum insulation thickness; thermal insulation; numerical simulation; moving grid method; redistribution schemes; adaptive numerical methods

MSC:

PACS:

∗ Corresponding author

Last modified:

1. 1 Introduction

2. 2 Mathematical formulation of the physical problem

3. 3 Method of solution

4. 3.1 Quasi-Uniform Nonlinear Transformation

5. 3.1.1 A finite-difference scheme on a fixed uniform mesh

6. 3.1.2 Finite-difference scheme on a moving grid

7. 3.1.3 Construction and motion of the grid

8. 3.2 Validation of the numerical solution

9. 4 Numerical Benchmarking

10. 4.1 Results and discussion

11. 4.2 Convergence study

12. 5 Optimum insulation in Brazil

13. 5.1 First Case: insulated vs. non-insulated buildings

14. 5.2 Second Case: inside insulation vs. outside insulation

15. 5.2.1 Parametric study

16. 5.3 Third Case: wall orientation influence

17. 5.3.1 The roof

18. 5.4 Optimal thermal insulation: economic analysis

19. 6 Conclusion

20. Nomenclature

Introduction

In Brazil, the electricity consumption in residential, commercial and public sectors accounted for $42.8\%$ of the Brazilian electricity consumption in $2016$ [17]. The consumption of electrical energy in buildings is related to gains or losses of heat through the building envelope that is also associated with the internal loads. In particular, the use of heating, ventilation, and air conditioning (HVAC) systems can represent $40\%$ of the total electrical energy consumption [2].

Due to current environmental issues, there is a great potential to reduce heating and cooling loads particularly through improvements of thermal performance of the envelope. One of the most efficient ways to reduce the transmission loads and, consequently, reduce the energy consumption is the use of an appropriated thermal insulation for the building envelope. However, in Brazil, the insulation costs are relatively high, due especially to the low demand and the lack of a compulsory regulation [32]. Moreover, the choice of the insulation thickness is a difficult task since it is a delicate compromise between the cost of insulation materials, which increases with the thickness, and the energy cost for heating/cooling the building, that decreases with the thickness [13, 21].

The determination of the optimum thickness is a difficult task. It requires the computation of the heating or cooling loads of a building along a whole year by solving the heat transfer equation through building walls and roofs. Several studies reported in literature enhance the importance of efficient modeling to assess wall thermal loads. The influence of important wall parameters on the energy efficiency has been investigated: insulation properties in [20], thermal inertia in [7, 39] or thermal conductivity in [26]. In [25], the effect of wall moisture content on the thermal response of the wall has been analysed. Some studies aimed at calibrating the results of the numerical models with observations in order to improve the reliability of the physical model, as for instance in [33, 24, 41, 42].

These studies highlight the importance of efficient numerical models to predict heat transfer through walls. However, the elaboration of such models is a difficult issue since the numerical complexity arises from several points [31, 12]. First, the numerical model has to deal with different time scales since the computation is carried out over one year, as illustrated in [11], while the characteristic time scale of the heat diffusion scales with $15\ \mathsf{min}\,$ or even less. In addition, the boundary conditions are defined according to non-periodic climatic data [36]. Third, the wall configuration involves several layers inducing variations of the thermal properties with space. These difficulties imply an important numerical effort of the numerical model to compute the optimal thermal insulation.

In literature, several works circumvent the numerical difficulties to reduce the computational efforts. For instance, in [12], a whole building model is replaced by lumped approaches, denoted as RC model to reduced the computation efforts. This approach is used in many algorithms to solve the partial differential equation of heat transfer as mentioned in [18, 37, 34]. For the study of optimal thermal insulation, studies also attempt to avoid the numerical difficulties. For instance, in [16], the simple degree-day method is used to compute the building energy consumption in four cities of Turkey. In [21, 13, 8] simulations are performed only on the hottest and coldest days/months of the year. In [35], the outside boundary condition is assumed as periodic by repeating a daily cycle of air temperature and solar radiation. In [4, 3], a representative day for each month of the year is used to compute the yearly loads. In [14], an analytical solution is proposed for simulations carried out only for the hottest and coldest days of the year. For the climate of Paris, [9] proposed a parametric solution of the transient heat transfer using a reduced order model to determine the annual loads and consequently the optimum insulation thickness.

In this way, we propose an innovative method to simulate the transient heat transfer to determine the optimum insulation thickness, which will be tested in Brazilian buildings context. Instead simulating the hottest and coolest months we simulate the whole year to get the transmission loads. The method used for simulations is called Quasi-Uniform Nonlinear Transformation (QUNT) which is adaptive regarding the spatial grid [23]. The grid points of the spatial mesh move following the monitor function which indicates their placement. The adapted grid is obtained according to the gradients of the problem indicated by solving a parabolic problem. As shown in [22], this approach upgrades the usual finite-differences method, with a semi-implicit scheme, to the fourth order by using the non-uniform grid [10]. Just by changing the distribution of nodes one can gain two extra orders of accuracy, which allows simulating with efficiency using less discrete points. This numerical model enables to compute the temperature distribution for the entire year at a low computational effort. From these results, the optimal thermal insulation can be easily computed.

Therefore, this paper first describes in Section 2 the employed physical model of heat transfer. In Section 3, basics of the moving grid method is detailed before exploring the features of the scheme applied to a nonlinear case, presented in Section 4. Then, in Section 5, the method is applied to investigate the optimum insulation thickness in Brazilian buildings. The main conclusions of this study are outlined in Section 6.

Mathematical formulation of the physical problem

The physical problem considers one-dimensional heat transfer through a porous material. The spatial domain is defined by $x\,\in\,[\,0,\,l\,]$, in which $x\ =\ 0$ corresponds to the surface in contact with the inside room and, $x\ =\ l$, corresponds to the outside surface. The heat transfer is governed by conduction, which formulation is derived from Fourier’s law:

[TABLE]

where $x$ is the spatial coordinate $[\mathsf{m}]$, $t$ is the time $[\mathsf{s}]$, $T$ is the temperature $[\mathsf{K}]$111The results are presented in $[\mathsf{{}^{\circ}C}]$ but simulations were performed considering $[\mathsf{K}]$. , $k$, the thermal conductivity $[\mathsf{W/(m\cdot K)}]$ and $c\,=\,\rho\,c_{\,p}\,$, the thermal storage, in which $c_{\,p}$ is the material specific heat $[\mathsf{J/(kg\cdot K)}]$ and $\rho$, the material density $[\mathsf{kg/m^{3}}]$.

The Dirichlet-type boundary conditions are written as:

[TABLE]

Moreover, the internal Robin-type boundary condition is expressed as:

[TABLE]

and on the external side, the boundary condition includes the convection, short- and long-wave radiation:

[TABLE]

where $T_{\,L}$ and $T_{\,R}$ are the temperature of the air $[\mathsf{K}]$ that varies over time, $h$ is the convective heat transfer coefficient $[\mathsf{W/(m^{2}\cdot K)}]$ and subscripts $L$ and $R$ represent the left and right boundaries. If the bounding surface is in contact with the outside building air, $q_{\,\infty}$ is the total solar radiation $[\mathsf{W/m^{2}}]\,$, which includes the direct, diffuse and reflexive radiations and $\alpha$ is the surface absorptivity $[-]$. For the long-wave radiation, $T_{\,\text{sky}}$ is the sky temperature $[\mathsf{K}]$, $\epsilon$ is the emissivity $[-]$ and $\sigma$ is the Stefan–Boltzmann constant $[\mathsf{W/(m^{2}\cdot K^{4})}]\,$.

For the initial condition, a linear temperature distribution in function of $x$ is considered:

[TABLE]

The heat flux density $q\ [\mathsf{W/m^{2}}]$ at the point $x_{\,j}$ is computed as:

[TABLE]

Finally, the integration of the resulting heat flux density during a time period will give the transmission loads $E\ [\mathsf{J/m^{2}}]$ [8, 9]:

[TABLE]

which can be evaluated on daily or monthly periods.

Method of solution

The problem previously described is solved by using the Quasi-Uniform Nonlinear Transformation (QUNT) [23]. However, the problem must be expressed in a dimensionless formulation before being solved numerically. For example, this process can reduce the number of parameters of the problem allowing direct comparisons of models. For a detailed explanation of why this practice is important one can consult [19].

In this way, we define the following dimensionless quantities:

[TABLE]

where the subscript [math] represents a reference value, chosen according to the application problem and the superscript $\star$ represents a dimensionless quantity of the same variable. Therefore, the governing Equation (2.1) can be written in a dimensionless form as:

[TABLE]

with the following dimensionless formulation of the Dirichlet and Robin-type boundary conditions:

[TABLE]

For the initial condition the dimensionless form is:

[TABLE]

Moreover, the dimensionless form of the heat flux density is written as:

[TABLE]

In the following, we drop $\star$ for the sake of notation compactness but without losing generality.

Quasi-Uniform Nonlinear Transformation

The Quasi-Uniform Nonlinear Transformation (QUNT) [23, 22, 10] is a numerical method with an adaptive spatial grid which provides a high resolution where it is required. The adaptivity is obtained by using a monitor function which indicate the location of the spatial grid points during the simulation. The adapted grid is obtained by means of the equidistribution principle. This moving grid nodes approach has the advantage of being conservative in space and also second-order accurate.

Before solving the problem on a moving mesh, one have to chose a robust and accurate scheme on a fixed grid. This scheme will be then generalized to incorporate the motion of the mesh. This scheme is described in the sequence.

3.1.1  A finite-difference scheme on a fixed uniform mesh

Consider the dimensionless form of the heat diffusion Equation (3.1) with a simplified notation, just for describing the numerical method:

[TABLE]

The spatial domain $\mathcal{I}\,=\,[\,0\,,l\,]$ is uniformly discretized, with $x_{\,j}\,=\,j\Delta x,\,j\in\mathbb{Z}$ being the spatial nodes, in which $\Delta x\,>\,0$ is the spatial step size. For the time domain $\mathcal{P}\,=\,[\,0,\,\tau\,]$, also consider a uniform discretization, where $t^{\,n}\,=\,n\,\Delta t\,,n\in\mathbb{Z}^{+}$ are the temporal nodes and $\Delta t\,>\,0$ being the time step. Thus, by using an IMplicit-EXplicit (IMEX) scheme, Equation (3.2) is written in a discrete form:

[TABLE]

This scheme approximates the continuous operator to order $\mathcal{O}\,\,(\Delta x^{\,2}\ +\ \Delta t)$. The advantage of this semi-implicit scheme over the fully implicit is to avoid sub-iterations in the solution procedure and, at the same time, being stable and consistent. For simplicity, the boundary conditions are not treated here. In what follows, the method based on this scheme with the moving grid is describe.

3.1.2  Finite-difference scheme on a moving grid

Consider the computational domain $\mathcal{I}\,=\,[\,0\,,l\,]\,$, the reference domain $\mathcal{Q}\,=\,[\,0\,,1\,]$ and a bijective time-depend mapping from $\mathcal{Q}$ to $\mathcal{I}\,$:

[TABLE]

This represents the quasi-uniform nonlinear transformation. In addition, the boundary points are required to map into each other:

[TABLE]

The reference domain $\mathcal{Q}$ is uniformly discretized into $N$ elements, with $q_{\,j}\,=\,j\,h,\,j\in\mathbb{Z}^{+}$ representing the nodes of the grid, which is equally spaced $h\,=\,1/N\,$. However, only the image of nodes $q_{\,j}$ under the map $x\,(q\,,t)$ are needed, since they constitute the nodes of the moving mesh:

[TABLE]

The diffusion Equation (3.2) is rewritten on the domain $\mathcal{Q}$, with the help of the composed function v\,(q\,,t)\mathop{\stackrel{{\scriptstyle\ \mathrm{def}}}{{:=}}\ }(u\circ x)\,(q\,,t)\,\equiv\,u\,\Bigl{(}x\,(q\,,t)\,,t\Bigr{)}\,:

[TABLE]

where $J\mathop{\stackrel{{\scriptstyle\ \mathrm{def}}}{{:=}}\ }\dfrac{\partial x}{\partial q}$ is the Jacobian of the transformation $x\,(q\,,t)\,$.

Therefore, the discrete form of the diffusion heat equation on a moving mesh grid is:

[TABLE]

In summary, this is the parametrization process of the mesh motion by a bijective mapping. In the sequence, the construction of this map is explained.

3.1.3  Construction and motion of the grid

The equidistribution method is used to construct the adaptive grid. To control the distribution of the nodes, a positive valued function $w\,(x,\,t)$ has to be chosen. This function is called monitor function. The choice of the monitor function and its parameters are important for the accuracy of the scheme. In particular, for unsteady nonlinear problems, there are several choice options. One can find more information, for example in [38]. In this work, the function chosen is based on [23] and can be written as:

[TABLE]

where $\alpha_{\,1}\,$, $\alpha_{\,2}\,$, $\beta_{\,1}$ and $\beta_{\,2}$ are positive real parameters. The quantity $w\,(x\,,t)$ assumes large positive values in areas where the solution has important gradients, inducing the elements of the grid mesh to approach. The discrete formulation of the monitor function presented in Eq. (3.3) is written as:

[TABLE]

A non-uniform grid $\mathcal{I}$ is given if we construct the mapping $x\,(q\,,t)\,:\ \mathcal{Q}\ \longmapsto\ \mathcal{I}$ and evaluate it in the nodes of the uniform grid. To this end, the equidistribution method proposes that the desired mapping $x\,(q\,,t)$ be obtained as a solution of a nonlinear parabolic problem, which is written as:

[TABLE]

with $x\,(0,\,t)\ =\ 0$ and $x\,(1,\,t)\ =\ l$ as the boundary conditions. Eq. (3.4) can be written in the discrete form as:

[TABLE]

with the same boundary conditions $x_{\,0}^{\,n\,+\,1}\,=\,0$ and $x_{\,N}^{\,n\,+\,1}\ =\ l\,$. The parameter $\beta\ >\ 0$ plays the role of the inverse diffusion coefficient and it controls the smoothness of nodes trajectory.

Initial grid generation.

As the problem depends also on time, it is of capital importance to obtain a high-quality initial mesh. The initial condition must be adapted to the new grid before starting the dynamical simulation. At $t\,=\,0$ we compute the monitor function $w\,(x,\,0)$ of the initial condition $u\,(x,\,0)\,$. Then, the mapping $x\,(q,\,0)$ is determined as the solution of Eq. (3.4), which is reduced to a simple second-order ordinary differential equation:

[TABLE]

supplemented with the Dirichlet boundary condition. Thus, the finite-difference approximation of this latter differential equation is:

[TABLE]

with the discrete boundary conditions $x_{\,0}^{\,0}\,=\,0\,$, $x^{\,0}_{\,N}\,=\,l\,$. This latter nonlinear system of equations is solved iteratively, with iterations initialized with a uniform grid as the first guess. Its solution satisfies the equidistribution principle: in areas where $w^{\,0}_{\,j\,+\,{\frac{1}{2}}}$ takes large values, the space between two neighboring nodes $x^{\,0}_{\,j+1}$ and $x^{\,0}_{\,j}$ has to be inversely proportionally small. One can find more information in [22].

Smoothing step.

To ensure the smoothness of the mesh motion the monitor function $w$ can be filtered. This process, called smoothing step enable to enlarge the set of acceptable values of the parameters $\alpha$ and $\beta$ of the monitor function. According to [23], the following implicit scheme can produce robust results:

[TABLE]

where $\sigma$ is a positive smoothing parameter. To complete Eq. (3.5), boundary conditions are taken as in the original problem:

[TABLE]

The smoothed discrete monitor function \Bigl{\{}\bar{w}_{\,j\,+\,{\frac{1}{2}}}\Bigr{\}}_{\,j\,=\,0}^{\,N\,-\,1} is then obtained by solving the linear System (3.5). Note that the smoothing operator is applied to the monitor function.

On the choice of the parameters.

The parameters $\alpha_{\,1}\,$, $\alpha_{\,2}\,$, $\beta_{\,1}\,$, $\beta_{\,2}\,$, $\beta$ and $\sigma$ are chosen according to the problem and the solution under consideration. To determine the optimal values one has to make some numerical experiments. Interested readers are invited to consult [23] for more details on this aspect.

Validation of the numerical solution

All the numerical results in this paper were computed using Matlab™ [28]. The reference solution is computed using the open source package Chebfun [15]. Using pde15s function, it enables to compute a numerical solution of a partial derivative equation with the Chebyshev polynomials adaptive spectral methods.

The Euclidean error (distance) for a determined field $u$ at time $t^{\,n}$ is computed as:

[TABLE]

where $u^{\,\text{ref}}$ is the reference solution computed with Chebfun and $u^{\,\text{num}}$ is the solution computed with the QUNT or classical schemes.

The error over the time is computed as a function of $x$ by the following formulation:

[TABLE]

where $N_{\,t}$ is the number of temporal steps. The global $\mathcal{L}_{\,\infty}$ error is given by the maximum value of $\varepsilon\,(x)\,$:

[TABLE]

For the heat flux density, the Euclidean error (distance) is computed as:

[TABLE]

where $q^{\,\text{ref}}$ is the reference flux density computed with Chebfun solution and $q^{\,\text{num}}$ is the flux density computed with the QUNT or classical solutions. The global $\mathcal{L}_{\,\infty}$ error of the flux density is also given by the maximum value but of the $\xi\,(t)\,$:

[TABLE]

Numerical Benchmarking

The purpose of this section is to validate the QUNT method by comparing it with the IMEX scheme and the reference solution. The problem is inspired by [19]. As the purpose is to evaluate the efficiency of the QUNT approach, and not analyzing the physical phenomena, the study is performed using the dimensionless problem described by Equation (3.1) with Dirichilet-type boundary conditions and the following dimensionless material properties:

[TABLE]

The dimensionless ambient temperature at the left and right boundaries are:

[TABLE]

which are graphically given in Figure 1. The initial temperature is considered uniform over the material and its dimensionless formulation is given by:

[TABLE]

The monitoring function is defined by Equation (3.3) with following values for its parameters:

[TABLE]

These parameters were determined by our numerical investigations and correspond to the best combination for this case study. For the grid diffusion parameter, the best value tested was $\beta\ =\ 100$ and for the smoothing parameter $\sigma\ =\ 10\,$. Simulations are performed up to the final time of $\tau\ =\ 48\,$, with $\Delta t^{\,\star}\ =\ 5\cdot 10^{\,-\,3}$ and a number of spatial points $N_{\,x}\ =\ 51\,$, at first fixed for both methods.

Results and discussion

The sample trajectory of the spatial nodes can be observed in Figure 2. The spatial points concentrate where the values of the gradient are high, particularly on the boundaries where important variations of the fields are imposed.

Figure 3(a) presents the temperature profile for $t^{\,\star}\,=\,16$ of the computed solutions. As the material has a low thermal conductivity, the temperature in the interior of the material does not change as fast as it does at near the boundaries. Consequently, high gradients can be observed on the boundaries, making the spatial nodes concentrate in these areas. The QUNT method handles much better these high gradients than the IMEX, providing more accurate results. As the QUNT method was built based on the IMEX scheme, we can conclude that the moving grid modifications improved the numerical scheme. The absolute error $\delta$ for this instant of time is presented in Figure 3(b). This figure shows that the QUNT has a spatial accuracy the order of $\mathcal{O}\,\,(10^{\,-\,4})\,$, while the IMEX scheme has an accuracy the order of $\mathcal{O}\,\,(10^{\,-\,3})\,$. Thus for the same quantity of spatial nodes $(N_{\,x}\,=\,51)\,$, the QUNT can be $10$ times more accurate than the numerical solution of the IMEX scheme. If we look at the boundaries, this ratio is even higher reaching $\approx\,100$ times.

The results for the error over the time are presented in Figure 4. The order of the error remains stable for the whole simulation, with a maximum error of $\varepsilon_{\,\infty}\,=\,4.78\cdot 10^{\,-\,4}$ for the QUNT and $\varepsilon_{\,\infty}\,=\,2.82\cdot 10^{\,-\,2}$ for the IMEX method, which confirms $\approx\,100$ ratio between two approaches.

After simulating the temperature field, it is of great importance to get the values of the heat flux density at the boundaries, to compute, for instance, the transmission loads. Therefore, in what follows, we study the influence of the number of spatial nodes on the heat flux density.

Figures 5(a) and 5(b) present the heat flux evaluated at the left and right boundaries, respectively, with the derivative approximated with a finite-difference scheme second-order accurate in space. The results show that there is a difference between the flux computed for the uniform grid with the one computed for the non-uniform grid.

Figures 6(a) and 6(b) represent the error of the fluxes $\xi\,$, computed with the QUNT and IMEX solutions for the left and right boundaries, respectively. In general, the flux computed with the non-uniform grid is more accurate than the flux computed with the uniform grid. By using $51$ spatial nodes, the maximum error of the flux computed with an uniform grid is $\xi_{\,\infty}\,=\,1.9\cdot 10^{\,-\,1}$ and the flux computed with a non-uniform grid is $\xi_{\,\infty}\,=\,8.5\cdot 10^{\,-\,2}\,$.

Convergence study

A convergence study was carried out to compare the QUNT with the IMEX scheme and to analyze the accuracy of both temperature distribution field and heat flux density. Thus, for different values of \Delta t^{\,\star}\,=\,\bigl{\{}\,10^{\,-\,1}\,;10^{\,-\,2}\,;10^{\,-\,3}\,;10^{\,-\,4}\,\bigr{\}}\,, the global error $\varepsilon_{\,\infty}$ was computed as a function of the number of spatial nodes $N_{\,x}\ \in\ [\,10\,,80\,]$. Figures 7(a), 7(b), 7(c), and 7(d) display the results.

It can be seen in Figure 7 that for $N_{x}\ <\ 20$ both methods provide solutions of the field with approximately the same order of accuracy, which is insufficient to get good approximations of the solution. For the IMEX method, less than $100$ nodes are not enough to give a solution with an accuracy better than $\mathcal{O}\,\,(10^{\,-\,2})\,$. Even if the $\Delta t^{\,\star}$ is decreased, the global error is still of the same order of accuracy due to insufficient spatial discretization. On the other hand, for the QUNT method, as the number of spatial nodes increase, the global error decreases faster than the IMEX method. For the QUNT adaptive grid, with only $40$ nodes, an accuracy of $\mathcal{O}\,\,(10^{\,-\,3})$ can be achieved for the three values of $\Delta t^{\,\star}\,$.

Figures 8(a), 8(b), 8(c) and 8(d) display the results for the global error of the flux density $\xi_{\,\infty}\,$, as a function of the number of spatial nodes $N_{\,x}\ \in\ [\,20\,,100\,]$ also for different values of $\Delta t^{\,\star}\,=\,\{\,10^{\,-1}\,;10^{\,-2}\,;10^{\,-3}\,;10^{\,-4}\,\}$. Results show that the computation of the flux density is less accurate than the computation of the temperature field. This happens since the error increases when the derivative has to be computed. The QUNT method also provides solutions for the flux more accurate than the IMEX method. As the analysis of the error for the field, the QUNT approach proved to be more sensible regarding variations of $\Delta t$, as the solution becomes smoother and the moving grid adapts better when $\Delta t$ decreases. Although, the same accuracy is observed for $\Delta t\,=\,10^{\,-\,3}$ and $\Delta t\,=\,10^{\,-\,4}$ because the order of the error is also related with the spatial discretization. To obtain higher accuracy, it is necessary to increase the derivative approximation order. In this case, the error of the flux is relatively high due to the nonlinearity of the case study. For a less nonlinear case, the error of the flux density can be considerably reduced.

Therefore, it has been clearly noticed the QUNT approach computes an accurate solution of the problem. Moreover, the moving grid method converges to an accurate solution with less spatial nodes than the method with an uniform grid. In a long-term simulation, it does provide significant gains in regard to the CPU usage. The computer used to assess the computational cost has a processor Intel®CoreTMi$5$ @ $2.80$GHz$\times 4$ with $15.6$GB of RAM.

To evaluate the gains in terms of computational cost, Table 1 provides the CPU time to compute the solution using the Crank–Nicolson, the QUNT and the IMEX methods, which have been evaluated using the Matlab™ platform. All methods computed a solution with the same order of accuracy $\varepsilon_{\,\infty}\cong\mathcal{O}\,\,(10^{\,-\,4})\,$. For this, it was used a time step of $\Delta t^{\,\star}\,=\,10^{\,-\,2}\,$, and $51$ spatial nodes for the QUNT and $1\,001$ and $501$ spatial nodes for the Crank–Nicolson and IMEX methods, respectively. The computational effort to perform the simulation increases linearly with the simulation time. The QUNT and IMEX methods can compute the solution $25\%$ faster than the Crank–Nicolson approach. It is preferable to focus on the ratio of computer run-time rather than on absolute values, since it is system-dependent. The Crank–Nicolson scheme was implemented with a fixed point approach to treat the non-linearity and it was used here only for comparison purposes. One should notice that the codes are not optimized, and the differences between the CPU time of QUNT and IMEX methods are low, even with the QUNT method performing two extra matrix inversions.

Optimum insulation in Brazil

In the previous section, the QUNT method was presented and validated. Now, the method is used to compute the solution of a transient one-dimensional case study of heat transfer through building walls for a period of one year. The results are obtained with less computational efforts and enable to estimate the optimum thermal insulation for large cities in Brazil.

Climatic zones.

Brazil has a total surface area of $8.5$ million $\mathsf{km^{2}}$, the fifth largest country in the world, occupying almost half of the entire South America continent. In such a large country a variety of climates exist. However, most of the country can be defined as being tropical and sub-tropical and energy consumption for space conditioning is dominated by cooling. Although, in the Southern region of Brazil, energy is also used for heating during the winter months (June – August).

For this study, four Brazilian cities — Curitiba, Rio de Janeiro, São Paulo and Salvador — have been selected for determining their optimal insulation thickness, chosen due to the diversity of the climate and due to their large population [32]. However, the method used in this study can be applied to any city desired. Based on ASHRAE Standard $90.1-2013$ [6], these cities are located in three different climatic zones: Very Hot-Humid (zone $1$), Hot-Humid (zone $2$) and Warm-Humid (zone $3$). Table 2 presents some characteristics of the chosen cities.

Simulation parameters.

Bricks are commonly used in typical wall structures in Brazilian constructions [32, 1]. As illustrated in Figure 9, we consider $3$ types of walls: Wall$-1$ is composed only of brick; Wall$-2$ includes also a layer of insulation on the internal side; Wall$-3$ has a layer of insulation on the external side. The brick is $15\ \mathsf{cm}$ thick, while the thickness of insulation material (extruded polystyrene) varies. The thermal properties of the materials are gathered from the IEA Annex$-24$ [27] and are reported in Table 3.

As the initial condition, a steady-state temperature distribution within the wall is assumed:

[TABLE]

Climatic conditions of temperature and solar radiation are considered on the outside part of the wall for a period of one year. The outside temperature and total radiation are given in Figures 10 and 11. The outside heat flux is computed using radiation data from Domus [30] for a South-oriented wall. The solar absorptivity of the outside surface is taken as $\alpha\,=\,0.5\,$. The indoor air temperature $T_{\,R}(t)$ has a sinusoidal variation as illustrated in Figure 10. During the summer period, the temperature is set to $25\mathsf{{}^{\circ}C}$ and during the winter to $20\mathsf{{}^{\circ}C}$ [29]. The convective heat transfer coefficients at the outer and the inner surfaces are set to constant values $h_{\,L}\,=\,25\ \mathsf{W/(m^{\,2}\cdot K)}$ and $h_{\,R}\,=\,10\ \mathsf{W/(m^{\,2}\cdot K)}$, respectively. The long-wave radiative transfer on the external surface is only considered for the roof simulation.

According to Figure 10, the inside temperature in all cities varies in the same way. What can be observed, is that Curitiba is the coldest city, which means it needs more heating than cooling. In Rio de Janeiro and Salvador, the external temperatures are higher than the inside temperature, which means that the cooling system predominates. In addition, Salvador has a higher mean external temperature than Rio de Janeiro. São Paulo has the mildest weather among the four cities, but the external temperatures are mostly lower than the internal ones, which makes the heating need prevails.

In next sections, the thermal behavior of buildings located in these cities is studied for different wall configurations, various thicknesses of insulation material and different wall orientations, with respect to the heat flux density $q$ and transmission loads $E$ on the inside surface of the building wall. Three cases are distinguished. First, the benefits of insulating buildings are demonstrated. The second case analyzes if the insulation should be placed on the inside or outside part of the building. Last, the effect of the wall orientation is studied.

First Case: insulated vs. non-insulated buildings

For the first simulation, two South-oriented walls in the city of Curitiba are considered. One has no insulation (Wall$-1$) and the other one has $10\ \mathsf{cm}$ insulation (Wall$-2$). Climatic conditions of the temperature and solar radiation are given in Figures 10(a) and 11.

With the numerical values provided previously, Problem (2.1) is solved using the QUNT method. The monitoring function is defined by Equation (3.3), which has the following values for its parameters: $\alpha_{1}\,=\,0.8\,$, $\beta_{\,1}\,=\,2\,$, $\alpha_{\,2}\,=\,0.2\,$, $\beta_{\,2}\,=\,3\,$. The grid diffusion and the smoothing parameters are, respectively, $\beta\,=\,50$ and $\sigma\,=\,5\,$. The solution is obtained for a number of $N_{\,x}\,=\,41$ spatial nodes with a time discretization of $\Delta t^{\,\star}\,=\,0.1\,$, which is equivalent to $\Delta t\,=\,6\,\mathsf{min}$.

The temperature is computed for all points in the single- and multi-layered wall configurations. Then, the heat flux on the inside surface of the building wall is estimated and presented in Figure 12(a). For the case with no insulation, the heat flux has an important variation, alternating between positive and negative values. With $10\ \mathsf{cm}$ insulation, it strongly modifies the thermal behavior of the wall and the heat flux is considerably reduced in this case.

Yearly transmission loads are calculated based on Equation (2.4). Hourly variation of inside surface heat flux density is integrated over a $24\ \mathsf{h}$ period to obtain a daily total load. Then, daily total loads are integrated to get the monthly heating/cooling transmission loads presented in Figure 12(b), for Wall$-1$ and Wall$-2\,$. The heating/cooling loads are ten times higher for the wall without insulation indicating that much more energy is spent in order to maintain the temperature inside of the building.

With these results, the first conclusion is that insulation may play an important role in the construction sector in Brazil. In what follows, two configurations are tested, one with the insulation inside of the building wall (Wall$-2$) and another one with the insulation outside (Wall$-3$).

Second Case: inside insulation vs. outside insulation

For this second case, it is also considered a South-oriented wall in the city of Curitiba, with $10\ \mathsf{cm}$ insulation. Two wall configurations are verified: one with the insulation inside the building (Wall$-2$) and another one with the insulation outside (Wall$-3$). The same climatic conditions and numerical values of the previous case (Section 5.1) are taken.

Figure 13(a) shows the effect of the insulation position on the heat flux density during one year of simulation. The wall with external insulation has the lowest heat flux values but with a small discrepancy, with the mean difference between the heat flux of Wall-$2$ and Wall-$3$ being of the order of $0.1\,\mathsf{W/m^{2}}\,$, for a $10\ \mathsf{cm}$ insulation. This difference can also be noted in the monthly heating/cooling transmission loads as shown in Figure 13(b). Depending on the month, the value of the energy load can be higher for the insulation outside of the building wall. However, the annual heating/cooling transmission loads are very close, with $33.83\,\mathsf{MJ/m^{2}}$ for Wall$-2$ and $34.81\,\mathsf{MJ/m^{2}}$ for Wall$-3\,$. In [8] they also have similar results for the annual loads regarding the position of the insulation (inside/outside), even with different assumptions and climatic conditions.

5.2.1  Parametric study

To improve this analysis, a parametric study is performed, in which the thickness of the insulation material varies on the interval $l_{\,i}\,=\,[0.01,\,0.3]\,\mathsf{m}\,$, for both wall configurations (Wall$-2$ and Wall$-3$).

Simulations are performed with the same numerical and physical configuration of the previous case studies. The annual heating/cooling transmission loads are calculated and presented as a function of $l_{\,i}$ in Figure 14. As expected the thicker the insulation, the lower the energy consumption for either heating or cooling. One can easily notice that transmission loads drop fast for small insulation thickness values. Besides, for thin insulations, the difference between the two configurations is more evident. Although, this difference is minimized for high values of thicknesses.

With those results, another conclusion from this case study is that the difference of configurations with the insulation inside or outside does not cause important changes on annual heating/cooling when the thickness is higher than $2\ \mathsf{cm}$, despite the thermal mass on the internal side provides important effects on a shorter timescale. For this reason, just for simulation purposes, in what follows, the configuration of Wall$-2$ is used to study the impact of wall-orientation.

Third Case: wall orientation influence

For the third case study, the same features of the previous case are considered but changing only the facing wall orientation (North, South, East, and West). As Brazil is located in the Southern hemisphere, the North façade wall of a building is the one that receives more short-wave radiation.

Figure 15 shows the effect of wall orientation on the monthly heating/cooling transmission loads of Wall$-2$ model with a $10\ \mathsf{cm}$ thick insulation. The monthly transmission loads values are mostly negative indicating a heating demand, which is compatible with the weather of Curitiba. In colder months, the South façade has the highest heating transmission loads because of the low solar heat gain, which is the opposite for the North-facing wall.

The effect of wall orientation on yearly heating and cooling transmission loads are presented in function of the insulation thickness in Figure 16. The South-facing wall has the highest transmission loads, then, East-, North- and West-facing walls with very similar values. As previously observed, transmission loads quickly drop is faster for small values of insulation thickness.

5.3.1  The roof

The roof is one of the building components whereby occurs the most significant heat gain/loss. An identical analysis of the wall structure is performed for a roofing system. For this analysis, a flat roof is taken into account and we consider configurations identical to the wall structure, with internal and external insulation. The roofing is made of a $15\ \mathsf{cm}$ concrete fixed layer while the thickness of the insulation material varies. The thermal properties of the materials are reported in Table 3.

Simulations were performed by considering radiative and convective heat transfer at the external roof surface, which exchanges heat with the air and the sky, and receives solar radiation. The solar absorptivity of the outside surface of the roof is taken as $0.5\,$. For the long-wave radiation exchange, a thermal emissivity of $0.9$ is considered and the sky temperature is provided by the software Domus [30]. Instead of computing the transmissions loads only for Curitiba, the four cities were comprised in this study.

Figure 17 presents the annual heating/cooling transmission loads as a function of insulation thickness, for the four Brazilian cities. Salvador and Rio de Janeiro are the cities that have the highest transmission loads, which means that they will need a thicker insulation. The differences for the transmission loads are higher for values of insulation thickness. It can also be noticed that the differences are reduced, while the insulation thickness increases.

Optimal thermal insulation: economic analysis

To assess the economic viability of the insulation, a simple economic analysis based on [40, 9] is performed. The total cost of the insulation $C_{\,T}$ is estimated as:

[TABLE]

where $C_{\,E}$ is the annual cost of the electric energy consumed and $C_{\,I}$ is the the cost of the insulation. The energy cost is calculated by using the annual conduction load $E\,$, the energy efficiency $\eta$ and the energy price $p_{\,e}\,$:

[TABLE]

The insulation cost $p_{\,i}$ comprises the price of the material as well as its installation and it varies according to the thickness $l_{\,i}\,$:

[TABLE]

The optimal insulation thickness is given when $C_{\,T}$ has minimum owning and operating costs. For this analysis, the parameters used are given in Table 4. The energy price is the average price of the electric power in Brazil for $2017$ [5].

Figure 18 presents the insulation cost, the energy consumption cost, and the total cost, as a function of the insulation thickness. As the insulation thickness increases, the energy cost decreases because less energy is used to maintain the temperature inside of the building. However, installation cost increases linearly with the insulation thickness. Thus, the compromise between the energy cost and the insulation cost is the minimum value of their sum. For the case study of Section 5.2.1, it corresponds to an optimum insulation thickness of $0.5\ \mathsf{m}$ for both wall configurations (Wall$-2$ and Wall$-3$).

This economy study is extended to other three cities, considering the four oriented façades and with the insulation on the inside part of the building wall. Table 5 summarizes the values of optimum insulation thickness for East, West, South, and North oriented walls for the four Brazilian cities. In the case of the wall, the thickness of the insulation is the same for the external and internal configurations. Results show that Salvador is the city that needs the largest insulation thickness, as the inside temperature has the largest difference regarding the outside temperature. Salvador is the hottest city among them with cooling demand along the whole year.

The Northern and Western façades of Rio de Janeiro and São Paulo are the ones receiving the most solar radiation. In Rio de Janeiro, due to the cooling demand, these façades require larger insulation thickness, as reported in Table 5. However, in São Paulo the heating need prevails, and both façades orientations have smaller thickness values since their transmission loads are lower. On the other hand, in Curitiba and Salvador, the heat gains are higher for the North, West and East façades, making the South-facing wall require thicker insulation in Curitiba and thinner in Salvador.

For the roof, the position of the insulation generate differences in the optimum thickness. As shown in Table 6, Curitiba and São Paulo have different values of the optimal thickness of external and internal roof insulation. For Rio de Janeiro and Salvador, the thickness of the insulation must be higher, as they exchange more heat through the roof. For all cities, the roof needs more insulation than the wall, as previously mentioned.

Conclusion

In this work, the Quasi-Uniform Nonlinear Transformation (QUNT) method has been presented as an innovative method to perform the transient heat transfer, which was applied to an energy optimization problem. The method is based on a non-uniform adaptive grid technique that identifies were the spatial nodes must be placed. It was first compared with a reference solution to verify its features. As a result, a satisfactory accuracy is shown by using much less spatial nodes than traditional methods, by moving its nodes where the temperature gradients are higher i.e., an average ratio of $100$ compared with a uniform grid approach. Another remarkable point of this method is that it is easy to be implemented, by adding extra differential equations without being computational costly. An expressive gain of $25\%$ is obtained.

The high accuracy of the solution instigated the analysis of the physical behavior. The optimum insulation thickness of buildings walls is determined by taking into account the wall orientation and the position of the insulation in Brazilian buildings, which is determined through a parametric study. The annual heating/cooling transmission loads have been calculated using weather data and transient heat flow, simulated over one whole year.

First, a case was focused on comparing a wall with insulation and another one with no insulation, applied to the climatic conditions of one Brazilian city. This case proved the importance of the insulation in reducing the transmission loads. Then, a second case study is simulated to verify if the insulation is better inside or outside of the building wall. The last case considered analyzing the effect of wall orientation and the roof assembly. As Brazil is mostly located in the Southern Hemisphere, the North façade is the one receiving more solar radiation. In consequence, for the hottest climate cities, the North-facing wall is the one that requires the thicker insulation layers, while, in the South region – Curitiba –, the North-facing wall requires thinner layers of insulation.

A basic economic analysis was presented, demonstrating that the implementation of an insulation layer appears to be a cost-effective measure to save energy consumption. The optimum insulation thickness for all the wall and roof configurations studied is to be between $2\ \mathsf{cm}$ and $7\ \mathsf{cm}$. The city that requires less insulation is São Paulo, due to its mild climate, and the city that needs more insulation is Salvador, due to its temperature to be high above the thermal comfort values.

As the proposed numerical method looked very promising, further work is recommended to integrate it to whole-building simulation tools so that more realistic cases, including internal gains and complex boundary conditions, can be taken into account in order to provide more conclusive results and enhance building energy simulation tools.

Acknowledgements

This study was financed in part by the Coordenação de Aperfeiçoamento de Pessoal de Nível Superior – Brasil (CAPES) – in the framework of the International Cooperation Program CAPES/COFECUB (Grant $\#774/13$). The Authors also acknowledge the support from CNRS/INSIS and Cellule Énergie under the grant MN4BAT–2017. Professor Mendes thanks the Laboratory LAMA UMR 5127 for the warm hospitality during his visits in 2018, which were supported by the project MN4BAT under the AAP Recherche 2018 programme of the University Savoie Mont Blanc. Finally, the authors acknowledge the Junior Chair Research program “Building performance assessment, evaluation and enhancement” from the University Savoie Mont Blanc in collaboration with The French Atomic and Alternative Energy Center (CEA) and Scientific and Technical Center for Buildings (CSTB).

Nomenclature