User:Tmhoang81/A thermal problem

= The transient heat transfer =

The dimensional thermal problem
The governing equation has the following form:
 * {| style="width:100%" border="0"

$$  \displaystyle k \frac{\partial^2 u}{\partial x^2} + f(x,t) = \rho c \frac{\partial u}{\partial t} $$ where $$\displaystyle k,\, c,\, \rho$$ are material thermal conductivity, specific heat and density proportionally.
 * style="width:95%" |
 * style="width:95%" |
 * style= | (1)
 * }

The initial condition:
 * {| style="width:100%" border="0"

$$\displaystyle u(x,t=0) = u_0 $$
 * style="width:95%" |
 * style="width:95%" |
 * style= | (2)
 * }

The boundary conditions:
 * {| style="width:100%" border="0"

$$\displaystyle q(x=0,t) = q(x=L,t) = 0 $$
 * style="width:95%" |
 * style="width:95%" |
 * style= | (3)
 * }

Using the Galerkin's Finite Element Method (FEM) as appendix 1 for a transient heat transfer problem, we obtain the global equations:
 * {| style="width:100%" border="0"

\displaystyle \left[ C \right]\left\{ {\dot U} \right\} + \left[ K \right]\left\{ U \right\} = \left\{ {F_Q } \right\} + \left\{ {F_g } \right\} $$ where vector $$\left\{ U \right\} = \left[ {\begin{array}{*{20}c} {u_1 } & \ldots  & {u_n }  \\ \end{array}} \right]^T $$ is temperature at n nodes, that obtains by solving the original differential equation (4) with the initial condition (2).
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (4)
 * }

We will demonstrate the above procedure through an example.

Example 1: Solving the 1D following heat equation in a cylindrical rod having diameter d and length L:
 * {| style="width:100%" border="0"

$$  \displaystyle k \frac{\partial^2 u}{\partial x^2} + f(x,t) = \rho c \frac{\partial u}{\partial t}  \qquad 0 \le x \le L,\,\,0 \le t \le t_{end} = 30 s $$ with material parameters
 * style="width:95%" |
 * style="width:95%" |
 * style= |
 * }
 * $$k = 230\,W/m{}^0C,\,\,\,c = 900\,J/kg{}^0C,\,\,\,\rho = 2700kg/m^3

$$ and geometric parameters
 * $$L = 0.1m,\,\,\,d = 0.012m\,

$$

Heat generation: $$f(x,t) = \begin{cases} Q = 10^7\,\,\,\,\,\,\,\,\,\,\,\left( {x,t} \right) \in \left[ {0,\frac{L}{2}} \right] \times \left[ {0,t_1 } \right] \\ 0\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,otherwise \end{cases} $$

Initial condition:
 * {| style="width:100%" border="0"

$$\displaystyle u(x,t=0) = u_0 = 30 ^0C $$
 * style="width:95%" |
 * style="width:95%" |
 * style= |
 * }

The boundary condition:
 * {| style="width:100%" border="0"

$$\displaystyle q(x=0,t) = q(x=L,t) = 0 $$
 * style="width:95%" |
 * style="width:95%" |
 * style= |
 * }

Solution

Dividing the rod into four elements, each has length $$L_e = \frac{L}{4}$$. Thus, we have five nodes with five freedom degrees: $$\left\{ U \right\} = \left[ {\begin{array}{*{20}c}{u_1 } & {u_2 } & {u_3 } & {u_4 } & {u_5 } \\ \end{array}} \right]^T $$

FEM equations for each element are given as appendix 2

Subsituting numerical values into the global equation (4) we have:
 * {| style="width:100%" border="0"

\displaystyle \left[ K \right] = \left[ {\begin{array}{*{20}c} {1.0403} & { - 1.0403} & 0 & 0 & 0 \\   { - 1.0403} & {2.0806} & { - 1.0403} & 0 & 0  \\   0 & { - 1.0403} & {2.0806} & { - 1.0403} & 0  \\   0 & 0 & { - 1.0403} & {2.0806} & { - 1.0403}  \\   0 & 0 & 0 & { - 1.0403} & {1.0403}  \\ \end{array}} \right] $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }


 * {| style="width:100%" border="0"

\displaystyle \left[ C \right] = \left[ {\begin{array}{*{20}c} {2.2898} & {1.1449} & 0 & 0 & 0 \\   {1.1449} & {4.5796} & {1.1449} & 0 & 0  \\   0 & {1.1449} & {4.5796} & {1.1449} & 0  \\   0 & 0 & {1.1449} & {4.5796} & {1.1449}  \\   0 & 0 & 0 & {1.1449} & {2.2898}  \\ \end{array}} \right] $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }


 * {| style="width:100%" border="0"

\displaystyle \left\{ {F_Q } \right\} = \left\{ {\begin{array}{*{20}c} {14.1345} \\   {28.2690}  \\   {14.1345}  \\   0  \\   0  \\ \end{array}} \right\} \left\{ {F_g } \right\} = \left\{ {\begin{array}{*{20}c} 0 \\   0  \\   0  \\   0  \\   0  \\ \end{array}} \right\} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }

