Variational Approach to 2D and 3D Heat Conduction Modeling

Slavko Đurić^1*, Ivan Aranđelović², Milan Milotić¹

¹Faculty of Traffic Engineering, University of East Sarajevo, Doboj, Bosnia and Herzegovina.
²Faculty of Mechanical Engineering, University of Belgrade, Belgrade, Republic of Serbia.

DOI: 10.4236/jamp.2024.124085                

Abstract

The paper proposes an approximate solution to the classical (parabolic) multidimensional 2D and 3D heat conduction equation for a 5 × 5 cm aluminium plate and a 5 × 5 × 5 cm aluminum cube. An approximate solution of the generalized (hyperbolic) 2D and 3D equation for the considered plate and cube is also proposed. Approximate solutions were obtained by applying calculus of variations and Euler-Lagrange equations. In order to verify the correctness of the proposed approximate solutions, they were compared with the exact solutions of parabolic and hyperbolic equations. The paper also presents the research on the influence of time parameters $\tau$ as well as the relaxation times ${\tau }^{\ast }$ to the variation of the profile of the temperature field for the considered aluminum plate and cube.

Keywords

Classical Equation of Heat Conduction, Generalized Equation of Heat Conduction, Calculus of Variations, Approximate Solution

1. Introduction

In many engineering applications (thermal engineering, process engineering, aero-engineering, etc.), the classical Fourier heat conduction equation is considered an efficient method for solving heat conduction problems. The equation that describes the distribution of the temperature field is the diffusion equation. Due to the importance of diffusion in engineering practice, solutions to the 1D, 2D and 3D diffusion equations have been studied by a large number of researchers [1] - [7] . One-dimensional heat conduction problems with Dirichlet boundary conditions can be directly solved by the Fourier method of separation of variables. The application of this method in engineering practice to solving 2D and 3D heat conduction problems is quite complex. Application of some other classical analytical methods to the solution of 2D and 3D heat conduction problems such as Laplace transform, Duhamel’s theorem, Green’s function and the Riemann method is not very suitable for engineering practice. A method that is widely used in practice is a numerical method, such as the finite difference method [8] [9] .

Very fast material heating processes, e.g. during the absorption of energy coming from ultra-short laser pulses cannot be satisfactorily explained by the Fourier equation of heat conduction. It predicts an infinite rate of expansion of the thermal disturbance, which is physically inadmissible and contradicts the existing recent theories that treat the phenomenon of heat transfer. In order to resolve the paradox of thermal disturbance propagation in the classical Cattaneo diffusion theory [10] [11] , Vernotte [12] introduces the relaxation time into the Fourier law of heat conduction (CV model).

The CV model (generalized Fourier equation of heat conduction) describes hyperbolic partial differential equations (telegraph equation [13] ) and assumes a finite rate of propagation of the thermal disturbance. Based on these assumptions, in the considered theories, thermal disturbance in a solid medium behaves like a wave (the so-called second sound).

In recent years, several non-classical theories in the theory of thermoelasticity have been presented [14] - [21] . Those theories represent also a modified Fourier equation of heat conduction, including the hyperbolic form of the partial differential equation of heat transfer and the finite velocity of propagation of the thermal disturbance.

The main goal of the exhibited manuscript is to underline the possibility of variational description of 2D and 3D heat conduction through solid bodies. This paper develops a simplified analytical method (approximate solution) for 2D and 3D parabolic and hyperbolic heat conduction in aluminum plate and aluminum cube with Dirichlet boundary conditions using calculus of variations and the Euler-Lagrange equation. In the end, a numerical example is given, and the accuracy of the obtained solution using the approximate solution was checked by comparing it with the exact solution of 2D and 3D heat conduction presented in the paper.

2. Mathematical Formulation

Generalized heat conduction equation:

${\tau }^{\ast }\cdot \frac{{\partial }^{2}T}{\partial {\tau }^2}+\frac{\partial T}{\partial \tau }=k\cdot \left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}+\frac{{\partial }^{2}T}{\partial z^2}\right),$ (1)

can be derived from the Langrangian:

$L=\left[\frac{{\tau }^{\ast }}{2}{\left(\frac{\partial T}{\partial \tau }\right)}^2-\frac{k}{2}\cdot \left({\left(\frac{\partial T}{\partial x}\right)}^2+{\left(\frac{\partial T}{\partial y}\right)}^2+{\left(\frac{\partial T}{\partial z}\right)}^2\right)\right]\cdot {\text{e}}^{\frac{\tau }{{\tau }^{\ast }}}$ , (2)

using the Euler-Lagrange differential equation:

$\frac{\partial L}{\partial T}-\frac{\partial }{\partial \tau }\left(\frac{\partial L}{\partial {T^{\prime }}_{\tau }}\right)-\frac{\partial }{\partial x}\left(\frac{\partial L}{\partial {T^{\prime }}_x}\right)-\frac{\partial }{\partial y}\left(\frac{\partial L}{\partial {T^{\prime }}_y}\right)-\frac{\partial }{\partial z}\left(\frac{\partial L}{\partial {T^{\prime }}_z}\right)=0,$ (3)

and the following:

$\left(x,y,z\right)\in V$ ( $V\subseteq R^3$ ) are coordinates of an arbitrary point in the volume V;

$\tau$ —time;

${\tau }^{\ast }$ —relaxation time;

$T\left(x,y,z,\tau \right)$ —temperature value (temperature field) in volume V.

Equation (1) is also valid and represents the variational task of finding the stationary value of the integral:

$F={\displaystyle \underset{0}{\overset{\tau }{\int }}\text{d}\tau {\displaystyle \underset{V}{\iiint }L\text{d}V}},$ (4)

where the L is Lagrangian under (2).

2.1. Two-Dimensional Heat Conduction

The two-dimensional generalized (hyperbolic) heat conduction equation is given:

