## Abstract

The concept of both penetration and deviation times for rectangular coordinates along with the principle of superposition for linear problems allows short-time solutions to be constructed for a one-dimensional (1D) rectangular finite body from the well-known solutions of a semi-infinite medium. Some adequate physical considerations due to thermal symmetries with respect to the middle plane of a slab to simulate homogeneous boundary conditions of the first and second kinds are also needed. These solutions can be applied at the level of accuracy desired (one part in 10^{A}, with A = 2, 3, …, 15) with respect to the maximum temperature variation (that always occurs at the active surface and at the time of evaluation) in place of the exact analytical solution to the problem of interest consisting of an infinite number of terms and, hence, unapplicable.

## 1 Introduction

The symbolic expression of an “*exact*” analytical solution to a transient, 1D heat conduction problem of Fourier type for a finite body, denoted by $Te(x,t)$, with $x\u2208[0,L]$ and $t\u22650$, is available only in an infinite series form, that is, consisting of an infinite number of terms. Though the solution to the above initial-boundary value problem is unique [1 (Chapter 1)], the symbolic expression of the related exact analytical solution is not unique. In fact, there exist two different symbolic forms of it in the specialized heat conduction literature [2–4]: (1) the short-time form coming from Laplace transform combined to an asymptotic series expansion and (2) the large-time form deriving from the application of the separation-of-variables method. Both short- and large-time expressions of the solution exhibit single summations consisting of an infinite number of terms.

As no calculator can account for an infinite number of terms, the exact analytical solution needs to be replaced with an approximate one consisting of a finite number of terms. Though this is well-known, only a few existing references address the required number of terms in the series evaluation for small and large values of the time. It seems that the mathematical expression is all the reader needs, when in fact, for engineering applications, considerable additional work is usually required to obtain reliable and accurate numerical values.

### 1.1 Literature Review.

In an earlier paper by Beck et al. [5], the possibility of using small-time approximate solutions at the heated surface of some basic one-dimensional geometries was investigated and discussed. The treatment was, however, limited to the computation of only surface temperature for only constant-heat-flux boundary conditions. Also, for a plate, the approximate solution was consisting of two terms and valid for dimensionless times $t\u02dc=\alpha t/L2$ less than 0.2. But its accuracy in this time range was not given. A related paper restricted to 1D simple geometries was presented by Lavine and Bergman [6]. The small-time approximations derived by the authors using the Euler-Maclaurin summation formula for only surface temperature coincide with the 1D semi-infinite solutions for plane walls and, hence, it consists of only one term. It is suggested to use it for dimensionless times less than 0.2 but no accuracy is given, and only constant surface heat flux condition was investigated.

Another related paper was presented later by de-Monte et al. [7] for only rectangular bodies. But the proposed small-time approximate solutions are given for both temperature and heat flux and are valid for any location of the body and for any kind of boundary condition and its time variation. In addition, their time range of applicability is given as a function of location of evaluation and accuracy desired there. In detail, for a 1D finite rectangular body, temperature and heat flux can be assumed to be zero at a location $x$ with errors less than 10* ^{-A}*, being

*A =*2, 3, …, 15, if there the time is less than the so-called penetration time $tp=x2/(10A\alpha )$. It is defined as the time that it takes for a point of a body to be just affected by the “active” heated/cooled boundary in terms of temperature and heat flux. By “just affected” it means at the level of one part in 10

*with respect to the maximum temperature rise/fall that occurs at $x=0$ at the same time, i.e., $\Delta Tmax=|T(0,tp)\u2212Tin|$, where $Tin$ is the uniform initial temperature of the plate. In addition, temperature and heat flux can be computed at a location $x$ of the plate using the solutions of the related semi-infinite problem with errors less than 10*

^{A}*, if there the time is less than the so-called (first) deviation time $td(1)=(2L\u2212x)2/(10A\alpha )\u2265tp$. It is defined as the time that it takes for a point of a plane wall heated/cooled at $x=0$ to be affected (in terms of temperature and heat flux) by the presence of the homogeneous (inactive) boundary at $x=L$ at the level of one part in 10*

^{-A}*with respect to $\Delta Tmax$.*

^{A}A further paper on the same topic was proposed by McMasters et al. [8]. It extends the concept of penetration time given in Ref. [7] to cylindrical and spherical geometries. In addition, it is proven that the penetration time corresponding to a constant-temperature boundary condition, defined through $tp=x2/(10A\alpha )$ for rectangular 1D bodies [7], is the most conservative among the $T0(t/t0)n/2$ time-varying boundary conditions, with *n = *0, 1,2, …, and $t0$ arbitrary reference time. By extension, the first kind boundary condition can also be considered the most conservative case for spherical and cylindrical geometries. Later Ostrogorsky and Mikic [9] extend the semi-infinite solutions to radial systems such as cylinders and spheres for short values of the time. Boundary conditions of first, second and third kind are considered. The time ranges of validity of the small-time approximations consisting of only one term and applicable at any location of the body are given. Also, their accuracy is provided in a graphical form through a relative error.

### 1.2 Motivation and Outline.

The objective of this paper is to construct very simple analytical solutions for 1D finite rectangular bodies applicable at any location and at the level of accuracy desired (i.e., one part in $10A$ with respect to $\Delta Tmax$ defined before) for small values of the time but larger than the (first) deviation time $td(1)$ defined in the previous subsection. The extension of the time range of applicability of the solution requires increasing its number of terms consisting currently of only one single semi-infinite term for times less than or equal to $td(1)$ [7]. For this purpose, a new short-time solution will be constructed superimposing two appropriate semi-infinite terms and, hence, will be valid for longer times till $t=td(2)>td(1)$. An algebraic conservative expression is also derived for $td(2)$ (termed as “second deviation time”) depending on both the location of evaluation and accuracy desired.