To obtain solutions for first-order ordinary differential equations, we use the finite difference method or specifically the forward difference method. This method approximate the time derivative of the nodal temperature matrix as:

\displaystyle \left\{ {\dot U} \right\} \cong \frac $$

Substituting the above equation, the global equation becomes

\displaystyle \left\{ {U(t_{i + 1} )} \right\} = \left\{ {U(t_i )} \right\} + \left[ C \right]^{ - 1} \left( {\left\{ {F_Q (t_i )} \right\} + \left\{ {F_q (t_i )} \right\} - \left[ K \right]\left\{ {U(t_i )} \right\}} \right)\Delta t,\,\,\,\left\{ {U(t_0 )} \right\} = \left\{ {U_0 } \right\} $$

In this example, this equation yields two cases:

a) From $$0 \le t \le t_1$$: heat generation is non-zero

\displaystyle \left\{ {U(t_{i + 1} )} \right\} = \left\{ {U(t_i )} \right\} + \left[ C \right]^{ - 1} \left( {\left\{ {F_Q (t_i )} \right\} - \left[ K \right]\left\{ {U(t_i )} \right\}} \right)\Delta t,\,\,\,\left\{ {U(t_0 )} \right\} = \left\{ {U_0 } \right\} $$

b) From $$t_1 \le t \le t_{end}$$: heat generation is zero

\displaystyle \left\{ {\tilde U(t_{i + 1} )} \right\} = \left\{ {\tilde U(t_i )} \right\} - \left[ C \right]^{ - 1} \left[ K \right]\left\{ {\tilde U(t_i )} \right\}\Delta t ,\,\,\,\left\{ {\tilde U(t_0 )} \right\} = \left\{ {U(t_1 )} \right\} $$

Therefore, nodal temperatures are:

\displaystyle \left\{ U \right\} = \left[ {\begin{array}{*{20}c} {U(t_0 )} & {U(\Delta t)} & \ldots  & {U(t_1 )} & {\tilde U(t_1  + \Delta t)} &  \ldots  & {\tilde U(t_{end} )}  \\ \end{array}} \right] $$

With $$t_{end} = 30\,s,\,\,t_1  = \frac{3} = 10\,s,\,\,\Delta t = 0.1\,s$$, we plot values $$\left\{ { U(t)} \right\}$$ versus $$t$$.



The nondimensional thermal problem
Now we nondimensionalize the heat transfer equation (1). For this equation, there are four fundamental units $$\left( {F_1 ,F_2 ,F_3 ,F_4 } \right) = \left( {m,s,K,W} \right)$$ and ten physical quantities which are listed below together with their corresponding dimensional formula for the unit in terms of the fundamental units:

Thus, the dimension matrix is:

\displaystyle A = \left[ {\begin{array}{*{20}c} 0 & 0 & 1 & 1 & 0 & 0 & { - 1} & { - 3} & { - 2} & { - 3} \\   0 & 0 & 0 & 0 & 1 & 1 & 0 & 1 & 0 & 0  \\   1 & 1 & 0 & 0 & 0 & 0 & { - 1} & { - 1} & 0 & 0  \\   0 & 0 & 0 & 0 & 0 & 0 & 1 & 1 & 1 & 1  \\ \end{array}} \right]_{4 \times 10} $$ It can find the nullspace of matrix $$A$$:

\displaystyle Null\left( A \right) = \left[ {\begin{array}{*{20}c} { - 1} & 0 & 0  & 0  & { - 1}  & { - 1}   \\   1 & 0  & 0  & 0  & 0  & 0   \\   0  & { - 1}  & 0  & 2  & 1  & 2   \\   0  & 1  & 0  & 0  & 0  & 0   \\   0  & 0  & { - 1}  & { - 1}  & 0  & 0   \\   0  & 0  & 1  & 0  & 0  & 0   \\   0  & 0  & 0  & { - 1}  & { - 1}  & { - 1}   \\   0  & 0  & 0  & 1  & 0  & 0   \\   0  & 0  & 0  & 0  & 1  & 0   \\   0  & 0  & 0  & 0  & 0  & 1   \\ \end{array}} \right]_{10 \times 6} $$

From the above result, we have six independent dimensionless groups of the original variables.

\displaystyle R_1 ^{ - 1} R_2^1 = \frac{U} = {\prod}_1 ,\,\,R_3 ^{ - 1} R_4^1 \, = \frac{x}{L} = {\prod}_2 ,\,\,\,R_5 ^{ - 1} R_6^1  = \frac{t} = {\prod}_3 ,\,\,\, $$

