A State-of-the-Art Review on Integral Transform Technique in Laser–Material Interaction: Fourier and Non-Fourier Heat Equations

Abstract

Heat equations can estimate the thermal distribution and phase transformation in real-time based on the operating conditions and material properties. Such wonderful features have enabled heat equations in various fields, including laser and electron beam processing. The integral transform technique (ITT) is a powerful general-purpose semi-analytical/numerical method that transforms partial differential equations into a coupled system of ordinary differential equations. Under this category, Fourier and non-Fourier heat equations can be implemented on both equilibrium and non-equilibrium thermo-dynamical processes, including a wide range of processes such as the Two-Temperature Model, ultra-fast laser irradiation, and biological processes. This review article focuses on heat equation models, including Fourier and non-Fourier heat equations. A comparison between Fourier and non-Fourier heat equations and their generalized solutions have been discussed. Various components of heat equations and their implementation in multiple processes have been illustrated. Besides, literature has been collected based on ITT implementation in various materials. Furthermore, a future outlook has been provided for Fourier and non-Fourier heat equations. It was found that the Fourier heat equation is simple to use but involves infinite speed heat propagation in comparison to the non-Fourier heat equation and can be linked with the Two-Temperature Model in a natural way. On the other hand, the non-Fourier heat equation is complex and involves various unknowns compared to the Fourier heat equation. Fourier and Non-Fourier heat equations have proved their reliability in the case of laser–metallic materials, electron beam–biological and –inorganic materials, laser–semiconducting materials, and laser–graphene material interactions. It has been identified that the material properties, electron–phonon relaxation time, and Eigen Values play an essential role in defining the precise results of Fourier and non-Fourier heat equations. In the case of laser–graphene interaction, a restriction has been identified from ITT. When computations are carried out for attosecond pulse durations, the laser wavelength approaches the nucleus-first electron separation distance, resulting in meaningless results.

1. Introduction

Various engineering problems can be modeled using partial differential equations (PDEs) with initial and boundary conditions. For this purpose, numerical approaches, including finite element, finite difference, boundary element, and spectral techniques, are usually applied [1]. The integral transform technique (ITT) is mainly applied to solve linear system problems [2]. A generalized ITT was recently used to resolve numerous linear and non-linear models [3]. All methods reduce the PDEs to a set of ordinary differential equations (ODEs) that can be solved via well-known techniques. For finite difference and element techniques, the field variables are considered in a limited number of points; however, in the boundary elements approach, the number of boundary elements points is minimal. Furthermore, the spectral technique yields a solution in terms of short series numbers. Besides, ITT utilizes “Eigen Function” expansion.

The first law of thermodynamics, also known as the law of conservation of energy, explains that the total energy of an isolated system is constant; energy can be transformed from one form to another but can be neither created nor destroyed [4]. The second law of thermodynamics is about the “quality” of energy. It states that as energy is transferred or transformed, more and more of it is wasted. The second law also states that any isolated system has a natural tendency to degenerate into a more disordered state [4]. It is worth mentioning that all valid heat equations, including Fourier and non-Fourier, obey the first and second laws of thermodynamics.

In this study, the authors have presented the integral transform technique on Fourier and non-Fourier heat equations. We presented the generalized solutions of both techniques in laser–material interaction for the experimentalists to perform simulations before carrying out the actual experiments. We also discussed the various sections of Fourier and non-Fourier heat equations. To further facilitate the experimentalists in this field, we developed friendly user-defined software in MATHEMATICA and provided it within the manuscript for the audience. To the best of the authors’ knowledge, this is the first study that compiles ITT implementation on Fourier and non-Fourier heat equations along with the user-defined codes in the field of experimental physics.

In this review paper, ITT is illustrated as a reliable tool for solving heat equations engineering problems. This technique converts non-linear PDEs to a coupled non-linear ODEs so that they can be solved numerically. ITT on heat equations has been classified into (a) Fourier and (b) non-Fourier heat equations. The generalized solutions of Fourier and non-Fourier heat equations have been derived and discussed. The components of heat equations have been described in detail. Various applications of ITT have been illustrated. Besides, user-defined codes have been provided with explanations in the case of “MATHEMATICA” software programs. After reading this article, one will be able to implement ITT, for Fourier and non-Fourier, in various engineering applications.

1.1. Generalities in Heat Equation: Fourier Heat Equation Formalism and Its Solution

Consider a body at temperature, T(x,y,z), at the position M(x,y,z), and time (t). The temperature rises and falls in various body areas, moving from higher to lower temperatures. A small area (ΔS) around M(x,y,z) can be considered. The section to be heated (ΔQ) is proportional to (ΔtΔS) and normal-derivative $\left(\frac{\partial T}{\partial n}\right)$, as:

$$\mathsf{\Delta }Q=-k\left(x,y,z\right)\mathsf{\Delta }t\mathsf{\Delta }S\frac{\partial T}{\partial n}$$

(1)

Here, k(x,y,z) and n are the thermal conductivity and vector perpendicular to the ΔS and temperature declination orientation, respectively. In terms of thermal conductivity, it has been assumed that the solid-body probe is functioning isotopically. Here, heat passing through a surface area per unit time is denoted as “q.” Now, Equation (1) can be expressed as:

$$q=-k\frac{\partial T}{\partial n}.$$

(2)

To solve the heat equation, consider a random volume (V) covered by a surface (S) and perceive the temperature change within the interval time interval (t₁–t₂). After taking into account the above-defined explanations, one will obtain the following expression:

$$Q_1=-{{\displaystyle \int }}_{t_2}^{t_1}dt{{\displaystyle \int }}_S^{ }{\displaystyle \int } k\left(x,y,z\right)\frac{\partial T}{\partial n}dS.$$

(3)

For a given element volume (ΔV), the heat within the element is given as:

$$\mathsf{\Delta }{Q}_2=\left[T\left(x,y,z,t+\mathsf{\Delta }t\right)-T\left(x,y,z,t\right)\right]c\left(x,y,z\right)\rho \left(x,y,z\right)\mathsf{\Delta }V,$$

(4)

where c(x,y,z) and $\rho \left(x,y,z\right)$ are specific heat and density, respectively, for an element volume, it is essential to change from t₁ to t₂, expressed as:

$$\Delta T=T\left(x,y,z,t_2\right)-T\left(x,y,z,t_1\right)$$

(5)

Therefore

$$Q_2={{\displaystyle \iiint }}_V\left[T\left(x,y,z,t_2\right)-T\left(x,y,z,t_1\right)\right]c\rho dV.$$

(6)

Consequently:

$$T\left(x,y,z,t_2\right)-T\left(x,y,z,t_1\right)={{\displaystyle \int }}_{t_1}^{t_2}\frac{\partial T}{\partial t}dt.$$

(7)

Now, let us assume that there is a heat source within the solid body. In a unit volume per time unit, the source causes the heat absorbance or release, denoted by A(x,y,z,t). The heat absorbance or release resulted from the volume (V) during the time duration $\left(t_1,t_2\right)$ is given as:

$$Q_3={{\displaystyle \int }}_{t_1}^{t_2}dt{{\displaystyle \iiint }}_{V}A\left(x,y,z,t\right)dV.$$

(8)

Noticeably, thermal equilibrium is assumed within the volume V, as:

$$Q_2=Q_1+Q_3.$$

(9)

It infers:

$${{\displaystyle \int }}_{t_2}^{t_1}dt{{\displaystyle \int }}^{ }{{\displaystyle \int }}_V{{\displaystyle \int }}^{ }c\rho \frac{\partial T}{\partial t}dV=-{{\displaystyle \int }}_{t_2}^{t_1}dt{{\displaystyle \iint }}_{S}k\frac{\partial T}{\partial n}dS+{{\displaystyle \int }}_{t_2}^{t_1}dt{{\displaystyle \iiint }}_{V}A\left(x,y,z,t\right)dV.$$