${\tau }^{\ast }\cdot \frac{{\partial }^{2}T}{\partial {\tau }^2}+\frac{\partial T}{\partial \tau }=k\cdot \left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}\right),\left(x,y\right)\in D,\tau >0,$ (5)

where the area is $D=\left\{\left(x,y\right):0\le x\le a,0\le y\le b,a,b\in R_{+}\right\}$ .

The following conditions will be assumed for Equation (5):

Initial condition:

$T\left(x,y,0\right)=f\left(x,y\right)=x\left(a-x\right)y\left(b-y\right),T_{\tau }\left(x,y,0\right)=0,\left(x,y\right)\in D$ , (6)

Boundary condition:

$T\left(0,y,\tau \right)=T\left(a,y,\tau \right)=0$ , (7)

$T\left(x,0,\tau \right)=T\left(x,b,\tau \right)=0$ ,

$0\le x\le a$ , $0\le y\le b$ , $\tau >0$ .

Equation (5) with the given initial (6) and boundary conditions (7) will be solved by the Fourier method of separation of variables, i.e. the analytical solution will be sought in the form:

$T\left(x,y,\tau \right)=X\left(x\right)\cdot Y\left(y\right)\cdot \Psi \left(\tau \right)$ , (8)

and by substituting into Equation (5), we get:

$\frac{X\text{'}\text{'}}{X}+\frac{Y\text{'}\text{'}}{Y}=\frac{{\tau }^{\ast }}{k}\frac{\Psi \text{'}\text{'}}{\psi }+\frac{1}{k}\frac{\Psi \text{'}}{\Psi },$

i.e.

$\frac{X\text{'}\text{'}}{X}=-{\alpha }^2$ , (9)

$\frac{Y\text{'}\text{'}}{Y}=-{\beta }^2$ , (10)

$\frac{{\tau }^{\ast }}{k}\frac{\Psi \text{'}\text{'}}{\psi }+\frac{1}{k}\frac{\Psi \text{'}}{\Psi }=-{\gamma }^2$ , (11)

whereas ${\gamma }^2={\alpha }^2+{\beta }^2$ .

From Equations (9)-(11), three equations are obtained:

$X\text{'}\text{'}+{\alpha }^{2}X=0$ , (12)

$Y\text{'}\text{'}+{\beta }^{2}Y=0$ , (13)

${\tau }^{\ast }\Psi \text{'}\text{'}+\Psi \text{'}+k{\gamma }^2\Psi =0$ . (14)

The solution of Equation (12) is:

$X\left(x\right)=C_1\cos \left(\alpha x\right)+C_2\sin \left(\alpha x\right)$ ,

and using the boundary conditions (7), the value of the constant $C_1=0$ , it follows that:

$X\left(x\right)=C_2\sin \left(\alpha x\right)$ ,

i.e. $C_2\sin \left(\alpha a\right)=0$ implying:

${\alpha }_m=\frac{m\pi }{a}\left(m=1,2,\cdots \right).$ (15)

In a similar way, Equation (13) is obtained:

$Y\left(y\right)=C_3\cos \left(\beta y\right)+C_4\sin \left(\beta y\right)$ ,

and using the boundary conditions (7), the value of the constant $C_3=0$ implying:

$Y\left(y\right)=C_4\sin \left(\beta y\right),$

i.e. $C_4\sin \left(\beta b\right)=0$ implying:

${\beta }_n=\frac{n\pi }{b}\left(n=1,2,\cdots \right)$ , (16)

${\gamma }_{mn}^2={\alpha }_m^2+{\beta }_n^2$ . (17)

Equation (14) is a second-order differential equation with constant coefficients, so by solving it in a known way, the solution is given:

${\Psi }_{mn}\left(\tau \right)=A_{mn}{\text{e}}^{r_{1,mn}\tau }+B_{mn}{\text{e}}^{r_{2,mn}\tau }$ , (18)

with:

$r_{1,mn}=\frac{-1+\sqrt{1-4k{\gamma }_{mn}^2{\tau }^{\ast }}}{2{\tau }^{\ast }}$ ,

$r_{2,mn}=\frac{-1-\sqrt{1-4k{\gamma }_{mn}^2{\tau }^{\ast }}}{2{\tau }^{\ast }}$ ,

where the relaxation time belongs to the interval ${\tau }^{\ast }\in \left(0,\frac{1}{4k{\gamma }_{mn}^2}\right]$ .

Where every function $T_{mn}\left(x,y,\tau \right)$ is defined by:

$T_{mn}\left(x,y,\tau \right)=\sin \left({\alpha }_{m}x\right)\cdot \sin \left({\beta }_{n}y\right)\cdot \left(A_{mn}\cdot {\text{e}}^{r_{1,mn}\tau }+B_{mn}\cdot {\text{e}}^{r_{2,mn}\tau }\right)$ , (19)

satisfies equalities (5) and (6) and based on the principle of linear superposition and function $T(x,y,\tau )$ , it is defined as follows:

$T\left(x,y,\tau \right)={\displaystyle \underset{m=1}{\overset{\infty }{\sum }}{\displaystyle \underset{n=1}{\overset{\infty }{\sum }}\sin \left({\alpha }_{m}x\right)\cdot \sin \left({\beta }_{n}y\right)\cdot \left(A_{mn}\cdot {\text{e}}^{r_{1,mn}\tau }+B_{mn}\cdot {\text{e}}^{r_{2,mn}\tau }\right)}}$ . (20)

It represents the solution to problems (5)-(7).

Constants $A_{mn}$ i $B_{mn}$ will be determined from the conditions $T_{\tau }\left(x,y,0\right)=0$ with $A_{mn}=C_{mn}\cdot r_{2,mn}$ and $B_{mn}=-C_{mn}\cdot r_{1,mn}$ .

Equation (20) now has the form:

$T\left(x,y,\tau \right)={\displaystyle \underset{m=1}{\overset{\infty }{\sum }}{\displaystyle \underset{n=1}{\overset{\infty }{\sum }}{C}_{mn}\left(r_{2,mn}\cdot {\text{e}}^{r_{1,mn}}-r_{1,mn}\cdot {\text{e}}^{r_{2,mn}}\right)\sin \left({\alpha }_{m}x\right)\cdot \sin \left({\beta }_{n}y\right)}}$ . (21)

From the initial condition:

$T\left(x,y,0\right)={\displaystyle \underset{m=1}{\overset{\infty }{\sum }}{\displaystyle \underset{n=1}{\overset{\infty }{\sum }}{C}_{mn}\left(r_{2,mn}-r_{1,mn}\right)}}\cdot \sin \left({\alpha }_{m}x\right)\cdot \sin \left({\beta }_{n}y\right)=f\left(x,y\right)$ , (22)

and putting that it is $E_{mn}=C_{mn}\left(r_{2,mn}-r_{1,mn}\right)$ . Equation (22) is equivalent to the equation:

$T\left(x,y,0\right)={\displaystyle \underset{m=1}{\overset{\infty }{\sum }}{\displaystyle \underset{n=1}{\overset{\infty }{\sum }}{E}_{mn}}}\cdot \sin \left({\alpha }_{m}x\right)\cdot \sin \left({\beta }_{n}y\right)=f\left(x,y\right)$ , (23)

whereby $E_{mn}$ is determined in the expression:

$E_{mn}=\frac{4}{ab}{\displaystyle \underset{0}{\overset{a}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}f\left(x,y\right)\cdot \sin \left({\alpha }_{m}x\right)\cdot \sin \left({\beta }_{n}y\right)\text{d}x\text{d}y}}$ , (24)

and according to the given initial conditions $T\left(x,y,0\right)=f\left(x,y\right)=x\left(a-x\right)y\left(b-y\right),0\le x\le a,0\le y\le b$ ,

we get:

$E_{mn}=\frac{64a^2{b}^2}{{\left(2m-1\right)}^3\cdot {\left(2n-1\right)}^3\cdot {\pi }^6},m,n,=1,2,3,\cdots$

so the solution to the given problems (5)-(7) is given by:

$T\left(x,y,\tau \right)={\displaystyle \underset{m=1}{\overset{\infty }{\sum }}{\displaystyle \underset{n=1}{\overset{\infty }{\sum }}\frac{E_{mn}}{r_{2,mn}-r_{1,mn}}\cdot \left(r_{2,mn}\cdot {\text{e}}^{r_{1,mn}\tau }-r_{1,mn}\cdot {\text{e}}^{r_{2,mn}\tau }\right)\cdot \sin \left({\alpha }_{m}x\right)\cdot \sin \left({\beta }_{n}y\right)}}$ . (25)

2.2. Determining the Approximate Solution of the Two-Dimensional Heat Conduction Equation

Let us now find an approximate solution to the generalized heat conduction equation (Equation (5)). An approximate solution of the form will be assumed:

$T\left(x,y,\tau \right)=f\left(x,y\right)\cdot \phi \left(\tau \right)=x\left(a-x\right)y\left(b-y\right)\cdot \phi \left(\tau \right)$ , (26)

where $\phi \left(\tau \right)$ is anunknown function that satisfies the condition:

$\phi \left(0\right)=1$ . (27)

If we start from the functional:

$\begin{array}{c}F\left(T\left(x,y,\tau \right)\right)={\displaystyle \underset{0}{\overset{\tau }{\int }}\text{d}\tau }{\displaystyle \underset{D}{\iint }F\left(x,y,\tau ,T\left(x,y,\tau \right),{T^{\prime }}_x,{T^{\prime }}_{\tau }\right)\text{d}x\text{d}y}\\ ={\displaystyle \underset{0}{\overset{\tau }{\int }}{\displaystyle \underset{0}{\overset{a}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}\left[\frac{{\tau }^{\ast }}{2}\cdot {\left(\frac{\partial T}{\partial \tau }\right)}^2-\frac{k}{2}\cdot \left({\left(\frac{\partial T}{\partial x}\right)}^2+{\left(\frac{\partial T}{\partial y}\right)}^2\right)\right]}}}\cdot {\text{e}}^{\tau /{\tau }^{\ast }}\text{d}\tau \text{d}x\text{d}y,\end{array}$ (28)

by inserting the trial solution (26) into Equation (28) and solving the double integral, we get:

$F=\frac{a^3{b}^3}{180}\cdot {\displaystyle \underset{0}{\overset{\tau }{\int }}\left(\frac{a^2{b}^2{\tau }^{\ast }}{10}\cdot {\stackrel{\dot{}}{\phi }}^2-\left(a^{2}k+b^{2}k\right){\phi }^2\right)\cdot {\text{e}}^{\tau /{\tau }^{\ast }}\text{d}\tau }$ . (29)

If we now introduce the Lagrange function:

$L\left(\phi ,\stackrel{\dot{}}{\phi }\right)=\left(\frac{a^2{b}^2{\tau }^{\ast }}{10}\cdot {\stackrel{\dot{}}{\phi }}^2-\left(a^{2}k+b^{2}k\right)\cdot {\phi }^2\right)\cdot {\text{e}}^{\tau /{\tau }^{\ast }}$ , (30)

and using the Euler-Lagrange equation:

$\frac{\partial L}{\partial \phi }-\frac{\partial }{\partial \tau }\left(\frac{\partial L}{\partial \stackrel{\dot{}}{\phi }}\right)=0,$ (31)

and dividing with $e^{\tau /\tau *}$ , the following differential equation is obtained:

${\tau }^{\ast }\cdot \ddot{\phi }\left(\tau \right)+\stackrel{\dot{}}{\phi }\left(\tau \right)+\frac{10k\left(a^2+b^2\right)}{a^2{b}^2}\cdot \phi \left(\tau \right)=0$ . (32)

After a short calculation, the solution to Equation (32) is obtained:

$\phi \left(\tau \right)=C_1\cdot {\text{e}}^{r_1\tau }+C_2\cdot {\text{e}}^{r_2\tau }$ , (33)

with:

$r_1=\frac{-1+\sqrt{1-\frac{40k\left(a^2+b^2\right){\tau }^{\ast }}{a^2{b}^2}}}{2{\tau }^{\ast }}$ , $r_2=\frac{-1-\sqrt{1-\frac{40k\left(a^2+b^2\right){\tau }^{\ast }}{a^2{b}^2}}}{2{\tau }^{\ast }}$ and ${\tau }^{\ast }\in \left(0,\frac{a^2{b}^2}{40k\left(a^2+b^2\right)}\right]$ .

Constants $C_1$ i $C_2$ are determined from the conditions of the conditions $\phi \left(0\right)=1$ and $T_{\tau }\left(x,y,0\right)=0$ . The first condition implies:

$C_1+C_2=1$ , (34)

and the second condition implies:

$C_1\cdot r_1+C_2\cdot r_2=0$ . (35)

The constants are determined from the system of Equations (34) and (35). $C_1$ i $C_2$ :

$C_1=\frac{r_2}{r_2-r_1}$ , (36)

$C_2=\frac{r_1}{r_1-r_2}$ , (37)

so an unknown function $\phi \left(\tau \right)$ is determined by:

$\phi \left(\tau \right)=\frac{r_2}{r_2-r_1}\cdot {\text{e}}^{r_1\tau }+\frac{r_1}{r_1-r_2}\cdot {\text{e}}^{r_2\tau }$ . (38)

Now, according to (26), the approximate solution of Equation (5) is given by:

$T\left(x,y,\tau \right)=x\left(a-x\right)y\left(b-y\right)\cdot \left(\frac{r_2}{r_2-r_1}\cdot {\text{e}}^{r_1\tau }+\frac{r_1}{r_1-r_2}\cdot {\text{e}}^{r_2\tau }\right)$ . (39)

If the boundary condition ${\tau }^{\ast }\to 0$ is applied to Equation (32) or Equation (39), an approximate solution to the classical two-dimensional heat conduction equation is obtained:

$\frac{\partial T}{\partial \tau }=k\cdot \left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}\right)$ , (40)

which is given by:

$T\left(x,y,\tau \right)=x\left(a-x\right)y\left(b-y\right)\cdot {\text{e}}^{-\frac{10k\left(a^2+b^2\right)}{a^2{b}^2}\cdot \tau }$ . (41)

2.3. Three-Dimensional Heat Conduction

Consider a three-dimensional body of dimensions $a\times b\times c$ . The temperature of a point in the body is defined as a function of the temperature field $T=T\left(x,y,z,\tau \right)$ and $\left(x,y,z\right)$ being coordinates of the position of a point in the body and the $\tau$ is the time. If the body is homogeneous, the temperature field is described by the equation:

$\frac{\partial T}{\partial \tau }=k\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}+\frac{{\partial }^{2}T}{\partial z^2}\right)$ , (42)