\displaystyle R_3^2 R_5 ^{ - 1} R_7 ^{ - 1} R_8^1 = \frac = {\prod}_4 ,\,\,\,R_1^{ - 1} R_3 ^1 R_7 ^{ - 1} R_9^1  = \frac = {\prod}_5 ,\,\,\,R_1^{ - 1} R_3 ^2 R_7 ^{ - 1} R_{10}^1  = \frac = {\prod}_6 $$

Alternatively we have 6 nondimensional variables which are provided as below:

\displaystyle \frac{U} = \Gamma \,\,\,\,\,\,(5) \,\,\,\,\,\,\,\, \frac{x}{L} = \xi \,\,\,\,\,\, (6) \,\,\,\,\,\,\,\, \frac{t} = \tau \,\,\,\,\,\, (7) $$



\frac = \beta \,\,\,\,\,\, (8)\,\,\,\,\,\,\,\, \frac = \varphi \left( t \right) \,\,\,\,\,\, (9) \,\,\,\,\,\, \frac = \Phi \left( {\xi ,\tau } \right) \,\,\,\,\,\, (10) $$

In terms of these nondimensional variables the simplified thermal problem can be written in the following nondimensional form (the details of nondimensionalizing the heat problem are found in appendix 3):

Heat conduction equation:

\beta \frac + \Phi \left( {\xi ,\tau } \right) = \frac \,\,\,\,\,\, (11) ,\,\,\left( {\xi ,\tau } \right) \in \left[ {0,1} \right] \times \left[ {0,1} \right] $$

Initial condition:

\displaystyle \Gamma \left( {\xi ,\tau = 0} \right) = 1 \,\,\,\,\,\, (12) $$

Boundary condition:

\displaystyle \left. { - \frac} \right|_{\xi = 0}  = \varphi _{in} ,\,\,\left. { - \frac} \right|_{\xi = 1}  = \varphi _{out} \,\,\,\,\,\, (13) $$

where

\varphi _{in} = \frac, \,\,\,\,\,\, \varphi _{out} = \frac \,\,\,\,\,\, (14) $$

The complete set of nondimensional variables needed for the problem is given in equations (5)-(10). The nondimensional temperature $$\Gamma$$ can be expressed function of the nondimensional distance $$\xi$$ and the nondimensional time $$\tau$$ as well as function of four other nondimensional parameters.

The physical interpretation of the remaining four nondimensional parameters in equations (8)-(10) is the following. $$\beta$$, the Fourier number or a nondimensional thermal diffusivity, is the ratio of the rate of heat conduction and the rate of heat storage(thermal energy storage). $$\varphi$$ is the ratio of the incident heat flux and the rate of heat conduction, or can be seen as effective heat flux. Finally, $$\Phi$$ is the ratio of the heat source and the rate of heat storage, can interpreted as effective heat source.

Example 2: Solving the problem in example 1 with the nondimensional form. Using Galerkin's FEM, similarly to dimensional form, we have the desired FEM equation for one element:

\displaystyle \left[ {K^{\left( e \right)} } \right]\left\{ {\begin{array}{*{20}c} {\Gamma _1 \left( \tau \right)}  \\ {\Gamma _2 \left( \tau \right)}  \\ \end{array}} \right\} + \left[ {C^{\left( e \right)} } \right]\left\{ {\begin{array}{*{20}c} {\dot \Gamma _1 \left( \tau \right)}  \\ {\dot \Gamma _2 \left( \tau \right)}  \\ \end{array}} \right\} = \left\{ {F_Q^{\left( e \right)} } \right\} + \left\{ {F_g^{\left( e \right)} } \right\} $$ where

\displaystyle \left[ {K^{\left( e \right)} } \right] = \frac{\beta }\left[ {\begin{array}{*{20}c} 1 & { - 1} \\   { - 1} & 1  \\ \end{array}} \right]\,\,\,\,\,\,\,\,\,\,\left[ {C^{\left( e \right)} } \right] = \frac{6}\left[ {\begin{array}{*{20}c} 2 & 1 \\   1 & 2  \\ \end{array}} \right] $$



\displaystyle \left\{ {F_Q^{\left( e \right)} } \right\} = \left\{ {\begin{array}{*{20}c} {\int\limits_0^{L_e } {\Phi \left( {\xi ,\tau } \right)N_1 \left( \xi \right)d\xi } }  \\ {\int\limits_0^{L_e } {\Phi \left( {\xi ,\tau } \right)N_2 \left( \xi \right)d\xi } }  \\ \end{array}} \right\} = \frac{2}\left\{ {\begin{array}{*{20}c} 1 \\   1  \\ \end{array}} \right\} \left( {if \,\,\, \Phi  = const} \right) \,\,\,\,\,\,\,\,\,\,\,\left\{ {F_g^{\left( e \right)} } \right\} = \beta \left\{ {\begin{array}{*{20}c} {\left. { - \frac} \right|_{\xi = 0} }  \\ {\left. {\frac} \right|_{\xi = 1} }  \\ \end{array}} \right\} = \beta \left\{ {\begin{array}{*{20}c} {\varphi _1 } \\ { - \varphi _2 } \\ \end{array}} \right\} $$