Therefore, the time interval extension stated before is the prime objective of the current work as it allows matching better from a computational viewpoint the short- and large-time forms of the solution, so facilitating time partitioning. For this reason, the second deviation time may also be considered as a partitioning time. A first paper dealing with the analysis and implementation of the time partitioning method was proposed by Yen et al. [10] where, however, an algebraic expression was not given. In fact, the optimum partitioning time was estimated only a posteriori through several numerical experiments. Then, it was used for intrinsic verification purposes of numerical methods in transient heat conduction [11] and for improving convergence of infinite series involved in transient temperature solutions [4 (Sec. 5.3)].

The construction of the two-term short-time solution may be obtained by using: (1) some adequate physical considerations (due to thermal symmetries with respect to the middle plane of a slab) to simulate homogeneous boundary conditions of the first and second kinds on the back side of the slab (inactive surface); (2) the principle of superposition for linear problems [12 (Chapter 3)]; and (3) the concepts of penetration time and (first) deviation time that allow the well-known solutions of a semi-infinite body to be utilized as basic “building blocks” for constructing the short-time solution in terms of temperature and heat flux. Before proceeding to apply this methodology, the concept of “*computational*” analytical solution versus “*exact*” analytical solution and related relative error will be defined in Sec. 2.

The method here proposed applies to transient problems with active boundary conditions of the first, second, and third kinds on one side of the slab, and inactive boundary conditions of only the first and second kinds on the other side. In particular, the proposed method is used to find short-time solutions applicable at different times and different locations within the body with the accuracy desired for six different flat plates cases (Secs. 3–8). They are denoted by X11, X12, X21, X22, X31, and X32 in the heat conduction numbering system devised by Cole et al. [4 (Chapter 2)], where 1, 2, and 3 indicate the boundary conditions of first, second, and third kinds, respectively. See also Ref. [12 (Chapter 2 and Appendix A)] by Woodbury et al. Then, two examples on the maximum number of terms required to obtain a temperature solution with three different accuracies matching the two-term short-time form here derived and the well-known large-time form are given in Sec. 9. Some plots and tables are also provided.

## 2 Computational Analytical Solution

*computational*” analytical temperature solution is an approximate solution consisting of a finite number of terms, say $M$. It is denoted by $Tc(x,t,M)$, and replaces the “

*exact*” one, $Te(x,t)$, consisting of an infinite number of terms. If the plate is initially at a uniform temperature $Tin$, there is no heat source/sink within it (i.e., $g=0$), the $x=0$ boundary is the only “active” surface of the plate, and there the temperature can only change monotonically with time (heating or cooling), the error between the exact and computational analytical solutions at a location $x$ and time $t$ can be defined normalizing the absolute error with respect to the maximum temperature rise/fall that always occurs at $x=0$ and at time $t$. It results in [12 (Chapter 2)]

where $E=E(x,t,M)$ is hence the relative error.

where $T0(x,t)$ and $Tm(x,t)$ depend on (1) the form chosen for the analytical solution (short-time form or large-time form), (2) boundary condition type, and (3) boundary function prescribed.

Equation (2) states that, when $M\u2192\u221e$, $Tc(x,t,M)\u2192Te(x,t)$. Also, it shows that the numerator of Eq. (1) is the absolute value of the “tail” of the summation. This tail consisting of an infinite number of terms may be calculated conservatively by using the Euler-Maclaurin formula [13 (Eq. (4.14.1))]. An example of its use is given in Ref. [4 (Sec. 5.1.3)] for large-time solutions. The above formula can also be applied to the finite summation appearing at the denominator of Eq. (2).

By setting $E=10\u2212A$ in Eq. (1), where *A = *2, 3, …, 15 is the counting integer for the accuracy desired, and solving the resulting equation, the maximum number of terms $MA$ to get an accuracy of one part in $10A$ (with respect to the maximum temperature change $\Delta Tmax=|Tc(0,t,M)\u2212Tin|$) may be calculated as a function of time $t$, location $x$, and accuracy *A*. In detail, *A = *2 (accuracy of 1%) is for visual comparison and is acceptable in many engineering applications, while *A = *10 or even 15 is mainly for verification purposes of large numerical codes [14–16]. Note that the floating-point relative accuracy in double precision within MATLAB software is of 2.2204 × 10^{−16}. When using *A = *15, the “*computational*” analytical solution may be considered in every practical respect as the “*exact*” one.

where $Tc(0,t,MA)=Te(0,t)$ when the flat plate is subject to a boundary condition of the first kind.

Similarly, an error between computational and exact analytical heat flux solutions can be defined.

## 3 Slab With Prescribed Temperature on One Side and Insulated on the Other Side

Consider a slab of thickness *L* with temperature-independent properties, initially at zero temperature (without any loss of generality), subject to a step change in temperature ($0\u2192T0$) at *x = *0 and thermally insulated at *x = L*, as shown in Fig. 1, where *P* denotes a generic location of abscissa *x*. By using the numbering system devised by Cole et al. [4 (Chapter 2)], this problem may be denoted by X12B10T0, where “1” in B10 denotes a constant boundary condition at *x = *0, while “0” indicates a homogeneous boundary condition at *x = L*. Also, “T0” refers to a zero initial temperature within the plate.

Very simple solutions consisting of only two terms can be constructed and, hence, used to compute temperature and heat flux at any location of the slab at early times with the accuracy desired in place of the exact analytical solutions consisting of infinite terms, as shown in next subsections.

### 3.1 Penetration Time and First Deviation Time.

with a relative error (as defined by Eq. (1) for temperature) less than 10* ^{−A}*.

*. Note that “0” in X10 indicates a boundary condition of zeroth kind or of Beck type at*