(10)

By utilizing the Gauss–Ostrogradsky formalism, one must achieve:

$${{\displaystyle \int }}_{t_2}^{t_1}dt{{\displaystyle \int }}^{ }{{\displaystyle \int }}_V{{\displaystyle \int }}^{ }[c\rho \frac{\partial T}{\partial t}dV-div(k\cdot grad T)-A\left(x,y,z,t\right)]dV=0.$$

(11)

Equation (11) describes the thermal-field process that evolved when radiations interact with a solid. Exact analytical results can be achieved if a series of generalizations concerning temporal and spatial distributions for a given radiations source and the geometry are made. If the configurations of the radiation source and the geometry to be irradiated become more severe, it is very tough to determine an analytical solution and the thermal distribution, $T\left(x,y,z,t\right)$, expression can be evaluated numerically. A three-dimensional (3D) heat equation is generally expressed as:

$$\rho C\frac{\partial T}{\partial t}=\frac{\partial }{\partial x}\left(k\frac{\partial T}{\partial x}\right)+\frac{\partial }{\partial y}\left(k\frac{\partial T}{\partial y}\right)+\frac{\partial }{\partial z}\left(k\frac{\partial T}{\partial z}\right)+A\left(x,y,z,t\right).$$

(12)

Here, x, y, and z are the universal Cartesian coordinates, T is the temperature, t is the time, k is the thermal conductivity, C is the material’s specific heat, ρ is the material’s density, and A(x,y,z,t) is the heat source. Once the solid body is assumed to be consistent and isotropous, Equation (12) reduces to:

$${\nabla }^{2}T-\frac{1}{\gamma }\frac{\partial T}{\partial t}=-\frac{A\left(x,y,z,t\right)}{k}$$

(13)

where $\gamma$ ($=k/\rho C$) is the thermal diffusivity. For the steady-state condition ($\frac{\partial T}{\partial t}=0$), Equation (13) reduces to:

$${\nabla }^{2}T=-\frac{A\left(x,y,z\right)}{k},$$

(14)

If no external heating source is applied, Equation (13) becomes:

$${\nabla }^{2}T=\frac{1}{\gamma }\frac{\partial T}{\partial t}.$$

(15)

In the case of stationary and no heating source, one achieves the following term:

$${\nabla }^{2}T=0$$

(16)

Equation (10) has an infinite heat-wave speed, and the formalism is known as the “Fourier heat” equation.

$${\nabla }^{2}T=0$$

(17)

In the current study, a cylinder has been considered. It is the reason why the authors used cylindrical coordinates. Equation (12) in terms of cylindrical coordinates (r,z,ϕ) can be written as:

$$\frac{{\partial }^{2}T\left(r,z,\varphi ,t\right)}{\partial r^2}+\frac{1}{r}\frac{\partial T\left(r,z,\varphi ,t\right)}{\partial r}+\frac{{\partial }^{2}T\left(r,z,\varphi ,t\right)}{\partial z^2}+\frac{1}{r^2}\frac{{\partial }^{2}T\left(r,z,\varphi ,t\right)}{\partial {\varphi }^2}-\frac{1}{\gamma }\frac{\partial T\left(r,z,\varphi ,t\right)}{\partial t}=-\frac{A\left(r,z,\varphi ,t\right)}{k}.$$

(18)

The boundary conditions are:

$$K{\frac{\partial T\left(r,z,\varphi ,t\right)}{\partial r}|}_{r=b}+hT\left(b,z,\varphi ,t\right)=0.$$

(19)

$$K{\frac{\partial T\left(r,z,\varphi ,t\right)}{\partial r}|}_{z=0}+hT\left(r,0,\varphi ,t\right)=0.$$

(20)

$$K{\frac{\partial T\left(r,z,\varphi ,t\right)}{\partial r}|}_{z=a}+hT\left(r,a,\varphi ,t\right)=0.$$

(21)

Furthermore, the periodic boundary condition is:

$$T\left(r,0,z,t\right)=T\left(r,2\pi ,z,t\right).$$

(22)

After supposing that laser beam is emitting a beam in the transversal electromagnetic mode (TEM_mn), the thermal distribution variation T(r, z, φ, t) can be determined as [5]:

$$\begin{gathered} T\left(r,\varphi ,z,t\right)=\sum _{m,n}\sum _{i=1}^{\infty }\sum _{l=0}^{\infty }\sum _{j=1}^{\infty }{\widehat{f}}_{2l}\left({\mu }_{il},{\lambda }_j,l\right)g\left({\mu }_{il},{\lambda }_j,t\right)K_r\left({\mu }_{il},r\right)K_{\varphi }\left(2l,\varphi \right)K_z\left({\lambda }_j,z\right) \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+\sum _{m,n}\sum _{i=1}^{\infty }\sum _{l=0}^{\infty }\sum _{j=1}^{\infty }{\widehat{f}}_{2l-1}\left({\mu }_{il},{\lambda }_j,l\right)g\left({\mu }_{il},{\lambda }_j,t\right)K_r\left({\mu }_{il},r\right)K_{\varphi }\left(2l-1,\varphi \right)K_z\left({\lambda }_j,z\right) \end{gathered}$$

(23)

Equation (23) is the solution based on the boundary conditions given in Equations (19)–(22). Furthermore:

$${\stackrel{\wedge }{f}}_{2l}\left({\mu }_{il},{\lambda }_j,l\right)=\frac{1}{k\pi C_{il}{C}_j}{{\displaystyle \int }}_0^a{\alpha }_s{e}^{-{\alpha }_{s}z}{K}_z\left({\lambda }_j,z\right)dz{{\displaystyle \int }}_0^b{{\displaystyle \int }}_0^{2\pi }{I}_{mn}\left(r,\varphi \right)r{K}_r\left({\mu }_{il},r\right)K_{\varphi }\left(2l,\varphi \right)drd\varphi .$$

(24)

$${\stackrel{\wedge }{f}}_{2l-1}\left({\mu }_{il},{\lambda }_j,l\right)=\frac{1}{k\pi C_{il}{C}_j}{{\displaystyle \int }}_0^a{\alpha }_s{e}^{-{\alpha }_{s}z}{K}_z\left({\lambda }_j,z\right)dz{{\displaystyle \int }}_0^b{{\displaystyle \int }}_0^{2\pi }{I}_{mn}\left(r,\varphi \right)r{K}_r\left({\mu }_{il},r\right)K_{\varphi }\left(2l-1,\varphi \right)drd\varphi .$$

(25)

Here, ${\stackrel{\wedge }{f}}_{2l}\left({\mu }_{il},{\lambda }_j,l\right)$ and ${\stackrel{\wedge }{f}}_{2l-1}\left({\mu }_{il},{\lambda }_j,l\right)$ are the source terms of the solution. The term “$A\left(r,z,\varphi ,t\right)$” is the laser beam interaction with the material defined in correspondence with the Lambert–Beer law. It can be identified in the final solution of the heat equation given by ${\stackrel{\wedge }{f}}_{2l}\left({\mu }_{il},{\lambda }_j,l\right)$ and ${\stackrel{\wedge }{f}}_{2l-1}\left({\mu }_{il},{\lambda }_j,l\right)$ source terms of the solution. We also have:

$$I_{mn}\left(x,y\right)=I_{0mn}{\left\{H_m\left(\frac{\sqrt{2x}}{w}\right)H_n\left(\frac{\sqrt{2x}}{w}\right)exp\left[-\left(\frac{x^2+y^2}{w^2}\right)\right]\right\}}^2.$$