Using assembly procedure to obtain the global equations:

\displaystyle \left[ C \right]\left\{ {\dot \Gamma (\tau)} \right\} + \left[ K \right]\left\{ \Gamma (\tau) \right\} = \left\{ {F_Q } \right\} + \left\{ {F_g } \right\} \,\,\,\,\,\, (15) $$

Dividing the rod into four element, each element has length $$L_e = \frac{1}{4}$$, we have 5 nodes with 5 freedom degree:

\left\{ \Gamma \right\} = \left[ {\Gamma _1 } \,\,\, {\Gamma _2 } \,\,\, {\Gamma _3 } \,\,\, {\Gamma _4 } \,\,\, {\Gamma _5 } \,\,\, \right]^T \,\,\,\,and\,\,\,\left\{ {\dot \Gamma } \right\} = \left[ {\dot \Gamma _1 } \,\,\, {\dot \Gamma _2 } \,\,\, {\dot \Gamma _3 } \,\,\, {\dot \Gamma _4 } \,\,\, {\dot \Gamma _5 } \,\,\, \right]^T $$

where



\left[ K \right] = \frac{\beta }\left[ {\begin{array}{*{20}c} 1 & { - 1} & 0 & 0 & 0 \\   { - 1} & 2 & { - 1} & 0 & 0  \\   0 & { - 1} & 2 & { - 1} & 0  \\   0 & 0 & { - 1} & 2 & { - 1}  \\   0 & 0 & 0 & { - 1} & 1  \\ \end{array}} \right] \,\,\,\,\,\, (16) $$



\left[ C \right] = \frac{6}\left[ {\begin{array}{*{20}c} 2 & 1 & 0 & 0 & 0 \\   1 & 4 & 1 & 0 & 0  \\   0 & 1 & 4 & 1 & 0  \\   0 & 0 & 1 & 4 & 1  \\   0 & 0 & 0 & 1 & 2  \\ \end{array}} \right] \,\,\,\,\,\, (17) $$



\left[ F_Q \right] = \frac{2}\left\{ {\begin{array}{*{20}c} 1 \\   2  \\   1  \\   0  \\   0  \\ \end{array}} \right\} \,\,\,\,\,\, (18) \,\,\,\,\,\, \left[ F_g \right] = \beta \left\{ {\begin{array}{*{20}c} {\varphi _1 } \\ 0 \\   0  \\   0  \\   { - \varphi _5 }  \\ \end{array}} \right\} \,\,\,\,\,\, (19) $$

Substituting the values of parameters, it can easily calculate conductance and capacitance matrix, forcing function vectors as

\left[ K \right] = \left[ {\begin{array}{*{20}c} & { - {\rm{3}}{\rm{.0288}}} & 0 & 0 & 0 \\ { - {\rm{3}}{\rm{.0288}}} & & { - {\rm{3}}{\rm{.0288}}} & 0 & 0  \\ 0 & { - {\rm{3}}{\rm{.0288}}} & & { - {\rm{3}}{\rm{.0288}}} & 0  \\ 0 & 0 & { - {\rm{3}}{\rm{.0288}}} & & { - {\rm{3}}{\rm{.0288}}}  \\ 0 & 0 & 0 & { - {\rm{3}}{\rm{.0288}}} &  \\ \end{array}} \right],\,\,\,\,\,\,\,\, $$



\left[ C \right] = \left[ {\begin{array}{*{20}c} & & 0 & 0 & 0  \\    &  &  & 0 & 0  \\   0 &  &  &  & 0  \\   0 & 0 &  &  & {1.1449}  \\   0 & 0 & 0 & {1.1449} &   \\ \end{array}} \right] $$



\left\{ {F_Q } \right\} = \left\{ {\begin{array}{*{20}c} \\    \\     \\   0  \\   0  \\ \end{array}} \right\},\,\,\,\left\{ {F_g } \right\} = \beta \left\{ {\begin{array}{*{20}c} {\varphi _1 } \\ {\rm{0}} \\ {\rm{0}} \\ 0 \\   { - \varphi _5 }  \\ \end{array}} \right\} = \left\{ {\begin{array}{*{20}c} 0 \\   0  \\   0  \\   0  \\   0  \\ \end{array}} \right\} \,\,\,\,\,\, \left( {\varphi _1 = \varphi _{in}  = 0,\,\,\varphi _5  = \varphi _{out}  = 0} \right) $$

The graph of nondimensional temperature versus nondimensional time is below:



Example 3: Now we want to generate a graph that will give us the maximum value of the temperature in a rod as the example 1(max temperature over both time and space) for any value of $$\displaystyle k,c,\rho ,L,U_0 ,Q,t_1 $$ . We do this this by solving the non-dimensional equation (11)-(13) for a range of these values. It can easily find that the maximum value of the nondimensional temperature $$\displaystyle \Gamma _{\max } $$ depends on three independent nondimensional variables:
 * $$\beta = \frac, \,\,\,\,\,\, \bar Q = \frac, \,\,\,\,\,\, \tau _1  = \frac