^{−A}*x → ∞*, where there is no physical boundary [4 (Chapter 2); 12 (Chapter 2 and Appendix A)]. Therefore, temperature and heat flux may be calculated at any point

*P*of abscissa $0\u2264x\u2264L$ as [7]

which are the well-known solutions of the X10B1T0 problem, being $erfc(.)$ the complementary error function [4 (Appendix E)].

### 3.2 Second Deviation Time.

The slab of Fig. 1 is equivalent to a slab *2 L* thick and boundary condition like that at *x = *0 imposed at *x = 2 L*, i.e., X12B10T0 = X11B11T0, as shown in Fig. 3. In fact, if the slab is subject to a second kind of homogeneous boundary condition at *x = L* (Fig. 1), it can be simulated by a second layer in perfect thermal contact with the first one and having its same thickness, and subject to a nonzero boundary term at *x = 2 L* that is equal to that at *x = *0, as shown in Fig. 3.

*a*) and (6

*b*) as well as the superposition defined by Eq. (7), temperature and heat flux at any location $0\u2264x\u2264L$ of the slab of Fig. 1 can be calculated as

where *2 L* has been replaced by *4 L*, the slab may be approximated to a semi-infinite rectangular body ($x\u22650$), i.e., X11B10T0 ≈ X10B1T0 (see Fig. 5(a)), with relative errors *E* less than 10* ^{−A}*.

*P*of abscissa

*x*as

which are the solutions of the X10B1T0 problem of Fig. 5(a).

where *2 L − x* has been replaced by *2 L + x*, the slab may be approximate to a semi-infinite rectangular body ($x\u22642L$), i.e., X11B01T0 ≈ X01B1T0 (see Fig. 5(b)), with relative errors *E* less than 10* ^{−A}*.

*a*) and (10

*b*), for $0\u2264x\u22642L$ and $0\u2264t\u2264td(1b)$, the temperature and heat flux may be calculated as

which are the solutions of the X01B1T0 problem of Fig. 5(b) related to the ones of the companion X10B1T0 case of Fig. 5(a).

Also, note that $td(1a)$ defined by Eq. (9) is always greater than $td(1b)$ given by Eq. (11) with the only exception of the inactive back side *x = L* of the slab of Fig. 1, where the above two deviation times are the same. As $td(1b)=min{td(1a),td(1b)}$, for conservative purposes it may be taken as second deviation time, that is, $td(2)=td(1b)$.

*a*) and (8

*b*), for $td(1)<t\u2264td(2)$, the thermal field and related heat flux may be derived at any location $0\u2264x\u2264L$ of the flat plate of Fig. 1 with relative errors

*E*less than 10

*. In detail, substituting Eqs. (10*

^{-A}*a*), (10

*b*), (12

*a*), and (12

*b*) in Eqs. (8

*a*) and (8

*b*), respectively, yields

## 4 Slab With Prescribed Temperature on One Side and Zero Temperature on the Other Side

Consider a slab of thickness *L* with temperature-independent properties, initially at zero temperature, subject to a step change in temperature ($0\u2192T0$) at *x = *0 and kept at zero temperature at *x = L*, as shown in Fig. 6 where *P* denotes a generic location of abscissa *x*. By using the numbering system [4 (Chapter 2)], this problem may be denoted by X11B10T0.

Very simple solutions consisting of only two terms can be constructed and, hence, used to compute temperature and heat flux at any location of the slab at early times with the accuracy desired in place of the exact analytical solutions consisting of infinite terms, as shown in next subsections.

### 4.1 Penetration Time and First Deviation Time.

The short-time computational solutions to the current problem, with a relative error $E\u226410\u2212A$, being *A = *2, 3, …, 15, are [7]

- $0\u2264t\u2264tp$, where $tp=x2/(10A\alpha )$$TX11B10T0(x,t)|0\u2264x\u2264L\u22480$(14a)$qX11B10T0(x,t)|0\u2264x\u2264L\u22480$(14b)
- $tp<t\u2264td(1)=(2L\u2212x)2/(10A\alpha )$:$TX11B10T0(x,t)|0\u2264x\u2264L\u2248T0erfc(x2\alpha t)$(15a)$qX11B10T0(x,t)|0\u2264x\u2264L\u2248kT0\pi \alpha te\u2212x24\alpha t$(15b)

### 4.2 Second Deviation Time.

The slab of Fig. 6 is equivalent to a slab *2 L* thick and boundary condition opposite to that at *x = *0 imposed at *x = 2 L*, i.e., X11B10T0 = X11B11T0, as shown in Fig. 7. In fact, if the slab is subject to a first kind of homogeneous boundary condition at *x = L* (Fig. 6), it can be simulated by a second layer in perfect thermal contact with the first one and having its same thickness, and subject to a nonzero boundary term at *x = 2 L* that is opposite to that at *x = *0, as shown in Fig. 7.

*a*) and (16

*b*), temperature and heat flux at any location $0\u2264x\u2264L$ of the slab of Fig. 6 can be calculated as

Consider now the X11B10T0 problem of Fig. 4(a). If $0\u2264t\u2264td(1a)$, where $td(1a)$ is the deviation time given by Eq. (9), the slab may be approximated to a semi-infinite rectangular body ($x\u22650$), i.e., X11B10T0 ≈ X10B1T0 (see Fig. 5(a)), with relative errors *E* less than 10* ^{−A}*.

Therefore, for $0\u2264x\u22642L$ and $0\u2264t\u2264td(1a)$, temperature and heat flux may be calculated at any point *P* of abscissa *x* using Eqs. (10*a*) and (10*b*), which are the solutions of the X10B1T0 problem of Fig. 5(a).