where is k thermal diffusivity of the material from which the body is made.

Equation (42) is a parabolic equation (classical equation) of heat conduction that predicts an infinite rate of propagation of a thermal disturbance in a body. This contradiction can be overcome by introducing a generalized heat conduction equation.

In this sense, let’s consider the three-dimensional generalized heat conduction equation:

${\tau }^{\ast }\frac{{\partial }^{2}T}{\partial {\tau }^2}+\frac{\partial T}{\partial \tau }=k\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}+\frac{{\partial }^{2}T}{\partial z^2}\right)$ . (43)

Equation (43) is a hyperbolic partial differential equation, the solution of which will be sought by the Fourier method of variable separation under given initial and boundary conditions.

Initial condition:

$T\left(x,y,z,0\right)=f\left(x,y,z\right)=x\left(a-x\right)y\left(b-y\right)z\left(c-z\right)$ , $T_{\tau }\left(x,y,z,0\right)=0$ . (44)

Boundary condition:

$T\left(0,y,z\right)=T\left(a,y,z\right)=0$ , (45)

$T\left(x,0,z\right)=T\left(x,b,z\right)=0$ , (46)

$T\left(x,y,0\right)=T\left(x,y,c\right)=0$ , (47)

$0\le x\le a,0\le y\le b,0\le z\le c,\tau >0$ .

Similar to the 2D Fourier method of separation of variables, the problems (43)-(47) are solved:

$\begin{array}{l}T\left(x,y,z,\tau \right)\\ ={\displaystyle \underset{m=1}{\overset{\infty }{\sum }}{\displaystyle \underset{n=1}{\overset{\infty }{\sum }}{\displaystyle \underset{l=1}{\overset{\infty }{\sum }}\frac{E_{mnl}}{r_{2,mnl}-r_{1,mnl}}\cdot \left(r_{2,mnl}\cdot {\text{e}}^{r_{1,mnl}}-r_{1,mnl}\cdot {\text{e}}^{r_{2,mnl}}\right)\cdot \sin \left({\alpha }_{m}x\right)\sin \left({\beta }_{n}y\right)\sin \left({\gamma }_{l}z\right)}}},\end{array}$ (48)

finding:

$E_{mnl}=\frac{8}{abc}\cdot {\displaystyle \underset{0}{\overset{a}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}{\displaystyle \underset{0}{\overset{c}{\int }}f\left(x,y,z\right)\sin \left({\alpha }_{m}x\right)\sin \left({\beta }_{n}y\right)\sin \left({\gamma }_{l}z\right)\text{d}x\text{d}y\text{d}z}}}$ ,

${\alpha }_m=\frac{m\pi }{a},\left(m=1,2,\cdots \right)$ ,

${\beta }_n=\frac{n\pi }{a},\left(n=1,2,\cdots \right)$ ,

${\gamma }_l=\frac{l\pi }{a},\left(l=1,2,\cdots \right)$ ,