$$. From equation (15), we can note that the maximum temperature $$\displaystyle \Gamma _{\max } $$ will be linear in $$\displaystyle \Phi = \bar Q = const$$ so that we do not need this as a parameter. To prove this conclusion, we scale the solution by changing variable:
 * $$\Lambda = \frac $$

At that time, the global equation (15) becomes:

\bar Q\left[ C \right]\left\{ {\dot \Lambda (\tau )} \right\} + \left[ K \right]\left\{ {\bar Q\Lambda (\tau ) + \left\{ 1 \right\}} \right\} = \left\{ F \right\} $$ or

\bar Q\left[ C \right]\left\{ {\dot \Lambda (\tau )} \right\} + \bar Q\left[ K \right]\left\{ {\Lambda (\tau )} \right\} = \left\{ {F_Q } \right\}\, - \left[ K \right]\left\{ 1 \right\} = \left\{ {F_Q } \right\} \,\,\,\,\,\, (because \left[ K \right]\left\{ 1 \right\} = \left\{ 0 \right\}) $$ Thus, dividing both sides of above equation to $$\bar Q \ne 0$$, we get

\left[ C \right]\left\{ {\dot \Lambda (\tau )} \right\} + \left[ K \right]\left\{ {\Lambda (\tau )} \right\} = \left\{ F \right\} $$ where $$K,\,C$$ is determined as equations (16),(17) and $$\displaystyle F = \frac{2}\left[ {1\,\,2\,\,1\,\,0\,\,0} \right]^T $$ with the initial condition $$\Lambda (0) = 0 $$

Summing up, the maximum nondimensional temperature $$\displaystyle \Gamma _{\max } $$ only depends on two independent nondimensional variables $$ \displaystyle \beta ,\tau _1 $$ The graph of $$\displaystyle \Gamma _{\max } $$ for a range of these two parameters is below:

Appendix 1: Finite element procedures for a thermal problem
Here the Galerkin's Finite Element Method (FEM) is applied for the heat transfer equation to obtain the element equations. A two-node element with linear interpolation functions is used and the temperature distribution in an element expressed as
 * {| style="width:100%" border="0"

$$  \displaystyle u(x,t) = N_1(x) u_1(t) + N_2(x) u_2(t) $$ where $$\displaystyle u_1(t) $$ and $$\displaystyle u_2(t) $$ are the temperature at nodes 1 and 2, which define the element, and the interpolation functions $$\displaystyle N_1(x) $$ and $$\displaystyle N_2(x) $$ are given by $$\displaystyle N_1(x) = (1- \frac{x}{L}) $$, $$\displaystyle N_2(x) = \frac{x}{L} $$
 * style="width:95%" |
 * style="width:95%" |
 * style= | (20)
 * }

Substitution of the discretized solution into the governing differential Equation results in the residual intergrals:
 * {| style="width:100%" border="0"

\displaystyle \int_{0}^{L} \left ( k \frac{\partial^2 u}{\partial x^2} + f(x,t) - \rho c \frac{\partial u}     {\partial t} \right ) N_i(x) A\, dx = 0 \quad i = 1,\,2 $$ where we note that the integration is over the volume of the element, that is, the domain of the problem, with $$\displaystyle dV = A dx $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (21)
 * }

or
 * {| style="width:100%" border="0"

\displaystyle k A \int_{0}^{L} N_i(x) \frac{\partial^2 u}{\partial x^2}\, dx + A \int_{0}^{L} N_i(x) f(x,t)\, dx - \rho c A \int_{0}^{L} N_i(x) \frac{\partial u}{\partial t}\, dx = 0 \quad i = 1,\,2 $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (22)
 * }

Intergratig the first term by parts and rearranging,
 * {| style="width:100%" border="0"

\displaystyle k A \int_{0}^{L} \frac{dN_i(x)}{dx} \frac{\partial u}{\partial x}\, dx  + \rho c   A \int_{0}^{L} N_i(x) \frac{\partial u}{\partial t}\, dx = A \int_{0}^{L} N_i(x) f(x,t)\, dx + \left. {kAN_i \left( x \right)\frac} \right|_{x = 0}^{x = L} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (23)
 * }

Substituting for $$\displaystyle u(x,t)$$ from (20) yields
 * {| style="width:100%" border="0"

\displaystyle k A \int_{0}^{L} \frac{dN_i(x)}{dx} \left ( \frac{dN_1(x)}{dx} u_1(t) + \frac{dN_2(x)}{dx} u_2(t)    \right )\, dx   + \rho c A \int_{0}^{L} N_i(x) \left ( N_1(x) \dot u_1(t) + N_2(x) \dot u_2(t) \right )\, dx $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }


 * {| style="width:100%" border="0"

\displaystyle = A \int_{0}^{L} N_i(x) f(x,t)\, dx + \left. {kAN_i \left( x \right)\frac} \right|_{x = 0}^{x = L} \quad i = 1,\,2 $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (24)
 * }