Consider now the X11B01T0 problem of Fig. 4(b). If $0\u2264t\u2264td(1b)$, where $td(1b)$ is the deviation time given by Eq. (11), the slab may be approximate to a semi-infinite rectangular body ($x\u22642L$), i.e., X11B01T0 ≈ X01B1T0 (see Fig. 5(b)), with relative errors *E* less than 10* ^{−A}*.

Therefore, by changing $x\u2192(2L\u2212x)$ in Eqs. (10*a*) and (10*b*), for $0\u2264x\u22642L$ and $0\u2264t\u2264td(1b)$ the temperature and heat flux may be calculated through Eqs. (12*a*) and (12*b*), respectively, which are the solutions of the X01B1T0 problem of Fig. 5(b).

*a*) and (18

*b*), for $td(1)<t\u2264td(2)$ where $td(2)=td(1b)$, the thermal field and related heat flux may be derived at any location $0\u2264x\u2264L$ of the flat plate with relative errors

*E*less than 10

*. In detail, substituting Eqs. (10*

^{−A}*a*),(10

*b*), (12

*a*), and (12

*b*) in Eqs. (18

*a*) and (18

*b*) yields

## 5 Slab With Prescribed Heat Flux on One Side and Insulated on the Other Side

Consider a slab of thickness *L* with temperature-independent properties, initially at zero temperature, subject to a step change in heat flux ($0\u2192q0$) at *x = *0 and thermally insulated at *x = L*, as shown in Fig. 1, where $T0$ has to be replaced with $q0$. By using the numbering system, this problem may be denoted by X22B10T0.

The short-time computational solutions to the current problem for $0\u2264x\u2264L$ and for a relative error $E\u226410\u2212A$, with *A = *2, 3, …, 15, are

- $0\u2264t\u2264tp$, where $tp=x2/(10A\alpha )$$TX22B10T0(x,t)|0\u2264x\u2264L\u22480$(20a)$qX22B10T0(x,t)|0\u2264x\u2264L\u22480$(20b)
- $tp<t\u2264td(1)$, where $td(1)=(2L\u2212x)2/(10A\alpha )$$TX22B10T0(x,t)|0\u2264x\u2264L\u2248TX20B1T0(x,t)|x\u22650$(21a)$qX22B10T0(x,t)|0\u2264x\u2264L\u2248qX20B1T0(x,t)|x\u22650$(21b)
- $td(1)<t\u2264td(2)$, where $td(2)=(2L+x)2/(10A\alpha )$$TX22B10T0(x,t)|0\u2264x\u2264L\u2248TX20B1T0(x,t)|x\u22650+TX20B1T0(2L\u2212x,t)|x\u22650$(22a)$qX22B10T0(x,t)|0\u2264x\u2264L\u2248qX20B1T0(x,t)|x\u22650\u2212qX20B1T0(2L\u2212x,t)|x\u22650$(22b)

where the $ierfc(.)$function is the complementary error function integral [4 (Appendix E)].

## 6 Slab With Prescribed Heat Flux on One Side and Zero Temperature on the Other Side

Consider a slab of thickness *L* with temperature-independent properties, initially at zero temperature, subject to a step change in heat flux ($0\u2192q0$) at *x = *0 and kept at zero temperature at *x = L*, as shown in Fig. 6, where $T0$ has to be replaced with $q0$. By using the numbering system, this problem may be denoted by X21B10T0.

The short-time computational solutions to the current problem for $0\u2264x\u2264L$ and for a relative error $E\u226410\u2212A$, with *A = *2, 3, …, 15, are (see Sec. 4.1)

- $0\u2264t\u2264tp$, where $tp=x2/(10A\alpha )$$TX21B10T0(x,t)|0\u2264x\u2264L\u22480$(24a)$qX21B10T0(x,t)|0\u2264x\u2264L\u22480$(24b)
- $tp<t\u2264td(1)$, where $td(1)=(2L\u2212x)2/(10A\alpha )$$TX21B10T0(x,t)|0\u2264x\u2264L\u2248TX20B1T0(x,t)|x\u22650$(25a)$qX21B10T0(x,t)|0\u2264x\u2264L\u2248qX20B1T0(x,t)|x\u22650$(25b)
- $td(1)<t\u2264td(2)$, where $td(2)=(2L+x)2/(10A\alpha )$$TX21B10T0(x,t)|0\u2264x\u2264L\u2248TX20B1T0(x,t)|x\u22650\u2212TX20B1T0(2L\u2212x,t)|x\u22650$(26a)$qX21B10T0(x,t)|0\u2264x\u2264L\u2248qX20B1T0(x,t)|x\u22650+qX20B1T0(2L\u2212x,t)|x\u22650$(26b)

## 7 Slab With Prescribed Fluid Temperature on One Side and Insulated on the Other Side

Consider a slab of thickness *L* with temperature-independent properties, initially at zero temperature, subject to a step change in fluid temperature ($0\u2192T\u221e,0$) at *x = *0 and thermally insulated at *x = L*, as shown in Fig. 1, where $T0$ has to be replaced with $T\u221e,0$. Also, the heat exchange by convection at that boundary is characterized by a coefficient $h0$. By using the numbering system, this problem may be denoted by X32B10T0.

The short-time computational solutions to the current problem for $0\u2264x\u2264L$ and for a relative error $E\u226410\u2212A$, with *A = *2, 3, …, 15, are