${\delta }_{mnl}^2={\alpha }_m^2+{\beta }_n^2+{\gamma }_l^2,$

$r_{1,mnl}=\frac{-1+\sqrt{1-4k{\delta }_{mnl}^2{\tau }^{\ast }}}{2{\tau }^{\ast }},$

$r_{2,mnl}=\frac{-1-\sqrt{1-4k{\delta }_{mnl}^2{\tau }^{\ast }}}{2{\tau }^{\ast }},$

${\tau }^{\ast }\in \left(0,\frac{1}{4k{\delta }_{mnl}^2}\right].$

$T\left(x,y,z,0\right)=f\left(x,y,z\right)=x\left(a-x\right)y\left(b-y\right)z\left(c-z\right)$ initial contitions given, it is resulting in:

$E_{mnl}=\frac{512a^2{b}^2{c}^2}{m^3{n}^3{l}^3{\pi }^9},\left(m,n,l=1,2,\cdots \right).$

2.4. Determining the Approximate Solution of the Three-Dimensional Heat Conduction Equation

Let us now find an approximate solution to the generalized (hyperbolic) heat conduction equation (Equation (43)). An approximate solution of the form will be assumed:

$T\left(x,y,z,\tau \right)=f\left(x,y,z\right)\cdot \phi \left(\tau \right)=x\left(a-x\right)y\left(b-y\right)z\left(c-z\right)\cdot \phi \left(\tau \right)$ , (49)

where $\phi \left(\tau \right)$ is anunknown function that satisfies the condition:

$\phi \left(0\right)=1$ . (50)

If we start from the functional:

$\begin{array}{l}F\left(T\left(x,y,z,\tau \right)\right)={\displaystyle \underset{0}{\overset{\tau }{\int }}\text{d}\tau }{\displaystyle \underset{V}{\iiint }F\left(x,y,z,\tau ,T,{T^{\prime }}_x,{T^{\prime }}_y,{T^{\prime }}_z,{T^{\prime }}_{\tau }\right)\text{d}x\text{d}y\text{d}z}\\ ={\displaystyle \underset{0}{\overset{\tau }{\int }}{\displaystyle \underset{0}{\overset{a}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}{\displaystyle \underset{0}{\overset{c}{\int }}\left[\frac{{\tau }^{\ast }}{2}\cdot {\left(\frac{\partial T}{\partial \tau }\right)}^2-\frac{k}{2}\cdot \left({\left(\frac{\partial T}{\partial x}\right)}^2+{\left(\frac{\partial T}{\partial y}\right)}^2+{\left(\frac{\partial T}{\partial z}\right)}^2\right)\right]}}}}\cdot {\text{e}}^{\frac{\tau }{{\tau }^{\ast }}}\text{d}\tau \text{d}x\text{d}y\text{d}z,\end{array}$ (51)

and by inserting the trial solution (49) into Equation (51) and solving the triple integral, we get:

$F=\frac{a^3{b}^3{c}^3}{5400}\cdot {\displaystyle \underset{0}{\overset{\tau }{\int }}\left(\frac{a^2{b}^2{c}^2{\tau }^{\ast }}{10}\cdot {\stackrel{\dot{}}{\phi }}^2-k\cdot \left(b^2{c}^2+a^2{c}^2+a^2{b}^2\right){\phi }^2\right)\cdot {\text{e}}^{\tau /{\tau }^{\ast }}\text{d}\tau }.$ (52)

If we now introduce the Lagrange function:

$L\left(\phi ,\stackrel{\dot{}}{\phi }\right)=\left(\frac{a^2{b}^2{c}^2{\tau }^{\ast }}{10}\cdot {\stackrel{\dot{}}{\phi }}^2-k\cdot \left(b^2{c}^2+a^2{c}^2+a^2{b}^2\right){\phi }^2\right)\cdot {\text{e}}^{\tau /{\tau }^{\ast }},$ (53)

and using the Euler-Lagrange equation:

$\frac{\partial L}{\partial \phi }-\frac{\partial }{\partial \tau }\left(\frac{\partial L}{\partial \stackrel{\dot{}}{\phi }}\right)=0,$ (54)

and sharing with ${\text{e}}^{\tau /{\tau }^{\ast }}$ , the following differential equation is obtained:

${\tau }^{\ast }\cdot \ddot{\phi }\left(\tau \right)+\stackrel{\dot{}}{\phi }\left(\tau \right)+\frac{10k\cdot \left(b^2{c}^2+a^2{c}^2+a^2{b}^2\right)}{a^2{b}^2{c}^2}\cdot \phi \left(\tau \right)=0$ . (55)

After a short calculation, the resulting solution to Equation (55):

$\phi \left(\tau \right)=C_1\cdot {\text{e}}^{r_1\tau }+C_2\cdot {\text{e}}^{r_2\tau }$ , (56)

given:

$r_1=\frac{-1+\sqrt{1-40kM{\tau }^{\ast }}}{2{\tau }^{\ast }}$ , $r_2=\frac{-1-\sqrt{1-40kM{\tau }^{\ast }}}{2{\tau }^{\ast }}$ , ${\tau }^{\ast }>0$ and ${\tau }^{\ast }\in \left(0,\frac{1}{40kM}\right]$ , $M=\frac{b^2{c}^2+a^2{c}^2+a^2{b}^2}{a^2{b}^2{c}^2}.$

Constants $C_1$ and $C_2$ are determined from the conditions of the conditions $\phi \left(0\right)=1$ and $T_{\tau }\left(x,y,z,0\right)=0$ . The first condition gives:

$C_1+C_2=1$ , (57)

and the second condition gives:

$C_1\cdot r_1+C_2\cdot r_2=0$ . (58)

The constants are determined from the system of Equations (57) and (58). $C_1$ and $C_2$ :

$C_1=\frac{r_2}{r_2-r_1}$ , (59)

$C_2=\frac{r_1}{r_1-r_2}$ , (60)