The two equations represented by equation (24) are conveniently combined into a matrix form by rewriting as
 * {| style="width:100%" border="0"

$$  \displaystyle u(x,t) = [N_1(x) N_2(x)] \begin{Bmatrix} u_1(t) \\ u_2(t) \end{Bmatrix} $$
 * style="width:95%" |
 * style="width:95%" |
 * style= | (25)
 * }

and substituting to obtain
 * {| style="width:100%" border="0"

\displaystyle k A \int_{0}^{L} \begin{bmatrix} \dot N_1(x) \\ \dot N_2(x) \end{bmatrix} \left [ \dot N_1(x) \dot N_2(x) \right ]dx\, \begin{Bmatrix} u_1(t) \\ u_2(t) \end{Bmatrix} + \rho c A \int_{0}^{L} \begin{bmatrix} N_1(x) \\ N_2(x) \end{bmatrix} \left [ N_1(x) N_2(x) \right ]dx\, \begin{Bmatrix} \dot u_1(t) \\ \dot u_2(t) \end{Bmatrix} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }


 * {| style="width:100%" border="0"

\displaystyle = A \int_{0}^{L} \begin{bmatrix} N_1(x) \\ N_2(x) \end{bmatrix} f(x,t)\, dx + k A \begin{Bmatrix} \left. { - \frac} \right|_{x = 0} \\ \left. {\frac} \right|_{^{x = L} } \end{Bmatrix} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (26)
 * }

Equation (26) is in the desired finite element form:
 * {| style="width:100%" border="0"

\displaystyle \left[ {k^{\left( e \right)} } \right]\left\{ {\begin{array}{*{20}c} {u_1 \left( t \right)} \\ {u_2 \left( t \right)} \\ \end{array}} \right\} + \left[ {c^{\left( e \right)} } \right]\left\{ {\begin{array}{*{20}c} {\dot u_1 \left( t \right)} \\ {\dot u_2 \left( t \right)} \\ \end{array}} \right\} = \left\{ {f_Q^{\left( e \right)} } \right\} + \left\{ {f_g^{\left( e \right)} } \right\} $$ where $$\left[ {k^{\left( e \right)} } \right]$$ is the conductance matrix defined as
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (27)
 * }
 * {| style="width:100%" border="0"

\displaystyle \begin{array}{l} \left[ {k^{\left( e \right)} } \right] = kA\int\limits_0^{L} {\left[ {\begin{array}{*{20}c} {\frac} \\ {\frac} \\ \end{array}} \right]} \left[ {\begin{array}{*{20}c} {\frac} & {\frac} \\ \end{array}} \right]dx = kA\int\limits_0^{L} {\left[ {\begin{array}{*{20}c} { - \frac{1}} \\ {\frac{1}} \\ \end{array}} \right]} \left[ {\begin{array}{*{20}c} { - \frac{1}} & {\frac{1}} \\ \end{array}} \right]dx \\ \,\,\,\,\,\,\,\,\,\,\,\, = \frac\left[ {\begin{array}{*{20}c} 1 & { - 1} \\   { - 1} & 1  \\ \end{array}} \right] \\ \end{array} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (28)
 * }

and $$\left[ {c^{\left( e \right)} } \right]$$ is the element capacitance matrix defined by
 * {| style="width:100%" border="0"

\displaystyle \begin{array}{l} \left[ {c^{\left( e \right)} } \right] = c\rho A\int_0^L {\left[ {\begin{array}{*{20}c} {N_1 } \\ {N_2 } \\ \end{array}} \right]\left[ {\begin{array}{*{20}c} {N_1 } & {N_2 } \\ \end{array}} \right]} dx = c\rho A\int_0^L {\left[ {\begin{array}{*{20}c} {1 - \frac{x}{L}} \\ {\frac{x}{L}} \\ \end{array}} \right]\left[ {\begin{array}{*{20}c} {1 - \frac{x}{L}} & {\frac{x}{L}} \\ \end{array}} \right]} dx \\ \,\,\,\,\,\,\,\,\,\,\, = c\rho AL\int_0^1 {\left[ {\begin{array}{*{20}c} {1 - \xi } \\ \xi  \\ \end{array}} \right]\left[ {\begin{array}{*{20}c} {1 - \xi } & \xi  \\ \end{array}} \right]} d\xi = \frac{6}\left[ {\begin{array}{*{20}c} 2 & 1 \\   1 & 2  \\ \end{array}} \right] \\ \end{array} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (29)
 * }