- $0\u2264t\u2264tp$, where $tp=x2/(10A\alpha )$:$TX32B10T0(x,t)|0\u2264x\u2264L\u22480$(27a)$qX32B10T0(x,t)|0\u2264x\u2264L\u22480$(27b)
- $tp<t\u2264td(1)$, where $td(1)=(2L\u2212x)2/(10A\alpha )$:$TX32B10T0(x,t)|0\u2264x\u2264L\u2248TX30B1T0(x,t)|x\u22650$(28a)$qX32B10T0(x,t)|0\u2264x\u2264L\u2248qX30B1T0(x,t)|x\u22650$(28b)
- $td(1)<t\u2264td(2)$, where $td(2)=(2L+x)2/(10A\alpha )$$TX32B10T0(x,t)|0\u2264x\u2264L\u2248TX30B1T0(x,t)|x\u22650+TX30B1T0(2L\u2212x,t)|x\u22650$(29a)$qX32B10T0(x,t)|0\u2264x\u2264L\u2248qX30B1T0(x,t)|x\u22650\u2212qX30B1T0(2L\u2212x,t)|x\u22650$(29b)

The $erfcx(z)=ez2erfc(z)$ is the scaled complementary error function available in Matlab ambient. It is recommended for large values of *x*, *t,* and *h*_{0} as difficulties might arise during the numerical evaluation of Eq. (30*c*) in its former expression because of the positive argument of the exponential term (depending on the variables listed before). On the contrary, the scaled complementary error function (appearing in the latter expression) approaches zero for large values of *x*, *t* and *h*_{0}. A proof of this statement may be derived using the asymptotic expansion of [4 (Eq. (E.6))].

## 8 Slab With Prescribed Fluid Temperature on One Side and Zero Temperature on the Other Side

Consider a slab of thickness *L* with temperature-independent properties, initially at zero temperature, subject to a step change in fluid temperature ($0\u2192T\u221e,0$) at *x = *0 and kept at zero temperature at *x = L*, as shown in Fig. 6, where $T0$ has to be replaced with $T\u221e,0$. Also, the heat exchange by convection is characterized by $h0$. By using the numbering system, this problem may be denoted by X31B10T0.

The short-time computational solutions to the current problem for $0\u2264x\u2264L$ and for a relative error $E\u226410\u2212A$, with *A = *2, 3, …, 15, are

- $0\u2264t\u2264tp$, where $tp=x2/(10A\alpha )$$TX31B10T0(x,t)|0\u2264x\u2264L\u22480$(31a)$qX31B10T0(x,t)|0\u2264x\u2264L\u22480$(31b)
- $tp<t\u2264td(1)$, where $td(1)=(2L\u2212x)2/(10A\alpha )$$TX31B10T0(x,t)|0\u2264x\u2264L\u2248TX30B1T0(x,t)|x\u22650$(32a)$qX31B10T0(x,t)|0\u2264x\u2264L\u2248qX30B1T0(x,t)|x\u22650$(32b)
- $td(1)<t\u2264td(2)$, where $td(2)=(2L+x)2/(10A\alpha )$:$TX31B10T0(x,t)|0\u2264x\u2264L\u2248TX30B1T0(x,t)|x\u22650\u2212TX30B1T0(2L\u2212x,t)|x\u22650$(33a)$qX31B10T0(x,t)|0\u2264x\u2264L\u2248qX30B1T0(x,t)|x\u22650+qX30B1T0(2L\u2212x,t)|x\u22650$(33b)

## 9 Results and Discussion

The penetration time and the first and second deviation times defined in the previous sections are given in dimensionless form by Table 1 for different values of the accuracy *A* and at two different locations, namely, the front side $x\u02dc=0$ of the slab and its back side $x\u02dc=1$, with $x\u02dc=x/L$. In detail, $t\u02dcp=x\u02dc2/(10A)$, $t\u02dcd(1)=(2\u2212x\u02dc)2/(10A)$, and $t\u02dcd(2)=(2+x\u02dc)2/(10A)$, where the dimensionless time may be defined as $t\u02dc=\alpha t/L2$.

$x\u02dc=0$ | $x\u02dc=1$ | |||||
---|---|---|---|---|---|---|

A | $t\u02dcp$ | $t\u02dcd(1)$ | $t\u02dcd(2)$ | $t\u02dcp$ | $t\u02dcd(1)$ | $t\u02dcd(2)$ |

2 | 0 | 0.200 | 0.200 | 0.050 | 0.050 | 0.450 |

3 | 0 | 0.134 | 0.134 | 0.033 | 0.033 | 0.300 |

4 | 0 | 0.100 | 0.100 | 0.025 | 0.025 | 0.225 |

5 | 0 | 0.080 | 0.080 | 0.020 | 0.020 | 0.180 |

7 | 0 | 0.056 | 0.056 | 0.014 | 0.014 | 0.128 |

10 | 0 | 0.040 | 0.040 | 0.010 | 0.010 | 0.090 |

15 | 0 | 0.026 | 0.026 | 0.007 | 0.007 | 0.060 |

$x\u02dc=0$ | $x\u02dc=1$ | |||||
---|---|---|---|---|---|---|

A | $t\u02dcp$ | $t\u02dcd(1)$ | $t\u02dcd(2)$ | $t\u02dcp$ | $t\u02dcd(1)$ | $t\u02dcd(2)$ |

2 | 0 | 0.200 | 0.200 | 0.050 | 0.050 | 0.450 |

3 | 0 | 0.134 | 0.134 | 0.033 | 0.033 | 0.300 |

4 | 0 | 0.100 | 0.100 | 0.025 | 0.025 | 0.225 |

5 | 0 | 0.080 | 0.080 | 0.020 | 0.020 | 0.180 |

7 | 0 | 0.056 | 0.056 | 0.014 | 0.014 | 0.128 |

10 | 0 | 0.040 | 0.040 | 0.010 | 0.010 | 0.090 |

15 | 0 | 0.026 | 0.026 | 0.007 | 0.007 | 0.060 |