so an unknown function $\phi \left(\tau \right)$ is determined by:

$\phi \left(\tau \right)=\frac{r_2}{r_2-r_1}\cdot {\text{e}}^{r_1\tau }+\frac{r_1}{r_1-r_2}\cdot {\text{e}}^{r_2\tau }$ . (61)

Now, according to (49), the approximate solution of Equation (43) is given by:

$T\left(x,y,z,\tau \right)=x\left(a-x\right)y\left(b-y\right)z\left(c-z\right)\cdot \left(\frac{r_2}{r_2-r_1}\cdot {\text{e}}^{r_1\tau }+\frac{r_1}{r_1-r_2}\cdot {\text{e}}^{r_2\tau }\right)$ . (62)

If a boundary condition is applied to Equation (55) or Equation (62), ${\tau }^{\ast }\to 0$ an approximate solution to the classic parabolic three-dimensional heat conduction equation is obtained:

$\frac{\partial T}{\partial \tau }=k\cdot \left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}+\frac{{\partial }^{2}T}{\partial z^2}\right)$ , (63)

which is given by:

$T\left(x,y,z,\tau \right)=x\left(a-x\right)y\left(b-y\right)z\left(c-z\right)\cdot {\text{e}}^{-\frac{10k\left(b^2{c}^2+a^2{c}^2+a^2{b}^2\right)}{a^2{b}^2{c}^2}\cdot \tau }$ . (64)

3. Numerical Example

In the previous section, the analytical and approximate solution of the 2D and 3D generalized (hyperbolic) heat conduction equation was presented. At the relaxation time, the ${\tau }^{\ast }\to 0$ approximate solutions of the classical (Fourier) 2D and 3D heat conduction equations were obtained. In order to illustrate the relevant physical effects of the obtained results, two examples are considered. One example refers to the 2D heat conduction equation, and the other example to the 3D heat conduction equation. The material of a 5 × 5 cm aluminum plate and a 5 × 5 × 5 cm aluminum cube were chose. Aluminum was chosen because of its good characteristics and wide application in various branches of process technology. The thermal properties of aluminum are: density $\rho =2707\text{kg}/{\text{m}}^{\text{3}}$ , thermal conductivity $\lambda =204\text{W}/\left(\text{m}\cdot \text{K}\right)$ , specific heat capacity $c_p=896\text{J}/\left(\text{kg}\cdot \text{K}\right)$ [21] . The thermal diffusivity of aluminum is determined using the expression $k=\lambda /\left(\rho \cdot c_p\right)$ . According to the given data, the thermal diffusivity of aluminum amounts to $k=0.8411{\text{cm}}^{\text{2}}/\text{s}$ . The relaxation time for the material aluminum ranges in the interval ${\tau }^{\ast }={10}^{-14}\text{-}{10}^{-11}\text{s}$ [22] .

Relaxation time $\tau$ represents the statistical value of the delay time required to establish steady-state heat flow conditions in a small elementary volume of material when a temperature gradient suddenly appears at the boundary of the elementary volume. The relaxation time is quoted in the literature [22] for different materials and was found to range from ${\tau }^{\ast }={10}^{-10}\text{s}$ for gases to ${\tau }^{\ast }={10}^{-14}\text{s}$ za for metals.

Example 1. 2D heat conduction

We will consider the generalized (hyperbolic) heat conduction equation (Equation (5)), i.e.

${\tau }^{\ast }\cdot \frac{{\partial }^{2}T}{\partial {\tau }^2}+\frac{\partial T}{\partial \tau }=k\cdot \left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}\right)$ ,

and according to the condition of the task, the initial and boundary conditions are:

Initial condition: $T\left(x,y,0\right)=x\left(5-x\right)y\left(5-y\right)$ ,

Boundery condition: $T\left(0,y,\tau \right)=T\left(5,y,\tau \right)=0$ ,

$T\left(x,0,\tau \right)=T\left(x,5,\tau \right)=0$ ,

$0\le x\le 5$ , $0\le y\le 5$ , $\tau >0$ .

The considered problem will be solved using the exact solution (Equation (25)) and the approximate solution (Equation (39)). The results are shown in Table 1. Note that these two solutions agree well. Taking more terms in Equation (25) would give better results.

It can be observed that the value of the temperature field increases until the middle of the plate, and then the value of the temperature field decreases. Obviously, the hottest spot is the spot in the plate $x=2.5\text{cm}$ , $y=2.5\text{cm}$ and it is 41.73˚C for the exact solution and 39.06˚C for the approximate solution. It is also observed that over time, the profile of the temperature field practically does not change (Figure 1). MATLAB software was used to display the figures. During

Table 1. Comparative results of the temperature distribution obtained by applying the exact and approximate solution of the generalized 2D heat conduction equation of the 5 × 5 cm aluminum plate (relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ ).

Figure 1. Hyperbolic 2D profile of the temperature field in an aluminum plate (relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ ).

engineering calculations of heat conduction through solid materials, the practical question of determining time τ arises, during which the entire plate will cool down to a specific temperature $T^{\ast }$ . For the considered example, the time required for the aluminum plate to cool, for example, to 1˚C, it is necessary to solve the equation $T\left(5/2,5/2,\tau \right)=1$ and based on Equation (62), the result is $\tau \approx 5.45\text{s}$ .

Figures 2-5 show the comparation of parabolic and hyperbolicof the heat conduction profile in the aluminum plate at constant time. On the time micro scale $\left({\tau }^{\ast },\tau \right)\to \left(0,0\right)$ , the hyperbolic temperature profile is practically indistinguishable from the parabolic temperature profile (Figure 2). Due to the longer time of heat conduction through the aluminum plate (time macroscopic process),

Figure 2. Comparison of parabolic and hyperbolic temperature profiles for an aluminum plate at constant time $\tau =0.01\text{s}$ and relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ . (a) Hyperbolic and parabolic profile of the temperature field of aluminum plate; (b) 2D parabolic profile of the temperature field of the aluminum plate.

