Novel evaluation method for non-Fourier effects in heat pulse experiments
NNOVEL EVALUATION METHOD FOR NON-FOURIER EFFECTS IN HEAT PULSEEXPERIMENTS
A. FEH´ER , R. KOV ´ACS Abstract.
The heat pulse (flash) experiment is a well-known and widely accepted method to measure the thermaldiffusivity of a material. In recent years, it is observed that the thermal behavior of heterogeneous materials canshow deviation from the classical Fourier equation, resulting in a different thermal diffusivity and requiring furtherthermal parameters to identify. Such heterogeneity can be inclusions in metal foams, layered structure in composites,or even cracks and porous parts in rocks. Furthermore, the next candidate, the so-called Guyer-Krumhansl equation,is tested on these experiments with success. However, these recent evaluations required a computationally intensivefitting procedure using countless numerical solutions, even when a good initial guess for the parameters is found byhand. This paper presents a Galerkin-type discretization for the Guyer-Krumhansl equation, which helped us finda reasonably simple analytical solution for time-dependent boundary conditions. Utilizing this analytical solution,we developed a new evaluation technique to immediately estimate all the necessary thermal parameters using themeasured temperature history. Introduction
The engineering practice requires reliable ways to determine the necessary parameters, which are enough tocharacterize the material behavior. In what follows, we place our focus on the thermal description of materials,especially on heterogeneous materials such as rocks and foams. In recent papers [1, 2], it is reported that thepresence of various heterogeneities can result in a non-Fourier heat conduction effect on macro-scale under roomtemperature conditions. A particular one is depicted in Fig. 1 for a capacitor sample having a periodic layeredstructure. Such effects are observed in a so-called flash (or heat pulse) experiment in which the front side of thespecimen is excited with a short heat pulse, and the temperature is measured at the rear side. That temperaturehistory is used to find the thermal diffusivity in order to characterize the transient material behavior.
Figure 1.
Measured rear side temperature history for the capacitor sample and the predictionprovided by Fourier’s theory.This non-Fourier effect occurs on a specific time interval as Fig. 1 shows for a typical outcome of the flashexperiments; this is called over-diffusion. After that interval, the Fourier equation appears to be a suitable choice formodeling, the influence of the heterogeneities vanishes (later, we show further examples). Also, there is no difference
Date : January 5, 2021. a r X i v : . [ phy s i c s . a pp - ph ] J a n A. FEH´ER , R. KOV´ACS between the steady-states of Fourier and non-Fourier heat equations. Based on our experimental experience, theexistence of over-diffusion depends on various factors, for instance, sample thickness, characteristic parallel timescales, and excitation (i.e., boundary conditions) [3].In the following sections, we organize the discussion as follows. First, we briefly introduce the two heat conductionmodels to model the heat pulse experiments and used for evaluations with a particular set of dimensionless quantities.Second, we shortly present how the complete evaluation with the Fourier heat equation can be conducted. Then,we move on the way of the evaluation procedure with the Guyer-Krumhansl equation. After, we demonstrate thebenefits of this fitting procedure and revisit some previous measurements. Furthermore, we decided to place thederivation of the analytical solutions to the end of the paper as an Appendix. According to our knowledge, theGalerkin method has not been used before for the Guyer-Krumhansl equation, and is a novel result in this respect,we want to keep the focus on its practical utilizations.2. Models for heat pulse experiments
Although numerous generalizations of Fourier’s law exist in the literature [4], there is solely one of them, whichindeed proved to be reasonable as the next candidate beyond Fourier’s theory, this is called Guyer-Krumhansl (GK)equation, this constitutive equation reads in one spatial dimension τ q ∂ t q + q + λ∂ x T − κ ∂ xx q = 0 . (1)Here, τ q is the relaxation time for the heat flux q and κ is a kind of ‘dissipation parameter’, usually related tothe mean free path. Whereas it was first derived on the basis of kinetic theory [5], this model also has a strongbackground in non-equilibrium thermodynamics with internal variables (NET-IV) [6,7]. While in the first case, oneassumes an underlying mechanism for phonons as the kinetic theory requires it, this is entirely neglected in the caseof NET-IV, leaving the coefficients to be free (however, their sign is restricted by the II. law of thermodynamics).Eq. (1) is a time evolution equation for the heat flux, and in order to have a mathematically and physically completesystem, we need the balance of internal energy e , too, ρc∂ t T + ∂ x q = 0 , (2)in which the equation of state e = cT is used with c being the specific heat and rho is the mass density. All thesecoefficients are constant, only rigid bodies are assumed with no heat source.At this point, we owe an explanation of why we leave the Maxwell-Cattaneo-Vernotte (MCV) equation out ofsight.(1) Hyperbolicity vs. parabolicity. It is usually claimed that a heat equation model should be hyperbolic suchas the MCV theory, describing finite propagation speed. Indeed, this seems reasonable, but it does nothelp in the practical applications under common conditions (room temperature, heterogeneous materials).The Fourier equation is still well-applicable in spite of its parabolic nature, therefore we do not see it as adecisive property.(2) In a low-temperature situation, the MCV model was useful, primarily due to the observed wave phenomenonin super-fluids, called second sound [8, 9]. Despite the GK equation’s parabolic nature, it also helped theresearchers find the second sound in solids as well [10].(3) There is a significant effort to find the trace of wave propagation at room temperature (in a macro-scaleobject, so nano-structures does not count now), sadly with no success [11, 12].(4) There are higher-order models as well, such as ballistic-diffusive models [13–15], but they are related to adifferent research program, and for this work, investigating macro-scale objects, they do not seem relevant.(5) On the analogy of the MCV model, the so-called dual-phase lag (DPL) equation [16] usually used in manyworks as the best candidate after Fourier’s law. Sadly, this model introduces two time constants in anad hoc manner, violating basic physical principles [17, 18], leading to mathematically ill-posed problems aswell [19, 20].Last but not least, we also must mention a relatively less-known model from the literature, the Ny´ıri equation [21], q + λ∂ x T − κ ∂ xx q = 0 , (3)which one is indeed similar to the Guyer-Krumhansl model but leaves the time lagging effects out of sight, henceit is purely a spatially nonlocal heat equation. Testing its solutions with the method presented in the Appendix, itturned out to be inaccurate for measurements, unfortunately. Consequently, the GK model is indeed the simplestbut necessary extension for the Fourier equation, neither the MCV nor the Ny´ıri models are capable of describingthese experiments accurately. In other words, the two new parameters ( τ q and κ ) are truly needed. OVEL EVALUATION METHOD FOR NON-FOURIER EFFECTS IN HEAT PULSE EXPERIMENTS 3
T and q-representations.
Depending on the purpose, it is useful to keep in mind that for such linear models,it is possible to chose a ‘primary’ field variable, which could ease the definition of boundary conditions in somecases. For the GK equation, the temperature T and the heat flux q are the candidates, and their forms areT-representation: τ q ∂ tt T + ∂ t T − α∂ xx T − κ ∂ txx T = 0 , (4)q-representation: τ q ∂ tt q + ∂ t q − α∂ xx q − κ ∂ txx q = 0 . (5)We note that in T -representation, it is unknown how to define boundary condition for q since it requires knowledgeon ∂ xx q . On the other hand, in q -representation, it becomes meaningless to speak about T -boundaries. In a previousanalytical solution for the GK equation [22], this difference was inevitable to realize. In the present work, we use thesystem (1)-(2). It is also interesting to notice that the GK model can recover the solution of the Fourier equationwhen κ /τ q = α , this is called Fourier resonance [1, 23]. Overall, the coefficients τ q , α , and κ must be fitted to thegiven temperature history.2.2. Dimensionless set of parameters.
Following [1], we introduce these definitions for the dimensionless pa-rameters (quantities with hat):time: ˆ t = tt p and ˆ x = xL ;thermal diffusivity: ˆ α = αt p L with α = λρc ;temperature: ˆ T = T − T T end − T with T end = T + ¯ q t p ρcL ;heat flux: ˆ q = q ¯ q with ¯ q = 1 t p (cid:90) t p q ( t )d t ;heat transfer coefficient: ˆ h = h t p ρc ; (6)together with ˆ τ q = τ q t p , ˆ κ = κ L , where ˆ t differs from the usual Fourier number in order to decouple the thermaldiffusivity from the time scale in the fitting procedure. Furthermore, t p denotes the constant heat pulse durationfor which interval ¯ q averages the heat transferred with the heat pulse defined by q ( t ). Here, L is equal with thesample thickness. T end represents the adiabatic steady-state, and T is the uniform initial temperature. In the restof the paper, we shall omit the hat notation, otherwise we add the unit for the corresponding quantity. Utilizingthis set of definitions, one obtains the dimensionless GK model: ∂ t T + ∂ x q = 0 ,τ q ∂ t q + q + α∂ x T − κ ∂ xx q = 0 . (7)The initial condition is zero for both fields. For further details, we refer to the Appendix in which we present theanalytical solution for the two heat equations. This set of dimensionless parameters does not change the definitionof the Fourier resonance condition, i.e., it remains ˆ κ / ˆ τ q = ˆ α .3. Evaluation with the Fourier theory
The analytical solution of the Fourier equation is found for the rear side in the form of T ( x = 1 , t ) = Y exp( − ht ) − Y exp( x F t ) , x F = − h − απ , t > , (8)where all the coefficients are expressed in detail in the Appendix. First, we must estimate the heat transfercoefficient h by choosing arbitrarily two temperature values at the decreasing part of temperature history. In thisregion, exp( x F t ) ≈
0, thus h = − ln( T /T ) t − t . (9)For the Fourier theory, it is possible to express the thermal diffusivity explicitly, i.e., α F = 1 . L π t / , (10)and after registering t / , it can be directly determined. This is the ratio of the thermal conductivity λ and thespecific heat capacity ρc . Then, the top of the temperature history ( T max ) follows by reading the time instant( t max ) when T max occurs. Figure 2 schematically summarizes this procedure. Overall, we obtained the heat transfercoefficients, the thermal diffusivity and T max , which all used for the Guyer-Krumhansl theory. A. FEH´ER , R. KOV´ACS Figure 2.
Schematically presenting the evaluation method using Fourier’s theory.4.
Evaluation with the Guyer-Krumhansl theory
The situation here becomes more difficult since this non-Fourier theory consists of two ‘time constants’ ( x and x ) instead of one ( x F in the Fourier theory). Consequently, it is not possible to find these exponents withoutmaking simplifications, in which one must be immensely careful. We prepared ‘parameter maps’ for all possible τ q and κ values that could be practically possible and beyond in order to check the effect of the simplifications madein the following. However, we still had to restrict ourselves to a domain, which is 3 > κ / ( ατ q ) ≥
1. Its lower limitexpresses the Fourier case, and any other combination falls on the over-diffusive region. The highest experimentallyobserved ratio so far is around 2 .
5, thus we expect 3 to be eligible. For κ , we consider 0 . < κ <
1. We wantto emphasize that the GK theory itself is not restricted on this domain, it would allow under-diffusive (‘wave-like’)propagation as well [10]. However, for the present situation, we consider it out of interest in the lack of experimentalobservation for room temperature experiments on macro-scale heterogeneous samples. In the GK theory, we canexpress the rear side temperature history as T ( x = 1 , t ) = Y exp( − ht ) − Z exp( x t ) − Z exp( x t ) , x , x < , (11)for the detailed calculation and parameter definitions, we refer to the Appendix again. This can be equivalentlyformulated realizing that Z = − P − Z , T ( x = 1 , t ) = Y exp( − ht ) − Z (cid:0) exp( x t ) − exp( x t ) (cid:1) + P exp( x t ) , (12)where merely one simplification becomes possible for all τ q and κ : exp( x t ) (cid:29) exp( x t ) when t >
60, i.e., T ( x = 1 , t >
60) = Y exp( − ht ) − Z exp( x t ) + P exp( x t ) . (13)This form is more advantageous because P remains practically constant for a given boundary condition, thus itsvalue can be assumed a-priori, this is exploited in the evaluation method. Now, let us present from step by stepthe determination of GK parameters, depicted on Fig. 3. • Step 1/A. We have to observe that the temperature predicted by Fourier’s theory always runs together withthe measured one at the beginning, after that, it rises faster at the top. In other words, in this region thesame temperature value (usually around 0 . − .
95) is reached sooner. Mathematically, we can express itby formally writing the equations for the Fourier and GK theories as follows, T F = Y exp( − ht ) − Y exp( x F t F ); T GK = Y exp( − ht m ) − Z exp( x t m ) + P exp( x t m ) , (14)where the t F time instant is smaller than the measured t m , also T F = T GK holds. Let us choose such twotemperatures arbitrarily and taking their ratio, it yieldsexp (cid:0) x F ( t F − t F ) (cid:1) = exp (cid:0) x ( t m − t m ) (cid:1) − Z + P exp (cid:0) ( x − x ) t m (cid:1) − Z + P exp (cid:0) ( x − x ) t m (cid:1) (15) OVEL EVALUATION METHOD FOR NON-FOURIER EFFECTS IN HEAT PULSE EXPERIMENTS 5 where the fraction on the right hand side is close to 1, mostly between 1 and 1 .
05 for ‘small’ time intervals.It could be possible to introduce it as a correction factor (denoted with c below) for x in an iterativeprocedure if more to be known about τ q and κ . After rearrangement, we obtain a closed form formula for x : x = ln(1 /c ) t m − t m + x F t F − t F t m − t m . (16)Taking c = 1 is equivalent with neglecting exp( x t ) from the beginning around reaching T max , and leadingto this same expression. Eventually, it introduces a correction for the Fourier exponent x F based onthe deviation from the measured data with the possibility to apply further corrections using c if needed.Practically, we take the 80% and 90% of T max and for the next 20 subsequent measurement points, then weconsider their mean value to be x . From mathematical point of view, closer data point pairs should performbetter, but it does not due to the uncertainty in the measurement data. According to our experience, itoffers a more consistent value for x . Figure 3.
The schematic representation of the evaluation method using the Guyer-Krumhansltheory. Here, the ’fitted curve’ belongs to the Fourier equation. • Step 1/B. In parallel with part A, we can determine the coefficient Z for each t m and for each corresponding x ,m , that is, Z ,m = exp( − x ,m t m ) (cid:16) T m − Y exp( − ht m ) (cid:17) (17)where the subscript m denotes the value related to one measurement point. Also, after 20 subsequentpoints, we take the mean value of the set { Z ,m } . • Step 2. At this point, we can exploit that P is ‘almost constant’, i.e., 2 < − P < .
03 holds. Here, 2 . − P , and also, it cannot be smallerthan 2. This property allows us to a-priori assume its value (such as P = − . P , we can obtain Z = − P − Z . Inorder to obtain x , we can rearrange the equation T = Y exp( − ht ) − Z exp( x t ) − Z exp( x t ) (18)for x , and calculate it as a mean value of the set { x ,m } filled with values related to each t m . When havingnoisy data, this approach can result in positive x values, unfortunately. These values must be excluded,otherwise, it leads to instability and a meaningless outcome. Careful data filtering can help to solve thisshortcoming, and in fact, we used it to ease the calculation (the details are plotted in the next section). A. FEH´ER , R. KOV´ACS • Step 3. Now, having both exponents and coefficients, it is possible to rearrange the analytical expressionsto the GK parameters explicitly and calculate α GK , τ q and κ : x , x ⇒ k , k ; Z ⇒ DP ⇒ M , M ⇒ τ q ⇒ α GK ⇒ κ . (19)For the detailed parameter definitions, we refer to the Appendix. • Step 4. As it is mentioned in Step 2, the overall outcome is sensitive to P . Therefore we choose to makea sweep on the possible interval with the step of 0 . R , the coefficient of determination. Lastly, we chose the bestset.Practically, this evaluation method reduces the number of ‘fitted’ parameters as only P has to be fine-tuned at theend. Besides, it is constrained into a relatively narrow range, consequently, the overall evaluation procedure takesonly a few seconds instead of hours to perform computationally intensive algorithms.5. Comparison with foregoing experiments
First, we revisit the experiments presented in [3] since that set of data on Basalt rock samples with thicknessesof 1 .
86, 2 .
75 and 3 .
84 mm, showed size dependence both on the thermal diffusivity and on the non-Fourier effects.However, the fitted parameters in [3] found by hand, thus not exactly precise. Here, we aim to specify the exactquantities for the GK model and establishing a more robust theoretical basis for the observations. Second, wereevaluate the data recorded on a metal foam sample with 5 . Figure 4.
The magnified view of the metal foam specimen (center) and the basalt rock sample (right).In some cases, the available data is too noisy for such an evaluation method, an example is presented in Fig. 5.That data is smoothed using the built-in Savitzky-Golay algorithm of Matlab. Besides, we paied much attentionto not smooth it overly in order to keep the physical content untouched as much as possible.
Figure 5.
The effect of data smoothening, using 10 neighboring points in the Savitzky-Golay algorithm.
OVEL EVALUATION METHOD FOR NON-FOURIER EFFECTS IN HEAT PULSE EXPERIMENTS 7
Basalt rock samples.
Regarding the exact details of measurements, we refer to [3]. Tables 1 and 2 consistof our findings using this evaluation algorithm. Comparing the outcomes of the two fitting procedures, we find thethermal diffusivities to be close to each other. However, this is not the case with the GK parameters τ q and κ ,they significantly differ from the previous values from [3]. Despite the huge difference, the size-dependence for boththe Fourier and non-Fourier behaviors is apparent nevertheless. The fitted temperature histories are depicted inFigs. 6, 7, 8 and 9 for each thickness, respectively. Each figure shows the R value for the fitted curve. For theFourier one, two of them are given: R t represents the one found without any fine-tuned thermal diffusivity, this ispurely theoretical. The other R stands for fine-tuned α F .In the first case ( L = 1 .
86 mm), although the difference for the non-Fourier samples seems negligible, it results in10% difference in the thermal diffusivity. It is more visible from Table 2, in which the Fourier resonance conditionspectacularly characterize the deviation from Fourier’s theory, it decreases for thicker samples. Regarding the thirdone ( L = 3 .
84 mm), Fourier’s theory seems to be ‘perfectly splendid’, and the GK model hardly improves it. Indeed,the 0 .
94 for the Fourier resonance is close enough to 1 to consider it to be a Fourier-like propagation.Basalt rocksamples Findings in [3] Refined results α F − [m /s] α GK − [m /s] τ q [s] κ − [m ] α F − [m /s] α GK − [m /s] τ q [s] κ − [m ]1 .
86 mm 0 .
62 0 .
55 0 .
738 0 .
509 0 .
68 0 .
61 0 .
211 0 . .
75 mm 0 .
67 0 .
604 0 .
955 0 .
67 0 .
66 0 .
61 0 .
344 0 . .
84 mm 0 .
685 0 .
68 0 .
664 0 .
48 0 .
70 0 .
68 1 0 . Table 1.
Summarizing the fitted thermal parameters. κ τ q α Findingsin [3] Refinedvalues1 .
86 mm 1 .
243 1 . .
75 mm 1 .
171 1 . .
84 mm 1 .
06 0 . Table 2.
Characterizing the non-Fourier behavior using the Fourier resonance condition.
Figure 6.
The rear side temperature history for the basalt rock sample with L = 1 .
86 mm.
A. FEH´ER , R. KOV´ACS Figure 7.
The rear side temperature history for the basalt rock sample with L = 2 .
75 mm.
Figure 8.
Demonstrating the complete fitting for the rear side temperature in case of the basaltrock sample with L = 2 .
75 mm.
Figure 9.
The rear side temperature history for the basalt rock sample with L = 3 .
84 mm.
OVEL EVALUATION METHOD FOR NON-FOURIER EFFECTS IN HEAT PULSE EXPERIMENTS 9
Metal foam.
Regarding the extent of the non-Fourier effect, the situation becomes remarkably differentfor the metal foam sample, presented first in [2]. The millimeter size inclusions can significantly influence thethermal behavior. The outcome is plotted in Fig. 10 together with the corresponding R values. Table 3 helps tocompare the fitted values found by Wolfram Mathematica to ours. Notwithstanding that the GK parameters arein correspondence, the most notable difference is on the thermal diffusivities, interestingly. The Fourier resonanceparameter is found to be 2 .
395 with our procedure, while on the contrary to 3 .
04 in [2]. Common in both cases,the ratio of α F and α GK is found to be 1 . .
29, which represents an indeed remarkable deviation from Fourier’stheory.Metal foamsample Findings in [2] with Wolfram Math Present algorithm α F − [m /s] α GK − [m /s] τ q [s] κ − [m ] α F − [m /s] α GK − [m /s] τ q [s] κ − [m ]5 . .
04 2 .
373 0 .
402 2 .
89 3 .
91 3 .
01 0 .
304 2 . Table 3.
Summarizing the fitted thermal parameters for the metal foam sample.
Figure 10.
The rear side temperature history for the metal foam sample with L = 5 . Discussion and summary
We developed an algorithm to efficiently evaluate room temperature heat pulse experiments in which a non-Fourier effect could exist. This is called over-diffusive propagation and detunes the thermal diffusivity, even whenthe deviation is seemingly small or negligible for the rear side temperature history. The presented method is basedon the analytical solution of the Guyer-Krumhansl equation, including temperature-dependent convection boundarycondition, thus the heat transfer to the environment can be immediately included in the analysis. The reevaluationof preceding experiments showed a real size-dependence for all thermal parameters, especially for the GK coefficients τ q and κ . Furthermore, it is in accordance with the result of the iterative ‘brute force’ iterative fitting procedureof Wolfram Math, basically, but using much less computational resource.We plan to improve this procedure by including the investigation of front side temperature history, too. When x is obtained from the rear side, it could be easily used to describe the front side’s thermal behavior. This is muchmore sensitive to the initial time evolution right after the excitation, therefore it could serve as a better candidateto achieve a more precise and robust estimation for the x exponent. Also, having two temperature histories wouldbe a remarkable step forward to ascertain the existence of non-Fourier heat conduction.We believe that this procedure lays the foundations for the more practical engineering applications of non-Fouriermodels, especially for the best candidate among all of them, the Guyer-Krumhansl equation. It sheds new lighton the classical and well-known flash experiments, and we provide the necessary tools to find additional thermalparameters to achieve a better description of heterogeneous materials. It becomes increasingly important with thespreading of composites and foams and helps characterize 3D printed samples with various inclusions. , R. KOV´ACS Acknowledgement
The authors thank Tam´as F¨ul¨op, P´eter V´an, and M´aty´as Sz¨ucs for the valuable discussions. We thank L´aszl´oKov´acs (K˝om´er˝o Kft., Hungary) and Tam´as B´arczy (Admatis Kft.) for producing the rock and metal foam samples.The research reported in this paper and carried out at BME has been supported by the grants National Research,Development and Innovation Office-NKFIH FK 134277, and by the NRDI Fund (TKP2020 NC, Grant No. BME-NC) based on the charter of bolster issued by the NRDI Office under the auspices of the Ministry for Innovationand Technology. This paper was supported by the J´anos Bolyai Research Scholarship of the Hungarian Academyof Sciences. 8.
Appendix: Galerkin-type solution of heat equations
Since the non-Fourier models are not well-known in the general literature, there are only a few available analyticaland numerical methods to solve such a spatially nonlocal equation like the Guyer-Krumhansl one. The nonlocalproperty is a cornerstone of these models since the usual boundary conditions do not work in the same way. Thatcould be a problem when the outcome seemingly violates the maximum principle, see for instance [24] in which theoperational approach is applied [25].Another particular candidate originates from the spectral methods, it is called Galerkin method, where both theweight and the trial functions are the same. Fortunately, following [22], we can surely apply sine and cosine trialfunctions in which terms the solution can be expressed. It is important to emphasize that we deal with a system ofpartial differential equations in our case. The physical (and mathematical) connection between the field variablesrestricts the family of trial functions. Namely, even in the simplest case of the Fourier heat equation, ∂ t T + ∂ x q = 0 , q + α∂ x T = 0 , (20) q and T are orthogonal to each other, and the trial functions must respect this property. Our choice is found bythe method called separation of variables, but it resulted in a too complicated outcome due to the time-dependentboundary condition. In [22], the heat pulse is modeled with a smooth q ( t ) = 1 − cos(2 πt/t p ) function on the0 < t ≤ t p interval. It is disadvantageous since the most interesting part falls beyond t p , and the solution in t p < t must account for the state at the time instant t p as an initial condition. Therefore, it results in cumbersomeexpressions for the coefficients.We overcome this difficulty by introducing a different function to model the heat pulse, i.e., we use q ( t ) = − (exp( − C t ) + exp( − C t )) /n , where C and C are chosen to have a sufficiently small q after t p , hence the valuesare C = 1 / .
075 and C = 6. The coefficient n is normalizing q to 1 from 0 to t p , so it is n = ( C − C ) / ( C C ).For larger time instants, the front side becomes adiabatic. Regarding the rear side boundary condition, we chooseto account heat convection for both models.(1) First, we restrict ourselves to the Fourier equation. According to our previous experiments [1–3], one cansafely use the Fourier equation where cooling effects become significant. This solution is used to estimatethe heat transfer coefficient, the maximum temperature and give a first approximation to the thermaldiffusivity.(2) In the second step, we repeat the calculation for the GK model with the same boundary conditions. Weuse the previously found Fourier parameters as the input to estimate the GK parameters and fine-tune thethermal diffusivity. The heat transfer coefficient and the temperature maximum can be kept the same.8.1. Step 1: solving the Fourier equation.
While there are several available solutions in the literature, we wantto see how the Galerkin approach performs on this model using our set of dimensionless parameters and boundaryconditions. Consequently, we can keep our findings to be consistent between the two heat equation models. Let usrecall the mathematical model for the sake of traceability. In the Fourier model, we have ∂ t T + ∂ x q = 0 , q + α∂ x T = 0 , (21)with q ( t ) = − n (cid:16) exp( − C t ) + exp( − C t ) (cid:17) , n = C − C C C , q ( t ) = hT ( x = 1 , t ) , (22)in which all parameters are dimensionless as presented in Sec. 2.2. The initial conditions are q ( x, t = 0) = 0 and T ( x, t = 0) = 0, the conducting medium is thermally relaxed. We emphasize that one does not need to separatelyspecify boundary conditions for the temperature field T as well. Regarding the heat flux field q , we must separatethe time-dependent part from the homogeneous one, q ( x, t ) = w ( x, t ) + ˜ q ( x, t ) , w ( x, t ) = q ( t ) + x ( q ( t ) − q ( t )) (23) OVEL EVALUATION METHOD FOR NON-FOURIER EFFECTS IN HEAT PULSE EXPERIMENTS 11 with ˜ q being the homogeneous field and w inherits the entire time-dependent part from the boundary (and x runsfrom 0 to 1). The spectral decomposition of ˜ q and T are˜ q ( x, t ) = N (cid:88) j =1 a j ( t ) φ j ( x ) , T ( x, t ) = N (cid:88) j =1 b j ( t ) ϕ j ( x ) , (24)with φ j ( x ) = sin( jπx ) and ϕ j ( x ) = cos( jπx ). Revisiting the boundary conditions, q is trivial, and q becomes: q ( t ) = hT ( x = 1 , t ) = h (cid:80) Nj =0 b j ( − j . Naturally, one has to represent also w ( x, t ) in the space spanned by φ j ( x ). Once one substituted these expressions into (21), multiplied them by the corresponding weight functions andintegrated respect to x from 0 to 1, one obtains a system of ordinary differential equations (ODE). Here, we exploitthat the square of the trial functions φ ( x ) and ϕ ( x ) are both integrable and after integration they are equal to1 /
2. Since the cos series have a non-zero part for j = 0, we handle it separately from the others corresponding to j > • For j = 0, we have ˙ b + ∂ x w = 0 , → ˙ b = − hb + q (25)with the upper dot denoting the time derivative, and a = 0 identically. • For j >
0, we obtained ˙ b j + jπa j = 0 , (26) a j = αjπb j + 2 jπ ( hb j − q ) , (27)where the term 2 / ( jπ ) comes from the sin series expansion of w . We note that for j > ∂ x w does notcontribute to the time evolution.Such ODE can be solved easily both numerically and analytically for suitable q functions. Figure 11 shows theanalytical solution programmed in Matlab in order to demonstrate the convergence to the right (physical) solution.In this respect, we refer to [26] in which a thorough analysis is presented on the analytical and numerical solutionof heat equations beyond Fourier. Our interest is to utilize as less terms as possible of the infinite series, which isable to properly describe the rear side temperature history (i.e., the measured one) from a particular time instant.In other words, we want to simplify the complete solution as much as possible but keeping its physical meaning.Starting with j = 0 case, we find that the terms in the particular solution with exp( − C t ) and exp( − C t ) extinctvery quickly, thus we can safely neglect them with keeping the exp( − ht ) as the leading term throughout the entiretime interval we investigate. Briefly, j = 0 yields b ( t ) = Y exp( − ht ) , Y = ( C − C ) / (cid:0) n ( C − h )( C − h ) (cid:1) . (28)Continuing with j = 1, we make the same simplifications and neglecting the same exponential terms as previouslyafter taking into account the initial condition, and we found b ( t ) = Y exp( − ht ) exp( − π t ) , Y = 2( C − C ) / (cid:0) n ( C + x F )( C + x F ) (cid:1) , x F = − h − απ . (29)Based on the convergence analysis (Fig. 11), we suppose that these terms are eligible to properly describe thetemperature history after t >
30 (which is equal to 0 . t p = 0 .
01 s). Finally, we can combine these solutions,thus T ( x = 1 , t ) = b − b (the alternating sign originates in cos( jπ )).8.2. Step 2: solving the Guyer-Krumhansl equation.
Here, we repeat the calculations using the same set oftrial and weight functions for the GK model, that is, we solve ∂ t T + ∂ x q = 0 , τ q ∂ t q + q + α∂ x T − κ ∂ xx q = 0 , (30)with q ( t ) = − (cid:16) exp( − C t/t p ) + exp( − C t/t p ) (cid:17) /n, q ( t ) = hT ( x = 1 , t ) , q ( x, t = 0) = 0 , T ( x, t = 0) = 0 . (31)Analogously with the Fourier case, we obtain a set of ODE as follows. • For j = 0, we have ˙ b + ∂ x w = 0 , → ˙ b = − hb + q , (32)which is the same as previously due to a = 0 identically. , R. KOV´ACS Figure 11.
Convergence analysis for the Fourier equation in two different cases on the rear sidetemperature history. The first one shows the adiabatic limit, and the second one presents the casewhen the heat transfer coefficient h is not zero. In both cases, we applied 1, 3 and 20 terms in thespectral decomposition (24), with α = 0 . • For j > a j changes ˙ b j + jπa j = 0 , (33) τ q ˙ a j + (cid:0) κ j π (cid:1) a j = αjπb j + 2 jπ (cid:104) ( hb j − q ) + τ q ( h ˙ b j − ˙ q ) (cid:105) . (34)Consequently, the zeroth term, b ( t ) remains the same with the particular solution being omitted, b ( t ) = Y exp( − ht ) , Y = ( C − C ) / (cid:0) n ( C − h )( C − h ) (cid:1) . (35)However, for b ( t ), the particular solution P ( t ) becomes more important, its initial value influences the temperaturehistory, i.e., P = P ( t = 0) and DP = d t P ( t = 0) appears in the coefficients Z and Z , and the P and DP quantities are important in the evaluation method, too. Thus b ( t ) reads b ( t ) = Z exp ( x t ) + Z exp ( x t ) + P ( t ) , Z = − DP − P x x − x , Z = − P + DP − P x x − x . (36)The exponents x and x depend on the GK parameters τ q and κ , and obtained as the roots of the quadraticequation x j + k j x + k j = 0: x , = x j , | j = 1 , x j , = 12 (cid:16) − k j ± (cid:113) k j − k j (cid:17) , k j = 1 + κ j π τ q + 2 h, k j = αj π τ q + 2 hτ q . (37)Furthermore, the particular solution reads as P j ( t ) = M j exp ( − C t ) + M j exp ( − C t ) , M j = (cid:18) C n − nτ q (cid:19) / (cid:2) k j − k j C + C (cid:3) ,M j = (cid:18) − C n + 2 nτ q (cid:19) / (cid:2) k j − k j C + C (cid:3) . (38) OVEL EVALUATION METHOD FOR NON-FOURIER EFFECTS IN HEAT PULSE EXPERIMENTS 13
Figure 12.
Convergence analysis for the Guyer-Krumhansl equation in two different cases on therear side temperature history. The first one shows the adiabatic limit, and the second one presentsthe case when the heat transfer coefficient h is not zero. In both cases, we applied 1, 3 and 20terms in the spectral decomposition (24), with α = 0 . τ q = 1, and κ = 10 ατ q .Hence P = M + M and DP = − M C − M C appears in b ( t ) with j = 1, too. After obtaining Z and Z , P ( t ) can be neglected since it becomes negligibly small at t > t p . Finally, we formulate the rear side temperaturehistory using b ( t ) and b ( t ) as T ( x = 1 , t ) = b − b = Y exp( − ht ) − Z exp ( x t ) − Z exp ( x t ) , (39)for which Figure 12 shows the convergence property. , R. KOV´ACS References [1] S. Both, B. Cz´el, T. F¨ul¨op, Gy. Gr´of, ´A. Gyenis, R. Kov´acs, P. V´an, and J. Verh´as. Deviation from the Fourier law in room-temperature heat pulse experiments.
Journal of Non-Equilibrium Thermodynamics , 41(1):41–48, 2016.[2] P. V´an, A. Berezovski, T. F¨ul¨op, Gy. Gr´of, R. Kov´acs, ´A. Lovas, and J. Verh´as. Guyer-Krumhansl-type heat conduction at roomtemperature.
EPL , 118(5):50005, 2017. arXiv:1704.00341v1.[3] T. F¨ul¨op, R. Kov´acs, ´A. Lovas, ´A. Rieth, T. Fodor, M. Sz¨ucs, P. V´an, and Gy. Gr´of. Emergence of non-Fourier hierarchies.
Entropy ,20(11):832, 2018. ArXiv: 1808.06858.[4] P. V´an. Theories and heat pulse experiments of non-Fourier heat conduction.
Communications in Applied and Industrial Mathe-matics , 7(2):150–166, 2016.[5] R. A. Guyer and J. A. Krumhansl. Solution of the linearized phonon Boltzmann equation.
Physical Review , 148(2):766–778, 1966.[6] P. V´an. Weakly nonlocal irreversible thermodynamics – the Guyer-Krumhansl and the Cahn-Hilliard equations.
Physic Letters A ,290(1-2):88–92, 2001.[7] P. V´an and T. F¨ul¨op. Universality in heat conduction theory – weakly nonlocal thermodynamics.
Annalen der Physik (Berlin) ,524(8):470–478, 2012.[8] L. Tisza. Transport phenomena in Helium II.
Nature , 141:913, 1938.[9] L. Landau. On the theory of superfluidity of Helium II.
Journal of Physics , 11(1):91–92, 1947.[10] R. A. Guyer and J. A. Krumhansl. Thermal Conductivity, Second Sound, and Phonon Hydrodynamic Phenomena in NonmetallicCrystals.
Physical Review , 148:778–788, 1966.[11] K. Mitra, S. Kumar, A. Vedevarz, and M. K. Moallemi. Experimental evidence of hyperbolic heat conduction in processed meat.
Journal of Heat Transfer , 117(3):568–573, 1995.[12] V. J´ozsa and R. Kov´acs.
Solving Problems in Thermal Engineering: A Toolbox for Engineers . Springer, 2020.[13] R. Kov´acs and P. V´an. Generalized heat conduction in heat pulse experiments.
International Journal of Heat and Mass Transfer ,83:613 – 620, 2015.[14] W. Dreyer and H. Struchtrup. Heat pulse experiments revisited.
Continuum Mechanics and Thermodynamics , 5:3–50, 1993.[15] I. M¨uller and T. Ruggeri.
Rational Extended Thermodynamics . Springer, 1998.[16] D. Y. Tzou. Longitudinal and transverse phonon transport in dielectric crystals.
Journal of Heat Transfer , 136(4):042401, 2014.[17] S. A. Rukolaine. Unphysical effects of the dual-phase-lag model of heat conduction.
International Journal of Heat and MassTransfer , 78:58–63, 2014.[18] S. A. Rukolaine. Unphysical effects of the dual-phase-lag model of heat conduction: higher-order approximations.
InternationalJournal of Thermal Sciences , 113:83–88, 2017.[19] M. Fabrizio, B. Lazzari, and V. Tibullo. Stability and thermodynamic restrictions for a dual-phase-lag thermal model.
Journal ofNon-Equilibrium Thermodynamics , 2017. Published Online:2017/01/10.[20] M. Fabrizio and F. Franchi. Delayed thermal models: stability and thermodynamics.
Journal of Thermal Stresses , 37(2):160–173,2014.[21] B. Ny´ıri. On the entropy current.
Journal of Non-Equilibrium Thermodynamics , 16(2):179–186, 1991.[22] R. Kov´acs. Analytic solution of Guyer-Krumhansl equation for laser flash experiments.
International Journal of Heat and MassTransfer , 127:631–636, 2018.[23] T. F¨ul¨op, R. Kov´acs, and P. V´an. Thermodynamic hierarchies of evolution equations.
Proceedings of the Estonian Academy ofSciences , 64(3):389–395, 2015.[24] K. Zhukovsky. Violation of the maximum principle and negative solutions for pulse propagation in Guyer–Krumhansl model.
International Journal of Heat and Mass Transfer , 98:523–529, 2016.[25] K. V. Zhukovsky. Exact solution of Guyer–Krumhansl type heat equation by operational method.
International Journal of Heatand Mass Transfer , 96:132–144, 2016.[26] ´A. Rieth, R. Kov´acs, and T. F¨ul¨op. Implicit numerical schemes for generalized heat conduction equations.
International Journalof Heat and Mass Transfer , 126:1177 – 1182, 2018. Department of Energy Engineering, Faculty of Mechanical Engineering, BME, Budapest, Hungary Department ofTheoretical Physics, Wigner Research Centre for Physics, Institute for Particle and Nuclear Physics, Budapest, Hun-gary3