Therefore, $t\u02dcd(2)=t\u02dcd(1)$ for $x\u02dc=0$, while $t\u02dcd(2)$ is nine times larger than $t\u02dcd(1)$ for $x\u02dc=1$. This indicates that the short-time computational solution consisting of two semi-infinite terms is applicable for relatively long times at the back side of the plate when dealing with smaller accuracies such as *A = *2 (1%) or 3 (0.1%). These are however acceptable in many engineering applications. Also, $t\u02dcd(1)=t\u02dcp$ for $x\u02dc=1$.

### 9.1 X12B10T0 Case.

where the dimensionless temperature is taken as $T\u02dc=T/T0$. Also, the dimensionless heat flux is $q\u02dcc(x\u02dc,t\u02dc,MA)=\u2212\u2202T\u02dcc/\u2202x\u02dc$.

where the function “$ceil(z)$” is a MATLAB function that rounds the number *z* to the nearest integer greater than or equal to *z*. When using *A = *15, the “*computational*” analytical solution may be considered in every practical respect as the “*exact*” one.

A numerical comparison between the two-term short-time computational solution Eq. (13*a*) and the corresponding exact (*A =* 15) large-time solution Eq. (34) is shown in Table 2. The comparison is performed at the back side of the slab, $x\u02dc=1$, for different values of the dimensionless time. At early times, less than or equal to 0.05, the two solutions have the same fifteen decimal digits as expected. In fact, the above time of 0.05 is less than $t\u02dcd(2)=0.06$ given in Table 1 for *A = *15 and $x\u02dc=1$. For times larger than 0.05, the two-term short-time solution starts to deviate from the exact one giving a higher numerical value, as pointed out by the decimal digits in bold type in Table 2. However, it is still applicable up to the time of nearly 0.50 that may be considered as the partitioning time at $x\u02dc=1$ with the large-time solution Eq. (34) but only for visual purposes. In fact, $t\u02dcd(2)=0.45$ for *A = *2 in Table 1.

$t\u02dc$ | $T\u02dcc(1,t\u02dc)$ Eq. (13a) | $T\u02dce(1,t\u02dc)$ Eq. (34) with Eq. (35) |
---|---|---|

0.01 | 0.000,000,000,003,075 | 0.000,000,000,003,075 |

0.02 | 0.000,001,146,606,288 | 0.000,001,146,606,288 |

0.03 | 0.000,089,114,181,208 | 0.000,089,114,181,208 |

0.05 | 0.003,130,804,516,005 | 0.003,130,804,516,005 |

0.07 | 0.015,052,630,332,916 | 0.015,052,630,332,914 |

0.10 | 0.050,694,637,354,937 | 0.050,694,637,315,530 |

0.20 | 0.227,692,596,013,316 | 0.227,688,393,141,409 |

0.30 | 0.393,411,204,917,894 | 0.393,196,182,780,912 |

0.50 | 0.634,621,015,725,828 | 0.629,222,570,200,476 |

0.70 | 0.796,049,439,013,876 | 0.773,637,283,867,688 |

1.00 | 0.959,000,244,373,907 | 0.892,022,955,555,891 |

$t\u02dc$ | $T\u02dcc(1,t\u02dc)$ Eq. (13a) | $T\u02dce(1,t\u02dc)$ Eq. (34) with Eq. (35) |
---|---|---|

0.01 | 0.000,000,000,003,075 | 0.000,000,000,003,075 |

0.02 | 0.000,001,146,606,288 | 0.000,001,146,606,288 |

0.03 | 0.000,089,114,181,208 | 0.000,089,114,181,208 |

0.05 | 0.003,130,804,516,005 | 0.003,130,804,516,005 |

0.07 | 0.015,052,630,332,916 | 0.015,052,630,332,914 |

0.10 | 0.050,694,637,354,937 | 0.050,694,637,315,530 |

0.20 | 0.227,692,596,013,316 | 0.227,688,393,141,409 |

0.30 | 0.393,411,204,917,894 | 0.393,196,182,780,912 |

0.50 | 0.634,621,015,725,828 | 0.629,222,570,200,476 |

0.70 | 0.796,049,439,013,876 | 0.773,637,283,867,688 |

1.00 | 0.959,000,244,373,907 | 0.892,022,955,555,891 |

Concerning this, when using $t\u02dcd(2)$ as partitioning time between the 2-term short-time computational solution and the $MA$-terms large-time one, the number of terms to be used can be optimized, as shown by Fig. 8 where the back side of the plate is considered. For *A = *3, the partitioning time $t\u02dcd(2)=0.30$ (Table 1) and the maximum number of required terms is 4 for $t\u02dc\u2208[0,\u221e)$. For *A = *15, when $t\u02dcd(2)=0.06$ (Table 1), it is amazing that the maximum number of terms is only 10 for the entire time domain.

In addition, temperature and heat flux plots can be obtained by using the zero-term, one-term, and two-term short-time computational solutions derived in Sec. 3. As an example, in Fig. 9 the temperature curves at times 0.01 and 0.05 were obtained using Eqs. (4*a*) and (5*a*) with a maximum relative error of 0.01 (1%) with respect to the maximum temperature rise at $x=0$ and at the same time, while the temperature profile at time 0.2 was derived using Eq. (13*a*) with the same error. As regards the curve at time 1.0, it was obtained using the computational (*A = *2) large-time solution Eq. (34) with the companion expression Eq. (35).

### 9.2 X22B10T0 Case.

As for the companion X12B10T0 case of Sec. 9.1, a numerical comparison between the two-term short-time computational solution defined through Eq. (22*a*) with Eq. (23*a*) and the corresponding exact (*A =* 15) large-time solution Eq. (36) is performed at $x\u02dc=1$ for different values of the time. The results are shown in Table 3. At early times, less than or equal to 0.07, the two solutions have the same fifteen decimal place digits. As a matter of fact, the above time of 0.07 is greater than $t\u02dcd(2)=0.06$ given in Table 1 for *A = *15 and $x\u02dc=1$, but it is here recalled that $t\u02dcd(2)=(2+x\u02dc)2/(10A)$ is a conservative equation. For times larger than 0.07, the two-term short-time solution starts to deviate from the exact one giving lower numerical values, as pointed out by the decimal digits in bold type in Table 3.