Figure 3. Comparison of parabolic and hyperbolic temperature profiles for an aluminum plate at constant time $\tau =0.1\text{s}$ and relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ . (a) Hyperbolic and parabolic profile of the temperature field of the aluminum plate; (b) 2D parabolic profile of the temperature field of the aluminum plate.

Figure 4. Comparison of parabolic and hyperbolic temperature profiles for an aluminum plate at constant time $\tau =1\text{s}$ and relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ . (a) Hyperbolic and parabolic profile of the temperature field of the aluminum plate; (b) 2D parabolic profile of the temperature field of the aluminum plate

Figure 5. Comparison of parabolic and hyperbolic temperature profiles for an aluminum plate at constant time $\tau =10\text{s}$ and relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ . (a) Hyperbolic and parabolic profile of the temperature field of the aluminum plate; (b) 2D parabolic profile of the temperature field of the aluminum plate.

greater differences in the temperature profile of the hyperbolic (generalized) equation than the parabolic (classical) equation are observed (Figures 3-5).

Example 2. 3D heat conduction

As an example, the temperature field will be considered $T\left(x,y,z,\tau \right)$ of a 5 × 5 × 5 cm aluminum cube. The physical properties of aluminum are shown in example 1. The plate has a constant temperature 0˚C at the edges and its initial temperature is $T\left(x,y,z\right)=x\left(5-x\right)y\left(5-y\right)z\left(5-z\right),˚\text{C}$ . The relaxation time of the aluminum plate is ${\tau }^{\ast }={10}^{-11}\text{s}$ .

A comparison of the distribution of the temperature field of the aluminum cube obtained by applying the exact solution (Equation (48)) and the approximate solution (Equation (62)) is shown in Figure 6. It can be observed that the values of the temperature field of the cube are uniform except in the middle of the cube where there are deviations of about 10%. Taking more terms in Equation (48) would give better results.

Figures 7-9 show comparison of parabolic and hyperbolic of the heat conduction profile in an aluminum cube at constant time. Due to the longer time of heat conduction through the aluminum cube (time macroscopic process), greater differences in the temperature profile of the hyperbolic (generalized) equation than the parabolic (classical) equation are observed.

4. Conclusions

The paper considers multidimensional 2D and 3D heat conduction. The classical (parabolic) and generalized (hyperbolic) equation of heat conduction is considered. The solutions of the two-dimensional and three-dimensional (2D and 3D) classical and generalized heat conduction equations were obtained using the Fourier method of separation of variables. The solutions are quite complex and

Figure 6. Profile of the temperature field of three-dimensional heat conduction depending on depth x, y and z in an aluminum cube at a fixed time $\tau =0.01\text{s}$ and a fixed relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ .

Figure 7. Profile of the temperature field of three-dimensional heat conduction depending on depth x, y and z in an aluminum cube at a fixed time $\tau =0.1\text{s}$ and a fixed relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ .

Figure 8. Profile of the temperature field of three-dimensional heat conduction depending on depth x, y and z in an aluminum cube at a fixed time $\tau =1\text{s}$ and a fixed relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ .

Figure 9. Profile of the temperature field of three-dimensional heat conduction depending on depth x, y and z in an aluminum cube at a fixed time $\tau =10\text{s}$ and a fixed relaxation time ${\tau }^{\ast }={10}^{-11}\text{s}$ .

impractical for wider engineering applications.

It is shown that by using calculus of variations by applying the Euler-Lagrange equation, approximate solutions of multidimensional heat conduction (2D and 3D) can be obtained. The model can also be adapted to the classical Fourier equation of heat conduction (classical diffusion) by letting the relaxation time ${\tau }^{\ast }\to 0$ . The main advantage of the proposed model (approximate solution) is that it avoids the solution for 2D and 3D, which is shown in the form of rows obtained by the Fourier method of variable separation or by derivation integration using the Green’s function or using the Riemann method.

Numerical example for a 5 × 5 cm aluminum plate and a 5 × 5 × 5 cm aluminum cube showed the effectiveness of the proposed approximate solution. In addition, the approximate method used is easy to find solutions for heat conduction and use them in engineering practice. The material used in the numerical example is aluminum due to its good thermal properties, its abundance and low price on the market. In the future, this research can be used to help solve 1D, 2D and 3D heat conduction phenomena by considering the boundary conditions used.

These studies also indicate that in practical engineering calculations, a hyperbolic heat conduction equation would be preferable when representing actual physical processes.

Since this manuscript is devoted to the application of calculus of variations to 2D and 3D heat conduction through solid bodies with constant thermophysical coefficients $\rho$ , c_p, $\lambda$ , and k, further considerations should be devoted to the study of heat conduction with variable thermophysical coefficients because they depend on changes in their temperature fields, which is required and dictated by engineering practice.

Acknowledgements

The research by I. Aranđelović is supported in part by the Serbian Ministry of Science, Technological Development and Innovations according to Contract 451-03-65/2024-03/200105, dated 5 February 2024.

Nomenclature

c_p: specific heat capacity of the material (kJ/kg∙K)

D: two-dimensional area

F: functional

k: thermal diffusivity of the material (cm²/s)

L: length (cm)

L: Lagrangian

R³: three-dimensional Euclidean space

R₊: the set of positive real numbers

T: temperature (◦C)

V: volume (cm³)

x: length (cm)

y: length (cm)

z: length (cm)

τ: time (s)

${\tau }^{\ast }$ : relaxation time (s)

ρ: material density (kg/m³)