(26)

$$g\left({\mu }_{il},{\lambda }_j,t\right)=\frac{1}{\left({\mu }_{il}^2+{\lambda }_j^2\right)\left[1-e^{-{\beta }_{ijl}^{2}t}-\left(1-e^{-{\beta }_{ijl}^2\left(t-t_0\right)}\right)h\left(t-t_0\right)\right]}$$

(27)

$${\beta }_{ilj}^2=\gamma \left({\mu }_{il}^2+{\lambda }_j^2\right).$$

(28)

Here, $g\left({\mu }_{il},{\lambda }_j,t\right)$ and ${\beta }_{ilj}^2$ are the part of the temporal part of the solution given by direct and inverse Laplace transform. Furthermore, h(t − t₀), t_0, and a_s are the Heaviside function, laser–material interaction time, and linear absorption coefficient. Besides, K_r(m_il,r), K_φ(2l,φ), K_φ(2l−1,φ) and K_z(λ_j,z) can be recognized as the “Eigen Functions” for the Eigen Values m_il, 2l, 2l−1 and λ_j. Therefore:

K_r(m_il,r) = J(m_ilr).

(29)

K_φ(2l,φ) = cos(lφ).

(30)

K_φ(2l−1,φ) = sin(lφ).

(31)

K_z(λ_j,z) = cos(λ_jz) + (h/kλ_j)sin(λ_jz).

(32)

Besides, J is the zero^th-order Bessel function, and C_il and C_j are the normalizing coefficients. It is worthy of mentioning that Eigen Values are dependent on the boundary conditions.

1.2. Non-Fourier Heat Equation Formalism and Its Solution

Consider the general form of the non-Fourier transform heat equation for laser–solid interaction. Keeping the standardized notation, the heat source term is expressed as:

$$A\left(x,y,z,t\right)=\{\begin{array}{l}\alpha I_0{e}^{-\alpha x}; 0\le t\le \tau \\ 0; t<0 and t>\tau \end{array}$$

(33)

Here, the Cartesian coordinates have been considered as a sample in parallelepiped shape is irradiated with a laser beam. The following equation shows the non-Fourier heat equation:

$$\frac{{\partial }^2\mathrm{T}}{\partial {\mathrm{x}}^2}+\frac{{\partial }^2\mathrm{T}}{\partial {\mathrm{y}}^2}+\frac{{\partial }^2\mathrm{T}}{\partial {\mathrm{z}}^2}-\frac{1}{\mathsf{\gamma }}\frac{\partial \mathrm{T}}{\partial \mathrm{t}}-\frac{{\mathsf{\tau }}_0}{\mathsf{\gamma }}\frac{{\partial }^2\mathrm{T}}{\partial {\mathrm{t}}^2}=-\frac{\mathrm{A}\left(\mathrm{x},\mathrm{y},\mathrm{z},\mathrm{t}\right)}{\mathrm{K}}\left(1+{\mathsf{\tau }}_0\mathsf{\delta }\left(\mathrm{t}\right)\right).$$

(34)

Here, $\gamma$ is the thermal diffusivity, and K is the thermal conductivity. After implementing ITT, one will obtain the following expression [5]:

$${\mu }_i^2\widehat{T}+{\beta }_l^2\widehat{T}+{\lambda }_j^2\widehat{T}+\frac{1}{\gamma }\frac{\partial \widehat{T}}{\partial t}+\frac{{\tau }_0}{\gamma }\frac{{\partial }^2\widehat{T}}{\partial t^2}=\frac{\widehat{A}\left(x, y, z, t\right)}{K}\left(1+{\tau }_0\delta \left(t\right)\right).$$

(35)

ITT is interpreted as a powerful general-purpose semi-analytical/numerical method. The method transforms partial differential equation models to a coupled system of ordinary differential equations. The implementation of this technique can be identified in-detail from Refs. [6,7]. Here, one should note that a, b and c are the length, width, and height of the substrate, with:

$$f\left({\mu }_i, {\beta }_l,{\lambda }_j, t\right)={{\displaystyle \int }}_0^a{{\displaystyle \int }}_0^b{{\displaystyle \int }}_0^c{\widehat{K}}_x\left({\mu }_i,x\right){\widehat{K}}_y\left({\beta }_l,y\right){\widehat{K}}_{Z }\left({\lambda }_j,z\right)\frac{\widehat{A}\left(x, y, z, t\right)}{K}\left(1+{\tau }_0\delta \left(t\right)\right).$$

(36)

Here, $f\left({\mu }_i, {\beta }_l,{\lambda }_j, t\right)$ and ${\widehat{K}}_x\left({\mu }_i, x\right)$ are the source terms of a solution, ${\widehat{K}}_x\left({\mu }_i,x\right), {\widehat{K}}_y\left({\beta }_l,y\right)$ and ${\widehat{K}}_{Z }\left({\lambda }_j,z\right)$ are the “Eigen Functions.” The first “Eigen Function (${\widehat{K}}_x\left({\mu }_i,x\right),$) is calculated as:

$${\widehat{K}}_x\left({\mu }_i, x\right)=\frac{1}{c_i}{K}_x\left({\mu }_i,x\right).$$

(37)

Furthermore, in the case of the ITT, we have the following expression [5]:

$$\frac{{\partial }^2{\widehat{K}}_x}{\partial x^2}+\lambda {\widehat{K}}_x=0.$$

(38)

The boundary conditions along with the x-axis are given as:

$${\left[k\left(\frac{\partial {\widehat{K}}_x}{\partial x}\right)-h{\widehat{K}}_x\right]}_{x=0}=0.$$

(39)

$${\left[k\left(\frac{\partial {\widehat{K}}_x}{\partial x}\right)+h{\widehat{K}}_x\right]}_{x=a}=0.$$

(40)

The “Eigen Function” can be defined as [5]:

$${\widehat{K}}_x\left({\mu }_i,x\right)=\frac{1}{C_i}\left(\cos \left({\mu }_{i}x\right)+\frac{h}{k{\mu }_i}\sin \left({\mu }_{i}x\right)\right).$$

(41)

where $C_i$ is the normalizing constant, calculated as [5]:

$$C_i={{\displaystyle \int }}_0^a{\widehat{K}}_x^2\left({\mu }_i,x\right)dx.$$

(42)

The expression for “Eigen Values”$\left({\mu }_i\right)$ is given as [5]:

$$2\cot \left({\mu }_{i}a\right)=\frac{{\mu }_{i}k}{h}-\frac{h}{k{\mu }_i}.$$

(43)

Note that:

$${\eta }_{ilj}^2=\frac{1}{{\gamma }^2}-4\left({\mu }_i^2+{\beta }_l^{2 }+{\lambda }_j^2\right)\frac{{\tau }_0}{\gamma }.$$

(44)

It is worthy to mention that the other two “Eigen Values” have also been calculated in the same manner. The temporal part of the generalized solution is obtained via direct and inverse Laplace transform and is expressed as [5]:

$$\begin{array}{ll}g\left({\mu }_i,{\beta }_l,{\lambda }_j,t\right)=& e^{\left(\frac{- 1 - \gamma {\eta }_{ilj}}{2{\tau }_0}\right)t}{C}_1+e^{\left(\frac{- 1+ \gamma {\eta }_{ilj}}{2{\tau }_0}\right)t}{C}_2\\ & -\left(4\frac{{\tau }_{0}A}{\gamma K}\right({\eta }_{ilj}-{\mu }_i^2{e}^{\left(\frac{- 1 - \gamma {\eta }_{ilj}}{2{\tau }_0}\right)t}{\tau }_{0}Unit Step \left[t\right]- {\beta }_l^2{e}^{\left(\frac{- 1 -\gamma {\eta }_{ilj}}{2{\tau }_0}\right)t}{\tau }_{0}Unit Step \left[t\right]\\ & - {\lambda }_j^2{e}^{\left(\frac{- 1- \gamma {\eta }_{ilj}}{2{\tau }_0}\right)t}{\tau }_{0}Unit Step \left[t\right]+ {\mu }_i^2{e}^{\left(\frac{- 1+ \gamma {\eta }_{ilj}}{2{\tau }_0}\right)t}{\tau }_{0}Unit Step \left[t\right]\\ & + {\beta }_l^2{e}^{\left(\frac{- 1+ \gamma {\eta }_{ilj}}{2{\tau }_0}\right)t}{\tau }_{0}Unit Step \left[t\right]\\ & + {\lambda }_j^2{e}^{\left(\frac{- 1+ \gamma {\eta }_{ilj}}{2{\tau }_0}\right)t}{\tau }_{0}Unit Step \left[t\right]\left)\right){\eta }_{ilj}^{-1 }{\left(\left(-\frac{1}{\gamma }+ {\eta }_{ilj}\right)\left(\frac{1}{\gamma }+ {\eta }_{ilj}\right)\right)}^{-1}.\end{array}$$

(45)

The unknown constants, $C_1$ and $C_2$, in the above solution, can be assessed using supplementary boundary conditions, as:

$$T\left(x,y,z,0\right)=0 \Rightarrow C_1+C_2=-\frac{4{\tau }_{0}A}{\gamma K\left({\eta }_{ilj}^2-\frac{1}{{\gamma }^2}\right)}.$$

(46)

$$T\left(\infty ,y,z,t\right)=0 \Rightarrow C_1=-C_2{e}^{\frac{ \gamma {\eta }_{ilj}t }{{\tau }_0}}.$$

(47)

One will obtain:

$$C_2=\frac{-4{\tau }_{0}A}{\gamma K\left(1-e^{\frac{ \gamma {\eta }_{ilj }t}{{\tau }_0}}\right)\left({\eta }_{ilj}^2-\frac{1}{{\gamma }^2}\right)}$$

(48)

$$C_1=\frac{4{\tau }_{0}A{e}^{\frac{\gamma {\eta }_{ilj }t}{{\tau }_0}}}{\gamma K\left(1- e^{\frac{ \gamma {\eta }_{ilj }t}{{\tau }_0}}\right)\left({\eta }_{ilj}^2-\frac{1}{{\gamma }^2}\right)}$$

(49)

Note that the boundary conditions provided in the Fourier section (Equations (19)–(22)) also imply in the non-Fourier section taking into consideration Cartesian Coordinates. Besides, Equations (39) and (40) are the boundary conditions along the x-axis during irradiation. The same type of boundary conditions applies for y- and z-axes. Furthermore, Equations (46) and (47) are the supplementary boundary conditions. It is worth mentioning that T represents the temperature variation, not the absolute temperature. The current study focuses on ITT’s implementation in the non-Fourier heat equation, also known as the Cattaneo–Vernotte equation. The current study’s authors made the first attempt to link the ITT with the non-Fourier heat equation, as presented in Ref. [8].

Fourier and non-Fourier heat equations are valid for equilibrium and non-equilibrium thermodynamics. However, in non-equilibrium thermodynamics, one can implement these solutions for the femtosecond time domain [9]. Besides these techniques, there are various other frameworks, including extended irreversible thermodynamics, generic thermodynamics, internal variables, and rational extended thermodynamics. Jou et al. [10] presented a new formulation of non-equilibrium thermodynamics, known as extended irreversible thermodynamics, which has fueled increasing attention. Evolution equations for these fluxes were obtained starting from a hypothesis and using methods similar to classical irreversible thermodynamics. These equations were reduced to the classical constitutive laws in the limit of slow phenomena. Still, they may also be applied to fast phenomena, such as second sound in solids, ultrasound propagation, or generalized hydrodynamics. In contrast with the classical theory, extended thermodynamics lead to hyperbolic equations with finite propagation speeds for thermal and viscous signals. The results of the macroscopic theory were confirmed by the kinetic theory of gases and non-equilibrium statistical mechanics. The presented theory is particularly useful for studying the thermodynamics of non-equilibrium steady states and systems with long relaxation times, such as viscoelastic media or systems at low temperatures. There is no difficulty in formulating the theory in the relativistic context. Szucs et al. [11] simultaneously applied the methodology of non-equilibrium thermodynamics with internal variables (NET-IV) and the framework of general equation for the non-equilibrium reversible–irreversible coupling (GENERIC). They demonstrated that, in heat conduction theories, entropy current multipliers could be interpreted as relaxed state variables. Fourier’s law and its various extensions–the Maxwell–Cattaneo–Vernotte, Guyer–Krumhansl, Jeffreys type, Ginzburg–Landau (Allen–Cahn) type, and ballistic-diffusive–heat conduction equations were derived in both formulations. They proved that the results might pave the way for microscopic/multiscale understanding of beyond-Fourier heat conduction and open new ways for numerical simulations of heat-conduction problems. Kovacs et al. [12] discussed two different theories: the kinetic-theory-based rational extended thermodynamics (RET) and non-equilibrium thermodynamics with internal variables (NET-IV). It was shown how NET-IV structure is related to RET and how the compatibility between them can be acquired. In another study by Van and Kovacs [13], a comparison of thermodynamic and variational techniques was presented. They found that the second law alone can effectively construct evolution equations for both dissipative and non-dissipative processes.

Furthermore, Jou [14] considered a few conceptual questions on extended thermodynamics to contribute to a higher contact between rational extended thermodynamics and extended irreversible thermodynamics. Both theories take several fluxes as independent variables. Still, they differ in the formalism dealing with the exploitation of the second principle (rational thermodynamics in the first one and classical irreversible thermodynamics in the second one). Rational extended thermodynamics is more restricted in the range of systems to be analyzed, but it is able to obtain a wider number of restrictions and deeper specifications from the second law. By contrast, extended irreversible thermodynamics is more phenomenological, its mathematical formalism is more elementary, but it may deal with a wider diversity of systems, although with less detail. Further comparison and dialogue between both branches of extended thermodynamics would be useful for a fuller deployment and deepening. Besides these two approaches, one should also consider the Hamiltonian approach, formalisms with internal variables, and microscopic approaches based on kinetic theory or non-equilibrium ensemble formalisms [15].

Besides several advantages, the accuracy of the classical heat conduction model, known as Fourier’s law, is highly questioned, dealing with the micro-/nano-systems and biological tissues. In other words, the results obtained from the classical equations deviate from the available experimental data. It means that the continuum heat diffusion equation is insufficient and inappropriate for modeling heat transport in these cases. Consequently, the development of novel models to improve the results of the classical equation while being less computationally expensive and more simple to use is always a topic of interest. There are several techniques for modeling non-Fourier heat conduction. The dual-phase-lag (DPL) model as an accurate modified constitutive equation replacing the Fourier law to simulate the heat transport in special cases such as micro/nanoscales, ultra-fast laser-pulsed processes, living tissues, and carbon nanotube has been trendy [16].

1.3. Comparison between Fourier and Non-Fourier Heat Equations

The second law of thermodynamics is a valuable and universal tool to derive the generalizations of Fourier’s law. In many cases, only linear relations are considered between the thermodynamic fluxes and forces, i.e., the conduction coefficients are independent of the temperature. Kovacs and Rogolino [17] investigated a non-linearity in which the thermal conductivity depends on the temperature linearly. Additionally, that assumption was extended to the relaxation time, which appears in the hyperbolic generalization of Fourier’s law, namely, the Maxwell–Cattaneo–Vernotte (MCV) equation. Although such non-linearity in the Fourier heat equation is well-known, its extension onto the MCV equation is rarely applied. Since these non-linearities have significance from an experimental point of view, an efficient way is needed to solve partial differential equations.

Table 1. Comparison between Fourier and non-Fourier heat equations.

Fourier Heat Equation	References	Non-Fourier Heat Equation	References
This equation is simple to use but involves infinite heat propagation in comparison to the non-Fourier heat equation.	[18]	This equation is complex and involves various unknowns in comparison to the Fourier heat equation.	[19]
It can be applied to both finite and infinite mediums; however, it yields excellent results with a finite medium.	[20]	It can deal with both finite and infinite mediums, but in the current scenario, the non-Fourier heat equation works well in the finite medium. However, it depends on the user’s choice to implement the non-Fourier heat equation in finite and infinite mediums. On the other hand, one should recognize that ITT works well on the finite target only.	[21]
It does not involve electron–phonon relaxation time.	[22]	It takes into account the electron–phonon relaxation time.	[22]
		To achieve an optimum solution, experimental data are always needed for normalizing coefficients.	[23]
The Fourier heat equation cannot be naturally linked with the Two-Temperature Model.	[8,9]	The non-Fourier heat equation can be linked with the Two-Temperature Model more efficiently.	[8,9]
Fourier and non-Fourier heat equations can follow equilibrium and non-equilibrium thermodynamics models for ultra-short laser heating.
The classical heat conduction model’s accuracy is highly questioned, dealing with the micro-/nano-systems and biological tissues. In simple words, the results obtained from the classical equations deviate from the available experimental data. It means that the continuum heat diffusion equation is insufficient and inappropriate for modeling heat transport in these cases.

summarizes a comparison between Fourier and Non-Fourier heat equations.

1.4. Results for Fourier and Non-Fourier Heat Equations via MATHEMATICA Software

In this section, the graphical results achieved via user-defined codes in “MATHEMATICA” software are provided. The sued-defined codes have been provided in Appendix A. The laser beam heats a sample in laser–target interaction, thus providing heat to the system (heat input). In the non-Fourier heat equation and Fourier equation, the heat transfer coefficient is the only parameter that provides heat output from the system, thus following the principle of the first law of thermodynamics. It is why the heat transfer coefficient has been considered in the non-Fourier heat equation, as mentioned in Appendix A. It is necessary to respect the law of conservation of energy. The parameters have been chosen based on the formalism provided for Fourier and non-Fourier heat equations with “SI” units. All the parameters are close to the experimental results [24].

Figure 1 shows the flow chart used to simulate in “MATHEMATICA” software for the Fourier heat equation. For each step, the codes are defined in Appendix A.1.

The above-defined flow chart was used to simulate the Fourier heat equation, and the result is shown in Figure 2. From this figure, it can be analyzed that the laser–material interaction lasts for 3 s, and after this time, the material starts cooling, which causes a sudden declination in temperature.

Figure 3 shows the flow chart used to carry out the simulation in “MATHEMATICA” software for the non-Fourier heat equation. For each step, the codes are defined in Appendix A.2.

Figure 4 shows the results after following the steps defined in Figure 3 when 500 s laser-material interaction time was selected. After material interaction, the thermal distribution intensity falls linearly, resulting in material cooling.

2. Various Applications of Integral Transform Technique (ITT): Fourier and Non-Fourier

In this section, efforts have been made to compile the processes in which ITT has been explored.

2.1. Metallic Materials

Metallic material is a category of materials that contain elemental metals and compounds and alloys [25,26]. There are 86 metals with different characteristic features among the 118 elements in the periodic table, and only a small fraction of these metals has engineering significance [27,28]. Over time, new processes for producing various materials with qualities superior to those of natural materials have been discovered [29,30]. Scientists are working to fully comprehend the links between the structural parts of materials and their qualities [31,32]. Changes in the relative proportions of the micro-constituents can produce drastic changes in terms of quality [33,34]. Phases are distinguished by their distinct crystal forms, compositions, and properties [35].

Metallic materials have been used for various applications. Currently, laser cladding is usually applied in various industrial applications. For this purpose, Oane et al. [36] developed a model in which the phase transition from solid to liquid formation is considered with an absorption coefficient that can highlight the liquid generation during heating. A semi-analytical model was proposed, which considered the melt pool a sphere and solved the heat equation in spherical coordinates. They identified that an increase in the laser scanning speed does not affect the thermal distribution profile significantly. Besides, the simulation results were found in good agreement with the experimental analyses with a 7–15% mean deviation [37]. El-Adawi et al. [38] used ITT for laser–two layers interaction. The thermal distribution was estimated within the thin film and the substrate via mathematical formalism. Furthermore, the front surface temperature was also attained. The simulations were carried for the laser to (a) aluminum–glass, (b) silver–glass, (c) copper–glass, and (d) gold–glass interactions. It was found that the amount of absorbed power determined the crucial time necessary to begin melting. The studied profiles were no longer linear functions of the thermal characteristics of the two-layer system’s materials. The computed values for the crucial time determined whether or not a single laser pulse can cause the damage. The temperature profiles in the case of silver thin film deposited on glass substrate at: (I) t = 5 ns; (II) t = 69.6 ns, and the aluminum thin film deposited on glass substrate at two different exposure times: (I) t = 8.5 ns; (II) t = 67.26 ns, and copper thin film deposited on glass substrate at exposure time t = 23.04 µs are shown in Figure 5a,b.

A thermal interaction model for laser–metal was developed by Oane et al. [39]. The formalism’s primary purpose was to calculate the thermal distribution of electrons and phonons. In laser beam–metal interaction, the Fourier heat equation was solved to provide three-dimensional temperature fields, surface temperature, and steady-state quantum effects temperature. ITT was used to solve the Fourier heat equation. Experimental results for iron irradiation with an Nd: YAG laser beam (=355 nm wavelength) were used to validate the model. The model was able to estimate results with an accuracy of 20% deviation.

Nicarel et al. [40] provided an analytical model for ultrafast thermal processes during femtosecond laser pulse–solid contact. An easy but powerful mathematical formalism was provided to analyze electron temperature’s spatial and temporal profiles under the irradiation with a single Gaussian femtoseconds laser pulse. A gold target was treated with 10¹⁵ W/cm² laser pulses to demonstrate the model’s reliability. The simulations revealed that the electron temperature rises rapidly to saturation in the first 50 fs and that the heat does not travel far beyond the incident laser beam’s focal point. However, when moving from the top surface to the middle of the Au target, the thermal field diminishes by four orders of magnitude. After 1 ps of interaction time, the temperatures of electrons and phonons achieve equilibrium. After 200 fs, when the maximum penetration depth of 7 m is reached, the temperature rises to its highest point. Figure 6a–c show the normalized electron thermal field in arbitrary units on the surface of an Au target versus time and distance at a distance of 4 µm with respect to the sample surface, inside the Au target, and at a depth of 4 μm versus time (from 0 to 1 ps), respectively.

Serban et al. [41] employed a novel approach based on the Fourier heat equation to characterize the laser–metal thermal interaction. The 3D thermal field, surface temperature, and steady-state quantum effects of laser irradiation of metals were measured. An ITT was applied to the Anisimov and Nolte models for this purpose. The laser–gold thin film interaction was the subject of the simulations. The model was able to estimate outcomes with an accuracy of 10⁻² K temperature after comparing with experimental data. In laser–metallic thin film interaction, Oane et al. [42] proposed an analytical technique for estimating thermal fields. A simplified model was presented for this purpose, which included the “global” heat equation. The major goal of this model was to take into account an absorption coefficient that was near to the real one. ITT was used to solve the heat equation for TEM₀₁–silicon copper thin films. The findings revealed that (a) the thin-film absorption coefficient determines the thermal field, and (b) the length of the contact influences thermal distributions.

Buca et al. [19] calculated the electron temperature fluctuation in metals due to the interaction of femtosecond laser pulses with metals. The classical Anisimov’s Two-Temperature Model was extended using the 3D telegraph Zhukovsky equation. The computational plots of electron thermal fields during the initial laser pulse interaction with a gold surface were deduced using this innovative approach. The interaction between the laser pulse and the metal sample during the initial picoseconds is governed by relaxation times and coupling factors over electron thermal conductivities. The higher the electron temperature, the lower the thermal conductivities. In contrast, the lower the electron temperature, the shorter the relaxation time. Figure 7a,b show the spatial–temporal distribution of the electrons’ thermal field generated by one laser pulse on an Au surface, when t = 100 fs and τ = 1 ps versus time (1 ps), and influence of g/K value on temperature intensity at g = 1.05 × 1016 W/m³K, K = 315 W/mK, and laser pulse duration of 100 fs.

Pelin et al. [43] implemented the ITT to solve Fourier and non-Fourier heat equations for the laser–nano-copper-particles interaction. The computations were carried out for 1, 2, and 4 particles-clusters. In comparison to the bulk material in pure form, placing groups or clusters of nano-particles-clusters on top of a layer exposed to irradiation results in a perceptible increase in temperature. Oane et al. [9] presented a new technique by combining the Two-Temperature Model with the Fourier heat equation to estimate the thermal distributions in electrons and phonons. The model was validated in the case of laser–gold, –copper, –silver, and –aluminum interactions. They determined that electron–phonon relaxation time plays a critical role in determining the precision of a given solution. In another study by Oane et al. [44], a semi-analytical model was presented for thermal fields induced in a small cylindrical sample made of tungsten. The model was validated with the experimental results for the tungsten irradiation by an electron beam with an energy of 6 MeV and average power of 62 W from a linear accelerator. The sample had diameter = 10 mm and length = 10 cm. The model predicted results with a mean absolute deviation of 10%.

2.2. Organic and Inorganic Materials

In modern chemistry, organic materials are described as carbon-based substances initially produced from living organisms but today include lab-synthesized counterparts [45]. The majority of them are made up of a few lightweight components, such as hydrogen, carbon, nitrogen, and oxygen [46]. Wood, feathers, leather, and synthetic materials are examples of organic materials [47]. Despite their diversity, they share some common qualities. For instance, prolonged exposure to light or other forms of radiation causes fading, yellowing, or embrittlement in many organic materials, which is caused by the breakdown of the covalent bonding structure shared by many carbon-containing molecules [48]. On the other hand, chemical substances that do not contain carbon are inorganic materials [49]. However, inorganic materials include elementary carbon and carbon compounds such as nitrogen, oxygen, and silicon [50].

Several pieces of research have been carried out on the ITT implementation in organic and inorganic materials. Using a hyperbolic heat conduction model, Vedavarz et al. [51] theoretically studied the transient temperature distributions in laser-irradiated materials. For temperature distributions, exact and limited mathematical solutions were generated, and significant factors are found. To account for the thermal wave’s finite speed, hyperbolic non-Fourier models were devised. The study considered laser surface interactions between two types of materials: (a) biological materials and (b) inorganic solids. They discovered that finite speed effects were substantial in short-pulse, where the laser input time is comparable to the material’s thermal characteristic time. Consequently, the resulting temperature variations were significantly different from those predicted by traditional infinite speed Fourier predictions. The ratio of thermal characteristic length to laser beam width was shown to be the most important parameter. Local temperature maxima or hot spots were detected at the beginning of the irradiation time with values > 0.2 in the range between 0.1 and 3, corresponding to different applications. The dual-phase lag (DPL) bio-heat transport equation was solved analytically by Kumar and Srivastava [52] using a finite-ITT (FITT). Three problems were formulated to show the applicability of the developed formalisms: (a) time-independent boundary conditions (persistent surface temperature), (b) time-dependent boundary conditions (sinusoidal surface heating), and (c) biological tissue irradiated with short-pulses. FITT-based analytical solutions of Fourier and non-Fourier heat conduction equations were linked with numerical solutions to calculate the thermal distribution. Thermal distribution intensity predicted by DPL was in between those achieved using the Fourier and hyperbolic heat conduction equations. When a comparison was made between DPL and Fourier heat models, the hyperbolic heat model generated more obvious wave characters in the anticipated temperature profiles. Talaee et al. [53] used separation of variables and Duhamel integral techniques to solve the three-dimensional hyperbolic heat conduction equation in a cubic media with rectangular cross-section under a pulsed heat flow on the top side. With both steady and pulsed fluxes, the closed-form solution of both Fourier and non-Fourier profiles was introduced. A model was created to simulate the interaction of a cubical tissue with a brief laser pulse. The Fourier and Non-Fourier temperature profiles showed a significant difference in the results. The findings can be used to treat biological tissues with lasers. To examine thermal damage in biological tissues, Zhou et al. [54] developed a thermal wave model of bioheat transmission and a seven-flux model for light propagation, and a rate process equation for tissue damage. They discovered that the thermal damage calculated using the thermal wave bioheat model differed significantly from that calculated using the traditional bioheat model. The assessment of thermal damage to biological tissue may not be reliable if the bioheat non-Fourier effect is not taken into account.

Brasoveanu et al. [55,56] devised a mathematical model to represent the thermal distribution based on the Cattaneo–Vernote formalism’s heat equation, solved by ITT in finite domains. The model was used to calculate the relaxation time, the hottest point, and the peak thermal distribution intensity. Experiments were carried out to record the temperature distribution in granular starch exposed to ionizing radiation to validate the created model. The accelerated electrons irradiated corn starch with a mean energy of 5.5 MeV (=435 s). A temperature sensor was inserted within the starch sample to record the temperature during and after the irradiation. The relaxation time was discovered to play a crucial influence in the cooling of the irradiated sample. The hottest spot was found at a depth of 14 mm in the sample (100 mm total depth). Figure 8 shows the temperature distribution in corn and starch after irradiation.

2.3. Semiconducting Materials

A semiconductor is a type of crystalline substance that is in between a conductor and an insulator in terms of electrical conductivity [57]. Semiconductors are usually used to produce diodes, transistors, and integrated circuits, among other electronic devices [58]. Such devices have found widespread use because of their compactness, reliability, power efficiency, and low cost [59]. They have been used in power devices, optical sensors, and light emitters, including solid-state lasers, as discrete components [60]. They can handle a wide range of current and voltage, and, more importantly, they can be easily integrated into complicated yet easily manufactured microelectronic circuits [61]. They serve communications, signal processing, computing, and control applications in both consumer and industrial markets and will continue to be used in the foreseeable future [62].

The ITT technique has also been explored in the semiconducting materials field. Mahdy et al. [63] identified the impact of three propagated waves: (a) elastic wave, (b) plasma wave, and (c) thermal waves using the hyperbolic generalized two-temperature theory. During the photothermal theory, the governing equations were investigated. They measured the effects of an external magnetic field and a laser pulse. The heat conductivity of semiconductor materials was examined. When the photothermal theory and the thermo-elasticity theory were connected, three different photo-thermo-elasticity models were developed. The principal equations were solved using the ITT in two-dimensional deformation. ITT was demonstrated using the double Fourier and Laplace transforms with appropriate conditions. The complete solutions were obtained numerically by inverting the double transforms with various thermal-elastic-mechanical-plasma boundary conditions. The photo-thermo-elasticity theory was used to compare three different models under an external magnetic field with variable thermal conductivity in the case of silicon material. Thermal memories with a negative thermal conductivity parameter constant, the effect of a magnetic field, and laser pulses in hyperbolic two temperatures were discovered to have a greater impact on the wave propagation of the primary fields. Figure 9 [63] shows the impact of three different cases according to the two-temperature parameter for the main physical fields against the horizontal distance. All numerical results in this category are made under the effect of the magnetic field and the laser pulses in the generalized Green and Lindsay (GL) model. The solid line curves in this category express the first case when the temperature (T) and angle (ϕ) are equal and when the heat supply is absent, which can be named the One-Temperature (OT) Model. The dotted lines curves refer to the Classical Two-Temperature Model (CTT), which is taken when the heat supply is also absent. The dashed lines curves show the general model, which is named the Hyperbolic Two-Temperature (HTT) Model. A clear significant effects in this figure are observed according to the three different cases of the hyperbolic two-temperature theory.

In another investigation by Mahdy et al. [64], the elastic semiconductor medium was irradiated with laser pulses. With the temporal fractional heat order, the laser pulses created vibrations in the medium’s inner structure. The effect of Hall current was observed when a strong external magnetic field was applied. In the case of silicon crystal material, the interaction between the strong magnetic field and the microstructure of the elastic medium was investigated. In the photothermal transport process, photo-excited electrons formed micro-temperature states. In the case of a semiconductor elastic rod, one-dimensional deformation was employed to describe the overlapping process between elastic-magneto-plasma-thermal distribution waves. In the field of microelements, the governing equations were solved utilizing the Laplace integral transform domain. Some external loads were applied to the medium’s outer free surface to achieve the physical quantities under examination. The experiment revealed that variations in laser pulse strength, fractional parameter, and Hall current significantly impacted physical quantities. When the propagated waves curves became similar with the increase in distance, all physical quantities reached equilibrium.

To simulate picosecond and femtosecond laser–silicon interactions, Xu and Wang [65] devised a lattice Boltzmann method (LBT). The LBT was used to calculate temperature fields compared to those derived using the parabolic heat conduction equation (PHCE) and the hyperbolic heat conduction equation (HHCE). Although the HHCE solved the PHCE’s infinite thermal propagation speed, it was unable to be applied to length scales comparable to the mean free path of energy carriers due to the failure of continuum techniques under severe non-equilibrium conditions. The LBT, which considers both impacts, might be used on both short and long-time scales. According to the LBT results, the speed of a thermal wave at the ballistic limit is equal to the speed of sound rather than the figure predicted by the HHCE, which is only true at the diffuse limit. It was further proved that utilizing the temperature gradient to calculate heat flow produces unreasonable results near the thermal wave front, whereas the LBM has no such drawback.

2.4. Graphene Material

A one-atom-thick layer of carbon atoms organized in a hexagonal lattice is known as “graphene” [66]. It is the building block of graphite, but graphene is a fascinating substance in and of itself, with a slew of astounding features that have earned it the moniker “wonder material” on numerous occasions [67]. With one-atom thickness, graphene is the thinnest material known to man and is roughly 200 times stronger than steel [68] at the same thickness scale. Furthermore, graphene is a good conductor of heat and electricity and has remarkable light absorption properties [69].

Scientists are working on ITT implementation in the case of graphene material. For this purpose, ITT was combined with the Anisimov–Nolte Two Temperature Model and the Cattaneo–Vernotte equation by Oane et al. [70]. For laser flash (TEM₀₀)–graphene interaction, they used ITT to develop straightforward formulas for electron and phonon temperatures. Longitudinal optical phonons, transverse optical phonons, out-of-plane acoustic phonons, longitudinal acoustic phonons, transverse acoustic phonons, and out-of-plane optical phonons were all able to generate 3D thermal fields using the model. With a 15% variation, the simulated results were in good accord with the literature [71]. Buca et al. [72] proposed a Multiple-Temperature model to explain a single sheet of graphene irradiation with a flash laser. The Fourier heat equations based on quantum principles, including heat operators, were solved using Zhukovsky’s mathematical technique. Concerning traditional mathematics, simple solutions were inferred. Simple equations for electron and phonon temperatures were established in the case of flash laser treatment of a single layer of graphene. The findings were in good accord with those described in [73,74]. Figure 10a shows the phonon’s temperature, based upon the current model, for the transverse optical (TO) phonons branch. It can be observed that where the laser beam intensity is zero, the temperature distribution is zero. This shows the perfect coupling between Zhukovsky’s formalism and Multiple-Temperature Model (MTM). However, thermal effects can be clearly observed beyond the laser beam spot size (=1 µm), as a very short laser–graphene interaction time has been chosen (of ps order). The same phenomena can be analyzed in-case of Figure 10b,c related to longitudinal optical phonons (LO) and longitudinal acoustic phonons (LA) braches, respectively. As the simulation time has been increased, from ps to ns, the heat-wave finds enough time to transfer from the irradiated spot to the rest of the target. This result has been presented in Figure 10d in the case of transverse acoustic phonons (TA) branch.

However, the model presented a few drawbacks, such as:

• The duration of the irradiation time/pulse should be greater than 1 fs. When simulations are for attosecond pulse durations, the laser wavelength approaches the nucleus-first electron separation distance, resulting in meaningless results.
• The target size is limited to a range of 20–100 nm. The Fourier law collapses when the particle size is less than 20 nm, and reliable results are impossible to obtain.
• The model was restricted to single-layer graphene alone. The study of multilayer graphene becomes complicated to calculate since the absorption law and heat-transfer coefficients must be evaluated in real-time for each layer, necessitating the use of quantum-field theory in solid-state treatment.

Chan et al. [75] used the Fourier transformation approach to determine the orientation of hexagonal graphene domains on a Cu substrate. They devised a hexagon function to characterize hexagonal graphene’s diffraction pattern. In the frequency domain, hexagonal graphene domains produced on Cu (111) had an average orientation value of around 3°. The optical and electrical properties of a large-area graphene film (2 cm²) were tested for transparent conducting electrode applications. The findings showed that graphene produced on Cu (111) outperformed graphene developed on polycrystalline Cu.

3. Criticism on Dual-Phase-Lag (DPL) Heat Conduction and Guyer-Rumhansl Equations

The DPL and also the Guyer-Krumhansl equations are criticized from several different points of view. Negative temperatures Rukolaine [76,77] obtained an analytical solution for the DPL equation assuming a Gaussian initial condition. According to these calculations, the solutions present an unphysical behavior of temperature history, and it goes into the negative domain. Zhukovsky [78,79] achieved a similar conclusion for the GK equation. Wang et al. [80] tested the Thermo mass and DPL among other different heat conduction models by calculating Taitel’s problem [81]. Inconsistent, unphysical behavior is shown, and the temperature achieves the negative domain again. It is an old problem, also mentioned several times for the Maxwell-Cattaneo-Vernotte equation.

Time shift paradox: the delayed form of the DPL model directly contradicts the requirement of the homogeneity of time. One can transform the equation to a single-phase-lag equation by shifting the zero point of t, the starting point of the time measurement. The proper mathematical representation of non-relativistic time considering all expected properties is a one-dimensional affine space, and with that model cannot have two relaxation times; only their difference plays a role [82].

Mathematical question: in the literature, the DPL-type constitutive equations are analyzed to prove the uniqueness and well-posedness of a process driven by such a constitutive equation. It is found by Fabrizio et al. [83,84,85] that there are mathematical conditions beyond the physical ones to obtain an exponentially stable equilibrium solution for the DPL equation. Such a condition requires a negative time delay (called a retarded effect) between the heat flux and temperature gradient. It is important because Tzou in Ref. [86] directly interprets both cases with the cause–effect concept, i.e., the quantity with the higher relaxation time is the effect caused by the other one. As Fabrizio states [83], the DPL model can be rewritten with the time delay difference, leading to a single-phase-lag model, but the temperature gradient has a relaxation time < 0. The opposite case is mathematically ill-posed, which enlightens the validity of the Maxwell–Cattaneo–Vernotte (MCV) equation but excludes equations based on arbitrary Taylor series expansion. Quintanilla et al. [87,88] obtain the same conclusion regarding the relaxation times and the ill-posedness. However, exponential stability seems to be a too-strict requirement, and the asymptotic stability of homogeneous equilibrium is satisfied only with non-negative coefficients. The condition of asymptotic stability is also reasonable from a physical point of view [89,90,91].

Second law: Fabrizio and co-workers also checked thermodynamical restriction. It was tested by using Clausius–Duhem inequality on cyclic histories. The conditions of exponential stability for the DPL equation turned out to be too strict [84]. On the other hand, the experimental pieces of evidence mostly show that this difference is negative. In the work of Liu and Chen [92], the DPL equation is fitted to experimental results, and in heterogeneous materials, there is a similar situation [93,94]. It draws attention to the practical aspects, which are emphasized further in the next section.

4. Future Outlook

Fourier and Non-Fourier heat equations’ solutions using the integral transform technique (ITT) have proved their viability in the interactions of laser–various materials. In the authors’ opinion, there are various domains where these equations can be still explored, as:

• Fourier and non-Fourier heat equations can be linked with Two- and Multi-Temperature Models to estimate the thermal distributions within the sample. The Two-Temperature Model can yield thermal distribution in electrons and phonons, while the Multi-Temperature Model can give thermal distributions in longitudinal optical phonons, transverse optical phonons, out-of-plane acoustic phonons, longitudinal acoustic phonons, transverse acoustic phonons, and out-of-plane optical phonons.
• These heat equations have not been explored in the case of laser–ceramic interaction. Ceramics are porous materials having voids and pores of various dimensions [95]. Furthermore, ceramics behave as insulators and present complex rheological properties [96]. Hence, there is a need to modify the Fourier and non-Fourier heat equations by taking into considerations all the aspects explained above. It will open a new area for laser–ceramic interaction heat equations’ solution using ITT. Furthermore, it is a long-lasting task to understand heat conduction phenomena beyond Fourier. Besides the low-temperature experiments on extremely pure crystals, it has been turned out recently that heterogeneous materials with macro-scale size can also show thermal effects that the Fourier equation cannot model. It is known as over-diffusive propagation, different from low-temperature observations, and is found in numerous samples made from metal foam, rocks, and composites. The measured temperature history is similar to what Fourier’s law predicts, but the usual evaluation cannot provide reliable thermal parameters. Feher et al. [97] reported experiments on several rock types, each having multiple samples with different thicknesses. They showed that the size-dependent thermal behavior could occur for both Fourier and non-Fourier situations. Moreover, based on the experimental data, they found an empirical relation between the Fourier and non-Fourier parameters, which could be helpful in later experiments to develop a more robust and reliable evaluation procedure.
• Another domain is to develop a solution for laser–metal matrix composites (MMCs) interaction [98]. MMCs are the composition of metals and ceramics mixed to achieve a combination of metallic and ceramic properties [99]. In this area, no Fourier and non-Fourier models are dealing with laser–MMCs interaction [100]. From the authors’ point of view, this problem can be solved by separately solving the heat equations for metals and ceramics. Following, “superposition” must be implemented based on the percentage of metals and ceramics within MMCs to attain the final solution.
• For the automation of Fourier and non-Fourier heat equations, these equations can be linked with artificial neural networking (ANN). A series of computations can be carried out based on inputs and corresponding outputs [101]. These values can be used to train ANN and further used to estimate output based on desired inputs [102]. Recently, Mahmood et al. used the same approach to link ANN with analytical modeling for process automation [103].
• This method has been explored for laser–metal, Two-Temperature Model for laser–metal interaction, ultra-short laser–material interaction, and electron beam–organic and –inorganic interaction. However, there are various applications that have not been studied in detail. For instance, the laser–graphene interaction, laser–ceramic reinforced metal matrix interaction, laser–nanoparticles interaction and laser–ceramic interaction are a few subjects of interest in the near future.

The authors believe that the Fourier and non-Fearer heat equations will exponentially gain attention from scientists and researchers after exploring the areas mentioned above.

5. Conclusions

Fourier and non-Fourier heat equations can estimate the thermal distribution and phase transformation based on the provided operating conditions. For this purpose, this review article discusses the heat equations, including Fourier and non-Fourier heat equations. An evaluation is presented between Fourier and non-Fourier heat equations. Besides, their generalized solutions have been presented. Various components of heat equations and their implementation in the various processes have been illustrated. For solving these equations, “MATHEMATICA” software codes have been provided for the scientific community. Furthermore, various applications of Fourier and non-Fourier techniques have been provided. Based on the current study, the following conclusions have been deduced:

• The Fourier heat equation is simple to use but involves infinite heat propagation compared to the non-Fourier heat equation. It can be applied to both finite and infinite mediums; however, it yields excellent results with a finite medium. It does not involve electron–phonon relaxation time. The Fourier heat equation cannot be linked with the Two-Temperature Model in a natural way.
• The non-Fourier heat equation is complex and involves various unknowns in comparison to the Fourier heat equation. It can deal with both finite and infinite mediums, but in the current scenario, the non-Fourier heat equation works well in the finite medium. However, it depends on the user’s choice to implement the non-Fourier heat equation in finite and infinite mediums. On the other hand, one should recognize that ITT works well only on the finite target. It takes into account the electron–phonon relaxation time. To achieve an optimum solution, experimental data are always needed for normalizing coefficients. The non-Fourier heat equation can be linked with the Two-Temperature Model in a more efficient way.
• Fourier and non-Fourier heat equations can follow equilibrium and non-equilibrium thermodynamics models for ultra-short laser heating. The classical heat conduction model’s accuracy is highly questioned, dealing with the micro-/nano-systems and biological tissues. In simple words, the results obtained from the classical equations deviate from the available experimental data. It means that the continuum heat diffusion equation is insufficient and inappropriate for modeling heat transport in these cases.
• Fourier and Non-Fourier heat equations have proved their reliability in the case of laser–metallic materials, electron beam–biological and –inorganic materials, laser–semiconducting materials, and laser–graphene material interactions. The presented simulation results showed a strong correspondence with experimental results. It is worth mentioning that the accuracy and transient behavior determine the accuracy of Fourier and non-Fourier heat equations. Especially in the case of biological systems, the Fourier model is not necessarily far from the measurements, but the transient behavior could be significantly different, and therefore the prediction is questionable.
• Through the literature survey and the authors’ research experience, it has been identified that the material properties play an essential role in defining the precise results of Fourier and non-Fourier heat equations.
• In the case of laser–graphene interaction, there is a certain restriction of the integral transform technique. The duration of the irradiation time/pulse should be greater than 1 fs. When computations are for attosecond pulse durations, the laser wavelength approaches the nucleus-first electron separation distance, resulting in meaningless results. The study of multilayer graphene becomes complicated to calculate since the absorption law and heat-transfer coefficients must be evaluated in real-time for each layer, necessitating the use of quantum-field theory in solid-state treatment.