$t\u0303$ | $T\u0303c(1,t\u0303)$ Eq. (22a) with Eq. (23a) | $T\u0303e(1,t\u0303)$ Eq. (36) with Eq. (37) |
---|---|---|

0.01 | 0.000,000,000,000,059 | 0.000,000,000,000,059 |

0.02 | 0.000,000,042,769,324 | 0.000,000,042,769,324 |

0.03 | 0.000,004,841,922,763 | 0.000,004,841,922,763 |

0.05 | 0.000,269,342,125,003 | 0.000,269,342,125,003 |

0.07 | 0.001,734,727,736,641 | 0.001,734,727,736,641 |

0.10 | 0.007,885,292,892,769 | 0.007,885,292,895,291 |

0.20 | 0.061,463,232,255,774 | 0.061,463,751,294,332 |

0.30 | 0.143,785,838,895,800 | 0.143,824,426,976,219 |

0.50 | 0.333,261,882,350,745 | 0.334,790,713,466,261 |

0.70 | 0.525,029,907,309,022 | 0.533,535,779,677,794 |

1.00 | 0.798,564,913,496,983 | 0.833,343,814,642,229 |

$t\u0303$ | $T\u0303c(1,t\u0303)$ Eq. (22a) with Eq. (23a) | $T\u0303e(1,t\u0303)$ Eq. (36) with Eq. (37) |
---|---|---|

0.01 | 0.000,000,000,000,059 | 0.000,000,000,000,059 |

0.02 | 0.000,000,042,769,324 | 0.000,000,042,769,324 |

0.03 | 0.000,004,841,922,763 | 0.000,004,841,922,763 |

0.05 | 0.000,269,342,125,003 | 0.000,269,342,125,003 |

0.07 | 0.001,734,727,736,641 | 0.001,734,727,736,641 |

0.10 | 0.007,885,292,892,769 | 0.007,885,292,895,291 |

0.20 | 0.061,463,232,255,774 | 0.061,463,751,294,332 |

0.30 | 0.143,785,838,895,800 | 0.143,824,426,976,219 |

0.50 | 0.333,261,882,350,745 | 0.334,790,713,466,261 |

0.70 | 0.525,029,907,309,022 | 0.533,535,779,677,794 |

1.00 | 0.798,564,913,496,983 | 0.833,343,814,642,229 |

If $t\u02dcd(2)$ is utilized as partitioning time between the 2-term short-time computational solution and the $MA$-term large-time one, the maximum number of terms required for the entire time domain can considerably be reduced, as shown by Fig. 10 at $x=1$. For *A = *3, the partitioning time $t\u02dcd(2)=0.30$ (Table 1) and the maximum number of required terms is 3 for $t\u02dc\u2208[0,\u221e)$. For *A = *15, when $t\u02dcd(2)=0.06$ (Table 1), it is amazing that the maximum number of terms is only 9 for the entire time domain.

## 10 Extension to Time-Dependent Boundary Conditions

The method illustrated in Secs. 3–8 for constant boundary conditions can also be extended to any time-dependent boundary function applied at the active surface $x=0$ as the most conservative boundary condition is the nonhomogeneous constant one of the first kind [8].

As an example, consider a slab initially at zero temperature with a power variation of the boundary condition (BC) on one side, say $T0(t/tref)n/2$ (first kind BC), $q0(t/tref)n/2$ (second kind BC) or $T\u221e,0(t/tref)n/2$ (third kind BC), with $tref$ arbitrary reference time and $n=\u22121,0,1,2,3,\u2026$. Also, the slab is at zero temperature or thermally insulated on the other side. This problem may be denoted by XIJB30T0, with *I* = 1, 2, 3 and *J* = 1 or 2, where “3” in B30 indicates some power other than 1 of $t$ [4 (Chapter 2)]. Alternatively, it can also be denoted by XIJB(3,n/2)0T0, as proposed later by Cole et al. [17]. When $n=2$, the problem is denoted by XIJB20T0, where “2” in B20 indicates a linear-in-time variation (unit power of $t$) of the surface temperature, heat flux or fluid temperature at $x=0$.

The short-time computational solutions to the current XIJB(3,*n*/2)0T0 problem, with *I* = 1, 2, 3 and *J* = 1, 2, for $0\u2264x\u2264L$ and for a relative error $E\u226410\u2212A$, with *A = *2, 3, …, 15, may be taken as

- $0\u2264t\u2264tp$, where $tp=x2/(10A\alpha )$$TXIJB(3,n/2)0T0(x,t)|0\u2264x\u2264L\u22480$(38a)$qXIJB(3,n/2)0T0(x,t)|0\u2264x\u2264L\u22480$(38b)
- $tp<t\u2264td(1)$, where $td(1)=(2L\u2212x)2/(10A\alpha )$$TXIJB(3,n/2)0T0(x,t)|0\u2264x\u2264L\u2248TXI0B(3,n/2)T0(x,t)|x\u22650$(39a)$qXIJB(3,n/2)0T0(x,t)|0\u2264x\u2264L\u2248qXI0B(3,n/2)T0(x,t)|x\u22650$(39b)
- $td(1)<t\u2264td(2)$, where $td(2)=(2L+x)2/(10A\alpha )$$TXIJB(3,n/2)0T0(x,t)|0\u2264x\u2264L\u2248TXI0B(3,n/2)T0(x,t)|x\u22650+(\u22121)JTXI0B(3,n/2)T0(2L\u2212x,t)|x\u22650$(40a)$qXIJB(3,n/2)0T0(x,t)|0\u2264x\u2264L\u2248qXI0B(3,n/2)T0(x,t)|x\u22650+(\u22121)J\u22121qXI0B(3,n/2)T0(2L\u2212x,t)|x\u22650$(40b)