The forcing function vectors on the right-hand side of equation (27) include the internal heat generation and boundary flux terms. These are given by
 * {| style="width:100%" border="0"

\displaystyle \left[ {f_Q ^{\left( e \right)} } \right] = A\left\{ {\begin{array}{*{20}c} {\int_0^L {f(x,t)N_1 (x)} dx} \\ {\int_0^L {f(x,t)N_2 (x)} dx} \\ \end{array}} \right\} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (30)
 * }
 * {| style="width:100%" border="0"

\displaystyle \left[ {f_g ^{\left( e \right)} } \right] = kA\left\{ {\begin{array}{*{20}c} {\left. { - \frac} \right|_{x = 0} } \\ {\left. {\frac} \right|_{x = L} } \\ \end{array}} \right\} = A\left\{ {\begin{array}{*{20}c} {\left. q \right|_{x = 0} } \\ {\left. { - q} \right|_{x = L} } \\ \end{array}} \right\} = A\left\{ {\begin{array}{*{20}c} {q_1 } \\ { - q_2 } \\ \end{array}} \right\} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (31)
 * }

Using assembly procedure for a transient heat transfer problem, we obtain the global equations:
 * {| style="width:100%" border="0"

\displaystyle \left[ C \right]\left\{ {\dot U} \right\} + \left[ K \right]\left\{ U \right\} = \left\{ {F_Q } \right\} + \left\{ {F_g } \right\} $$ where vector $$\left\{ U \right\} = \left[ {\begin{array}{*{20}c} {u_1 } & \ldots  & {u_n }  \\ \end{array}} \right]^T $$ is temperature at n nodes, that obtains by solving the differential equation (32) with the initial condition (2).
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= | (32)
 * }

Appendix 2: Finite element implementation of a dimensional problem
FEM equations for each element are given as (27)

Element 1:


 * {| style="width:100%" border="0"

\displaystyle \frac\left[ {\begin{array}{*{20}c} 1 & { - 1} \\   { - 1} & 1  \\ \end{array}} \right]\left\{ {\begin{array}{*{20}c} {u_1 }  \\  {u_2 }  \\ \end{array}} \right\} + \frac{6}\left[ {\begin{array}{*{20}c} 2 & 1 \\   1 & 2  \\ \end{array}} \right] \left\{ {\begin{array}{*{20}c} {\dot u_1 } \\ {\dot u_2 } \\ \end{array}} \right\} = QA\left\{ {\begin{array}{*{20}c} {\int\limits_0^{L_e } {N_1 (x)dx} } \\ {\int\limits_0^{L_e } {N_2 (x)dx} } \\ \end{array}} \right\} + A\left\{ {\begin{array}{*{20}c} {q_1 } \\ { - q_2 } \\ \end{array}} \right\} $$ or
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }
 * {| style="width:100%" border="0"

\displaystyle \frac\left[ {\begin{array}{*{20}c} 1 & { - 1} \\   { - 1} & 1  \\ \end{array}} \right]\left\{ {\begin{array}{*{20}c} {u_1 }  \\  {u_2 }  \\ \end{array}} \right\} + \frac{6}\left[ {\begin{array}{*{20}c} 2 & 1 \\   1 & 2  \\ \end{array}} \right] \left\{ {\begin{array}{*{20}c} {\dot u_1 } \\ {\dot u_2 } \\ \end{array}} \right\} = \frac{Q A L_e}{2} \left\{ {\begin{array}{*{20}c} {1}  \\   {1}   \\ \end{array}} \right\} + A\left\{ {\begin{array}{*{20}c} {q_1 } \\ { - q_2 } \\ \end{array}} \right\} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }

Element 2:
 * {| style="width:100%" border="0"

\displaystyle \frac\left[ {\begin{array}{*{20}c} 1 & { - 1} \\   { - 1} & 1  \\ \end{array}} \right]\left\{ {\begin{array}{*{20}c} {u_2 }  \\  {u_3 }  \\ \end{array}} \right\} + \frac{6}\left[ {\begin{array}{*{20}c} 2 & 1 \\   1 & 2  \\ \end{array}} \right] \left\{ {\begin{array}{*{20}c} {\dot u_2 } \\ {\dot u_3 } \\ \end{array}} \right\} = \frac{Q A L_e}{2} \left\{ {\begin{array}{*{20}c} {1}  \\   {1}   \\ \end{array}} \right\} + A\left\{ {\begin{array}{*{20}c} {q_2 } \\ { - q_3 } \\ \end{array}} \right\} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }

Element 3:
 * {| style="width:100%" border="0"

\displaystyle \frac\left[ {\begin{array}{*{20}c} 1 & { - 1} \\   { - 1} & 1  \\ \end{array}} \right]\left\{ {\begin{array}{*{20}c} {u_3 }  \\  {u_4 }  \\ \end{array}} \right\} + \frac{6}\left[ {\begin{array}{*{20}c} 2 & 1 \\   1 & 2  \\ \end{array}} \right] \left\{ {\begin{array}{*{20}c} {\dot u_3 } \\ {\dot u_4 } \\ \end{array}} \right\} = \frac{Q A L_e}{2} \left\{ {\begin{array}{*{20}c} {1}  \\   {1}   \\ \end{array}} \right\} + A\left\{ {\begin{array}{*{20}c} {q_3 } \\ { - q_4 } \\ \end{array}} \right\} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }

Element 4:
 * {| style="width:100%" border="0"

\displaystyle \frac\left[ {\begin{array}{*{20}c} 1 & { - 1} \\   { - 1} & 1  \\ \end{array}} \right]\left\{ {\begin{array}{*{20}c} {u_4 }  \\  {u_5 }  \\ \end{array}} \right\} + \frac{6}\left[ {\begin{array}{*{20}c} 2 & 1 \\   1 & 2  \\ \end{array}} \right] \left\{ {\begin{array}{*{20}c} {\dot u_4 } \\ {\dot u_5 } \\ \end{array}} \right\} = \frac{Q A L_e}{2} \left\{ {\begin{array}{*{20}c} {1}  \\   {1}   \\ \end{array}} \right\} + A\left\{ {\begin{array}{*{20}c} {q_4 } \\ { - q_5 } \\ \end{array}} \right\} $$
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }

The global equation (4) on the system assembly is:
 * {| style="width:100%" border="0"

\displaystyle \left[ C \right]\left\{ {\dot U} \right\} + \left[ K \right]\left\{ U \right\} = \left\{ {F_Q } \right\} + \left\{ {F_g } \right\} $$ where
 * style="width:95%" |$$
 * style="width:95%" |$$
 * style= |
 * }
 * {| style="width:100%" border="0"

\displaystyle \left[ K \right] = \frac\left[ {\begin{array}{*{20}c} 1 & { - 1} & 0 & 0 & 0 \\ { - 1} & 2 & { - 1} & 0 & 0 \\ 0 & { - 1} & 2 & { - 1} & 0 \\ 0 & 0 & { - 1} & 2 & { - 1} \\ 0 & 0 & 0 & { - 1} & 1 \\ \end{array}} \right]
 * style="width:95%" |$$
 * style="width:95%" |$$

\displaystyle \left[ C \right] = \frac{6}\left[ {\begin{array}{*{20}c} 2 & 1 & 0 & 0 & 0 \\ 1 & 4 & 1 & 0 & 0 \\ 0 & 1 & 4 & 1 & 0 \\ 0 & 0 & 1 & 4 & 1 \\ 0 & 0 & 0 & 1 & 2 \\ \end{array}} \right] $$
 * style= |
 * }


 * {| style="width:100%" border="0"

\displaystyle \left\{ {F_Q } \right\} = \frac{2}\left\{ {\begin{array}{*{20}c} 1 \\ 2 \\ 1 \\ 0 \\ 0 \\ \end{array}} \right\}
 * style="width:95%" |$$
 * style="width:95%" |$$

\displaystyle \left\{ {F_g } \right\} = A\left\{ {\begin{array}{*{20}c} {q_1 } \\ 0 \\ 0 \\ 0 \\ { - q_5 } \\ \end{array}} \right\} $$
 * style= |
 * }

Appendix 3: Nondimensionalize the thermal problem
To clarify the nondimensionalizing process, we change variables:

\displaystyle U\left( {x,t} \right) = U_0 \Gamma \left( {\xi ,\tau } \right),\,\,\frac{x}{L} = \xi ,\,\,\,\frac{t} = \tau $$

Therefore:

\displaystyle \frac = U_0 \frac = U_0 \frac\frac = \frac{L}\frac $$

\displaystyle \,\frac = \frac{\partial }\left( {\frac} \right) = \frac{\partial }\left( {\frac{L}\frac} \right) = \frac{L}\frac{\partial }\left( {\frac} \right) = \frac{L}\frac{\partial }\left( {\frac} \right)\frac = \frac\frac $$

\displaystyle \frac = U_0 \frac = U_0 \frac\frac = \frac\frac $$

Substituting the above quantities into the heat transfer problem (1):

\displaystyle k\frac\frac + f\left( {x,t} \right) = \rho c\frac\frac $$ or

\displaystyle \frac\frac + \frac = \frac $$

and the final form is:

\displaystyle \beta \frac + \Phi \left( {\xi,\tau} \right) = \frac $$

where

\displaystyle \beta = \frac, \,\,\,\,\,\, \Phi \left( {\xi,\tau} \right) = \frac $$

The initial condition:

\displaystyle U\left( {x,t = 0} \right) = U_0 \Rightarrow \Gamma \left( {\xi ,\tau  = 0} \right) = 1 $$ From

q = - k\frac =  - \frac{L}\frac \Rightarrow  - \frac = \frac = \varphi $$

we have boundary conditions:

\displaystyle \left. { - \frac} \right|_{\xi = 0}  = \frac = \frac = \varphi_{in}, \,\,\,\,\,\, \left. { - \frac} \right|_{\xi = 1}  = \frac = \frac = \varphi_{out} $$