λ: coefficient of thermal conductivity of the material (W/m∙K)

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Xu, M. (2019) A High-Order Finite Volume Scheme for Unsteady Convection-Dominated Convection-Diffusion Equations. Numerical Heat Transfer Part B: Fundamentals, 76, 253-272. https://doi.org/10.1080/10407790.2019.1665421
[2]	Andrianov, I.V., Topol, H. and Danishevskyy, V.V. (2017) Asymptotic Analysis of Heat Transfer in Composite Materials with Nonlinear Thermal Properties. International Journal of Heat and Mass Transfer, 111, 736-754. https://doi.org/10.1016/j.ijheatmasstransfer.2017.03.124
[3]	Mironov, L.V. (2021) Elliptic Equations of Heat Transfer and Diffusion in Solids. Institute for Physik of Microstructures, Russian Academy of Sciences, GSP-105, Nizhny Novgorod. https://www.semanticscholar.org
[4]	Lin, J., Zhang, Y. and Sergiy, R. (2021) A Semi-Analitycal Method 1D, 2D and 3D Time Fractional Second Order Dual-Phase-Lag Model of the Heat Transfer. Alexandria Engineering Journal, 60, 5879-5896. https://doi.org/10.1016/j.aej.2021.03.071
[5]	Wang, C., et al. (2021) Analysis of 2D Heat Conduction in Nonlinear Functionally Graded Materials Using a Local Semi-Analytical Meshless Method. AIMS Mathematics, 6, 12599-12618. https://doi.org/10.3934/math.2021726
[6]	Nimona, K.K. and Gizaw, D.H. (2022) Analysis of Two-Dimensional Heat Transfer Problem Using the Boundary Integral Equation. Advances in Mathematical Physics, 2022, Article ID: 1889774. https://doi.org/10.1155/2022/1889774
[7]	Jannot, Y. (2023) Heat Transfer, Volume 1: Conduction and Convection. John Wiley & Sons, Inc., Hoboken. https://www.wiley.com>en-cn
[8]	Drissa, O., Serge, W.I., Gael, L.S., Abdoulaye, C., Belkacem, Z. and Xavier, C. (2020) Modeling and Numerical Simulation of Heat Transfer in a Metallic Pressure Cooker Isolated with Kapok Wool. Modeling and Numerical Simulation of Material Science, 10, 15-30. https://www.scirp.org/journal/mnsms https://doi.org/10.4236/mnsms.2020.102002
[9]	Dalal, A.M., Nujud, M.A. and Eman, S.A. (2020) Finite Difference Approximation for Solving Transiet Heat Conduction Equation for Copper. Advances in Pure Mathematics, 10, 350-358. https://www.scirp.org/journal/apm https://doi.org/10.4236/apm.2020.105021
[10]	Cattaneo, C. (1958) Sur une forme de l’équation de la chaleuréliminant le paradoxe d’une propagation instantanée. Comptes Rendus de l’Académie des Sciences, 247, 431-433.
[11]	Cattaneo, C. (1949) Sul calcolo di alcuni potenziali e sul loro intervento nella risoluzione di particolari problemi Armonici. Atti del Seminario Matematico e Fisico dell Universita di Modena e Reggio Emilia, 3, 29-45.
[12]	Vernotte, P. (1958) La veritable equation de la chaleur. Comptes Rendus de l’Academie des Sciences de Paris, France, 247, 2103-2105.
[13]	Đurić, S., Aranđelović, I. and Milotić, M. (2024) Variational Approach to Heat Conduction Modeling. Journal of Applied Mathematics and Phisics, 12, 234-248. https://www.scirp.org/journal/jamp https://doi.org/10.4236/jamp.2024.121018
[14]	Iqbal, K., Kulvinder, S. and Eduard-Marius, C. (2023) New Modified Couple Stress Theory of Thermoelasticity with Hyperbolic Two Temperature. Mathematics, 11, Article 432. https://doi.org/10.3390/math11020432
[15]	Zarei, A. and Pilla, S. (2023) An Improved Theory of Thermoelasticity for Ultrafast Heating of Materials Using Short and Ultrashort Laser Pulses. International Journal of Heat and Mass Transfer, 215, Article ID: 124510. https://doi.org/10.1016/j.ijheatmasstransfer.2023.124510
[16]	Eman, M.H. (2023) Fractional Model in the Theory of Generalized Thermoelastic Diffusion. Special Topics & Reviews in Porous Media: An International Journal, 14, 1-16. https://www.dl.begellhouse.com https://doi.org/10.1615/SpecialTopicsRevPorousMedia.2022044574
[17]	Ya, J. (2018) Analytical Solutions to Hyperbolic Heat Conductive Models Using Green’s Function Method. Journal of Thermal Science and Technology, 13, JTST0012. https://doi.org/10.1299/jtst.2018jtst0012
[18]	Povstenko, Y. (2016) Generalized Theory of Diffusive Stresses Associated with the Fractional Diffusion Equation and Nonlocal Constitutive Equations for the Stress Tensor. Computers & Mathematics with Applications, 78, 1819-1825. https://doi.org/10.1016/j.camwa.2016.02.034
[19]	Zhang, Z.W., et al. (2022) Heat Conduction Theory Including Phonon Coherence. Physical Review Letters, 128, Article ID: 015901. https://doi.org/10.1103/PhysRevLett.128.015901
[20]	Arber, T.D., Goffrey, T. and Ridgers, C. (2023) Models of Thermal Conduction and Non-Local Transport of Relevance to Space Physics with Insights from Laser-Plasma Theory. Frontiers in Astronomy and Space Sciences, 10, Article ID: 1155124. https://doi.org/10.3389/fspas.2023.1155124
[21]	Farshad, S., Maryam, G., Juan, E.-D. and Masud, B. (2021) Recent Advances in Generalized Thermoelasticity Theory and the Modified Models: A Review. Journal of Computational Design and Engineering, 8, 15-35. https://doi.org/10.1093/jcde/qwaa082
[22]	Tamma, K.K. and Zhou, X.M. (2010) Macroscale and Microscale Thermal Transport and Thermo-Mechanical Interactions: Some Note Worthy Perspectives. Journal of Thermal Stresses, 21, 405-449. https://doi.org/10.1080/01495739808956154