where *J* = 1 or 2, and $TXI0B(3,n/2)T0(x,t)|x\u22650$ and $qXI0B(3,n/2)T0(x,t)|x\u22650$ are the analytical solutions of the associated semi-infinite problem denoted by XI0B(3,*n*/2)T0, with *I* = 1, 2 or 3.

They are [17]:

- X10(B3,
*n*/2)T0 solution [17 (Eqs. (15), (17*a*), and (17*b*))]$TX10B(3,n/2)T0(x,t)|x\u22650=T0(ttref)n/2\Gamma (1+n2)\xd72ninerfc(x2\alpha t)$(41a)$qX10B(3,n/2)T0(x,t)|x\u22650=kT02\alpha t(ttref)n/2\Gamma (1+n2)\xd72nin\u22121erfc(x2\alpha t)$(41b)

In Appendix A of Ref. [17], the temperature equation, Eq. (41*a*), is written out for specific values of $n=\u22121,0,1,2,3,4$.

- X20(B3,
*n*/2)T0 solution [17 (Eqs. (21*a*) and (21*b*))]$TX20B(3,n/2)T0(x,t)|x\u22650=q0k4\alpha t(ttref)n/2\Gamma (1+n2)\xd72nin+1erfc(x2\alpha t)$(42a)In Appendix B of Ref. [17], the temperature solution Eq. (42$qX20B(3,n/2)T0(x,t)|x\u22650=q0(ttref)n/2\Gamma (1+n2)\xd72ninerfc(x2\alpha t)$(42b)*a*) is written out for specific values of $n=\u22121,0,1,2,\u2026,6$. - X30(B3,
*n*/2)T0 solution [17 (Eqs. (23*a*) and (23*b*))]$TX30B(3,n/2)T0(x,t)|x\u22650=T\u221e,0(\u22121)n+1(\alpha tref)n/2(kh0)n\Gamma (1+n2)\xd7[Am(x,t,h0)\u2212\u2211r=0n(\u2212h0k4\alpha t)rirerfc(x2\alpha t)]$(43a)$qX30B(3,n/2)T0(x,t)|x\u22650=T\u221e,0(\u22121)n(\alpha tref)n/2(kh0)n\Gamma (1+n2)\xd7h0[Am(x,t,h0)\u2212k\pi \alpha te\u2212x24\alpha t\u2212\u2211r=0n(\u2212h0k4\alpha t)r\u22121ir\u22121erfc(x2\alpha t)]$(43b)

where $Am(x,t,h0)$ is defined by Eq. (30*c*). Also, in Appendix C of Ref. [17] the temperature expression Eq. (43*a*) is written out with specific values of $n=\u22121,0,1,2$. As regards the heat flux solution Eq. (43*b*), it is not given in Ref. [17] and has been derived here for the current work.

In Eqs. (41*a*) and (43*b*)$in\u22121erfc(.)$, $inerfc(.)$, and $in+1erfc(.)$ are the repeated integrals of the error function complement [4 (Appendix E)], and $\Gamma (.)$ is the gamma function [13]. In addition, $inerfc(0)=[2n\Gamma (1+n/2)]\u22121$, as shown in Ref. [4 (Fig. 6.1)].

*a*)–(40

*b*) are still applicable, provided $TXI0B(3,n/2)T0(x,t)|x\u22650$ and $qXI0B(3,n/2)T0(x,t)|x\u22650$ be replaced with the sum of the effects of each term in the function Eq. (44) using the principle of superposition for linear problems. Therefore

## 11 Conclusions

Very simple computational solutions consisting of two semi-infinite terms were constructed using the concept of penetration time, deviation time, principle of superposition, and adequate physical considerations to simulate homogeneous boundary conditions of the first and second kinds on the back side of a slab. These solutions can be used to compute temperature and heat flux at any location of the slab at early times in place of the exact analytical solutions (consisting of an infinite number of terms and, hence, unusable) with the accuracy desired. In fact, the large improvements in accuracy can be obtained without excessive computational costs. This contrasts significantly with the additional effort required to reduce errors using the numerical finite element method though this method allows solving problems with irregular geometry and nonlinear material properties.

Increasing the number of terms of the above short-time computational solution should extend its time range of applicability. However, this is a subject for future research. Another future topic concerns the possibility to investigate how to simulate the homogenous boundary condition of the third kind at the inactive back side of the plane wall.

## Acknowledgment

The authors would like to thank Professor Em. James V. Beck for his relevant contribution to the topic of this paper. Sadly, he passed away on July 28, 2022 [18].

## Data Availability Statement

The authors attest that all data for this study are included in the paper.

## Nomenclature

- $A$ =
counting integer for accuracy

- $E$ =
relative error based on $\Delta Tmax$

- $h$ =
heat transfer coefficient, $Wm\u22122K\u22121$

- $k$ =
thermal conductivity, $Wm\u22121K\u22121$

- $L$ =
slab thickness, $m$

- $M$ =
maximum number of terms

- $q$ =
heat flux, $Wm\u22122$

- $t$ =
time, $s$

- $td$ =
deviation time, $s$

- $tp$ =
penetration time, $s$

- $T$ =
temperature, $K$

- $x$ =
space coordinate, $m$

- $\alpha $ =
thermal diffusivity, $m2s\u22121$

- $\Delta Tmax$ =
maximum temperature change defined in the text