Analytical Determination of the Transfer Function of a Resistance Furnace

V. N. Chudinov

(Presented by the scientific seminar of the Department of Automation and Telemechanics)

One of the main problems that arise in creating automatic control systems is the problem of the mathematical description of the dynamics of the processes occurring in the object being automated. The most accurate information about the dynamic properties of an object is given by a system of differential equations. Partial differential equations describing objects with distributed parameters can often be reduced to transcendental transfer functions by selecting the corresponding boundary conditions [1].

In the analysis and synthesis of automatic control systems with such objects, it is of interest not so much to solve a complicated differential equation as to obtain the transfer function in fractional-rational form. Hence there arises the problem of approximating transcendental transfer functions in an appropriate way.

In the present work an attempt is made to determine analytically the transfer function of a resistance furnace by solving the general heat-conduction equation with the aid of Laplace transforms under boundary conditions of the third kind.

The indicated approach may be recommended for those cases in which the assumptions made below are valid.

As a result of approximation by means of Padé interpolation fractions, a simplified expression for the transfer function has been obtained.

The problem is solved as applied to indirectly heated resistance furnaces having a rectangular cross section; heaters are located in the furnace walls, Fig. 1. The energy supplied to the furnace provides heating of the air entering the working space, of the furnace walls, and of the product; part of the heat escapes into the surrounding atmosphere from the heated outer walls.

Fig. 1. Furnace diagram.
1 — furnace working space; 2 — heater; 3 — furnace wall; 4 — thermal-insulation casing

In the general case the product is heated owing to convective heat transfer from the heated air and to the radiant energy of the heated inner walls of the furnace.

Consider a furnace wall whose length is considerably greater than its width. We shall make the following assumptions. The heaters are uniformly distributed over the volume of the wall, which on the outside is assumed to be insulated (heat losses to the surrounding atmosphere are usually easy to take into account); the temperature gradient of the furnace wall depends only on one coordinate—the wall thickness; the workpiece is assumed, in the thermal sense, to be “thin.” We shall also assume that the gaseous medium inside the furnace is transparent to the radiant energy flux.

Introduce the following notation:

\(Q\)—furnace power, kcal/h;

\(T_{\mathrm{v}}\)—temperature of the air inside the furnace, \({}^{\circ}\mathrm{K}\);

\(T_{\mathrm{o}}\)—ambient temperature (outside temperature), \({}^{\circ}\mathrm{K}\);

\(T_{\mathrm{c}}\)—temperature of the furnace wall, \({}^{\circ}\mathrm{K}\);

\(T_{\mathrm{n}}\)—temperature of the workpiece, \({}^{\circ}\mathrm{K}\);

\(a\)—coefficient of thermal diffusivity, m;

\(\alpha_{\mathrm{c}}\)—coefficient of heat transfer from the inner walls of the furnace to the air passing inside the furnace,
\[ \frac{\mathrm{kcal}}{\mathrm{m}^{2}\cdot \mathrm{h}\cdot \mathrm{deg}}; \]

\(\alpha_{\mathrm{n}}\)—coefficient of heat transfer from the heated air to the workpiece,
\[ \frac{\mathrm{kcal}}{\mathrm{m}^{2}\cdot \mathrm{h}\cdot \mathrm{deg}}; \]

\(\lambda\)—coefficient of thermal conductivity of the furnace walls,
\[ \frac{\mathrm{kcal}}{\mathrm{m}\cdot \mathrm{h}\cdot \mathrm{deg}}; \]

\(C_{0}\)—emissive coefficient of an absolutely black body,
\[ \frac{\mathrm{kcal}}{\mathrm{m}^{2}\cdot \mathrm{h}\cdot \mathrm{K}^{4}}; \]

\(\varphi\)—irradiation coefficient;

\(\varepsilon_{n}^{\mathrm{pr}}\)—reduced blackness coefficient;

\(F_{\mathrm{c}}, F_{\mathrm{n}}\)—surface area of the wall and of the workpiece, respectively, \(\mathrm{m}^{2}\);

\(C_{\mathrm{n}}, C_{\mathrm{v}}\)—heat capacity of the workpiece material and of the air,
\[ \frac{\mathrm{kcal}}{\mathrm{kg}\cdot \mathrm{deg}}; \]

\(G_{0}\)—furnace throughput;

\(x\)—coordinate along the wall thickness;

\(\tau\)—current time;

\(\delta\)—wall thickness, m.

We regard the furnace wall as an infinite plate with internal heat sources, insulated on one side. Since a one-dimensional problem is considered, the heat-conduction equation for the wall has the form [2]

\[ \frac{\partial T_{\mathrm{c}}(x,\tau)}{\partial \tau} = a\frac{\partial^{2}T_{\mathrm{c}}(x,\tau)}{\partial x^{2}} + \frac{Q_{0}}{c\cdot \gamma}; \tag{1} \]

\[ T_{\mathrm{c}}(x,0)=T_{0};\qquad \frac{\partial T_{\mathrm{c}}(0,\tau)}{\partial x}=0; \qquad (\tau>0;\;0<x<\delta). \]

Let us write the boundary conditions for the steady state, bearing in mind that the quantity of heat supplied by thermal conductivity is expended in heating the air in the furnace chamber. In addition, the workpiece located in the furnace is heated by convective heat transfer from the heated air and by radiant energy from the heated inner walls of the furnace. Thus, the following equations hold:

\[ -\lambda F_{\mathrm{c}}\frac{\partial T_{\mathrm{c}}(x,\tau)}{\partial x} = \alpha_{\mathrm{c}}F_{\mathrm{c}}\,[T_{\mathrm{c}}(x,\tau)-T_{\mathrm{v}}(\tau)] + \alpha_{\mathrm{n}}F_{\mathrm{n}}\,[T_{\mathrm{v}}(\tau)-T_{\mathrm{n}}(\tau)] + \]

\[ + \varphi \varepsilon_n C_0 F_i \left\{ \left[\frac{T_c(x,\tau)}{100}\right]^4 - \left[\frac{T_i(\tau)}{100}\right]^4 \right\}, \]

\[ C_i G_0 [T_i(\tau)-T_0] = \alpha_i F_i [T_v(\tau)-T_i(\tau)] + \varphi \varepsilon_n C_0 F_i \left\{ \left[\frac{T_c(x,\tau)}{100}\right]^4 - \left[\frac{T_i(\tau)}{100}\right]^4 \right\}, \]

\[ C_v m_v [T_v(\tau)-T_0] = \chi_c F_c [T_c(x,\tau)-T_v(\tau)] - \alpha_i F_i [T_v(\tau)-T_i(\tau)]. \]

If the power supplied to the furnace is increased, then the temperature of the bodies located in the furnace will begin to rise, i.e., for \(Q=Q_0+\Delta Q\) the equations hold:

\[ \begin{aligned} T_c(x,\tau) &= T_{c0}+\Delta T_c(x,\Delta),\\ T_i(\tau) &= T_{i0}+\Delta T_i(\tau),\\ T_v(\tau) &= T_{v0}+\Delta T_v(\tau), \end{aligned} \tag{2} \]

where \(T_{c0}, T_{i0}, T_{v0}\) are the temperature values corresponding to the steady state.

Substituting (2) into equation (1) and into the expressions for the boundary conditions, discarding terms containing increments of degree higher than the first, and subtracting the statics equations from the equations obtained, we obtain:

for the wall

\[ \frac{\partial \Delta T_c(x,\tau)}{\partial \tau} = a\frac{\partial^2 \Delta T_c(x,\tau)}{\partial x^2} + \frac{\Delta Q}{c\gamma}. \]

Denoting

\[ \frac{\Delta T_c(x,\tau)}{T_{c0}}=\Theta_c(x,\tau), \]

we obtain

\[ \frac{\partial \Theta_c(x,\tau)}{\partial \tau} = a\frac{\partial^2 \Theta_c(x,\tau)}{\partial x^2} + \frac{\Delta Q}{T_{c0}\,c\gamma}. \]

After applying the Laplace transform to the last equation, we obtain

\[ \overline{\Theta}_c''(x,p) - \frac{p}{a}\overline{\Theta}_c(x,p) + \frac{T_0}{a} + \frac{\Delta Q}{pT_{c0}c\gamma a} = 0, \]

where \(\Theta_c(x,\tau)\doteq \overline{\Theta}_c(x,p)\), and primes denote differentiation with respect to the coordinate.

The solution for the transform may be written [2] in the form

\[ \overline{\Theta}_c(x,p) - \frac{T_0}{p} = \frac{\Delta Q}{p^2 T_{c0}\,c\gamma} + A\,\operatorname{ch}\sqrt{\frac{p}{a}}\,x. \tag{3} \]

Boundary condition \((x=\delta)\).

\[ - \frac{ \partial \dfrac{\Delta T_c(x,\tau)}{T_{c0}} }{ \partial x } = \frac{\Delta T_c(x,\tau)}{T_{c0}} \left[ \frac{\chi_c F_e}{\lambda F_c} + \frac{4T_{c0}^3\varphi\varepsilon_n C_0 F_i}{100^4\lambda F_c} \right] + \]

\[ + \frac{\Delta T_v(\tau)}{T_{v0}} \left[ \frac{\chi_i F_i T_{v0}}{T_{c0}\lambda F_c} - \frac{T_{v0}\chi_c}{T_{c0}\lambda} \right] - \]

\[ - \frac{\Delta T_i(\tau)}{T_{i0}} \left[ \frac{\alpha_i F_i T_{i0}}{T_{c0}\lambda F_c} + \frac{4T_{i0}^4\varphi\varepsilon_n C_0 F_i}{\lambda F_c T_{c0}\cdot 100^4} \right]. \tag{4} \]

Introduce the notation:

\[ \frac{\Delta T_v(z)}{T_{v0}}=\Theta_v(z),\qquad \frac{\Delta T_n(z)}{T_{n0}}=\Theta_n(z); \]

\[ \beta=\frac{100^4\alpha_c F_c+4T_{c0}^3\varphi\varepsilon_n C_0F_n}{100^4\lambda F_c}; \]

\[ \beta_1=\frac{T_{v0}}{T_{c0}}\cdot \frac{\alpha_n F_n-\alpha_cF_c}{\lambda F_c}; \]

\[ \beta_2=\frac{T_{n0}}{T_{c0}}\cdot \frac{100^4\alpha_nF_n+4T_{n0}^3\varphi\varepsilon_nC_0F_n}{100^4\lambda F_c}. \]

Taking the new notation into account, equation (4), after the Laplace transform, takes the form

\[ \frac{\partial \overline{\Theta}_c(x,p)}{\partial x} +\beta\overline{\Theta}_c(x,p) +\beta_1\overline{\Theta}_v(p) +\beta_2\overline{\Theta}_n(p)=0. \tag{5} \]

The equation for heating the product:

\[ \Delta T_n(z)\cdot \left( C_nG_0+\alpha_nF_n+ \frac{4\varphi\varepsilon_nC_0F_nT_{n0}^3}{100^4} \right) - \]

\[ -\alpha_nF_n\Delta T_v(z) -\frac{\Delta T_c(x,z)}{T_{c0}}\cdot \frac{4\varphi\varepsilon_nC_0F_nT_{c0}^4}{100^4} =0. \tag{6} \]

Introduce the notation:

\[ \beta'=\frac{4\varphi\varepsilon_nC_0F_nT_{c0}^3}{100^4}, \qquad \beta'_1=\frac{T_{v0}}{T_{c0}}\cdot\alpha_nF_n, \]

\[ \beta'_2=\frac{T_{n0}}{100^4T_{c0}} \left(C_nG_0+100^4\alpha_nF_n+4\varphi\varepsilon_nC_0F_nT_{n0}^3\right). \]

Taking into account the introduced notation, after applying the Laplace transform to equation (6), we obtain

\[ \beta'\overline{\Theta}_c(x,p) +\beta'_1\overline{\Theta}_v(p) -\beta'_2\overline{\Theta}_n(p)=0. \tag{7} \]

The equation for heating the air:

\[ \alpha_cF_c\Delta T_c(x,z) -\Delta T_v(z)\left(C_vm_v+\alpha_cF_c+\alpha_nF_n\right) +\Delta T_n(z)\alpha_nF_n=0. \]

Taking into account the notation:

\[ \beta''_1=\frac{T_{v0}}{T_{c0}}\cdot \frac{C_vm_v+\alpha_cF_c+\alpha_nF_n}{\alpha_cF_c}, \]

\[ \beta''_2=\frac{T_{n0}}{T_{c0}}\cdot \frac{\alpha_nF_n}{\alpha_cF_c}, \]

applying the Laplace transform, we rewrite equation (8) as

\[ \overline{\Theta}_c(x,p) -\beta''_1\overline{\Theta}_v(p) +\beta''_2\overline{\Theta}_n(p)=0. \tag{9} \]

Solving (3), (5), (7), and (9) simultaneously with respect to the variables of interest, we obtain

\[ \overline{\Theta}_n(p) \left( \beta^{**}\beta +\beta^{**}\sqrt{\frac{p}{a}}\operatorname{th}\sqrt{\frac{p}{a}}\cdot\delta -\beta^{**}\frac{\beta_2\beta'}{\beta'_2} + \right. \]

\[ + \frac{\beta_1 \beta_2' - \beta_2 - \beta_1'}{\beta_2'} \right) - \frac{T_0}{p}\sqrt{\frac{p}{a}}\,\operatorname{th}\sqrt{\frac{p}{a}}\,\delta - \frac{\Delta Q}{p^2 T_{c0} c \gamma} \sqrt{\frac{p}{a}}\,\operatorname{th}\sqrt{\frac{p}{a}}\,\delta = 0, \tag{10} \]

where

\[ \beta^* = \frac{\beta_1''\beta_2' - \beta_2''\beta_1'} {\beta_1' + \beta_2'\beta'} . \]

The transfer function \(W(p)\) is the operator of the dynamic relation between the input quantity (the power supplied to the furnace) and the output quantity (the air temperature in the furnace chamber). According to (10), we write:

\[ W(p) = \frac{ A\,\operatorname{th}\sqrt{\frac{p}{a}}\,\delta }{ p\,\operatorname{th}\sqrt{\frac{p}{a}}\,\delta + B\sqrt{p} } = \]

\[ = \frac{ A\left(\exp 2\sqrt{\frac{p}{a}}\,\delta - 1\right) }{ p\left(\exp 2\sqrt{\frac{p}{a}}\,\delta - 1\right) + B\sqrt{p}\left(\exp 2\sqrt{\frac{p}{a}}\,\delta + 1\right) }, \tag{11} \]

where

\[ A = \frac{\Delta Q}{T_{c0} c \gamma \beta^*}, \qquad B = \frac{ \sqrt{a}\,\beta\beta_2'\beta^* - \beta_2\beta'\beta^* + \beta_1\beta_2' - \beta_2\beta_1' }{ \beta^*\beta_2' }. \]

The transfer function (11) of the furnace is a transcendental expression. Further operations with such an expression are rather difficult. Therefore it is necessary to replace the exact expression (11) by an approximating transfer function in the form of polynomials in the operator \(p\).

In [3] a method is described for determining an approximating transfer function from the known coefficients of the expansion of the original transfer function into an infinite series. However, the indicated method does not always lead to a satisfactory result.

Another method of approximation is based on the use of Padé fractions [4] for monotone bounded processes. In [4] the concept of the inertance \(S\) of a transient monotone process is used, introduced for determining the coefficients of transfer functions and being a generalization of the concept of the process time constant. The expression for \(S_n\) has the form [4]:

\[ S_n = (-1)^n \frac{1}{n!}\lim_{p\to 0}\frac{d^n W(p)}{dp^n}. \]

The approximation of the furnace transfer function by Padé fractions, expressed through inertances of the corresponding orders, has the form

\[ W(p) = \frac{S_0}{ 1 + \frac{S_1}{S_0}p + \left[ \left(\frac{S_1}{S_0}\right)^2 - \frac{S_2}{S_0} \right]p^2 }, \tag{12} \]

where

\[ S_0 = \frac{A\delta}{Ba}, \qquad S_1 = \frac{A\delta^2(4+B\delta)}{4B^2a^3}, \qquad S_2 = \frac{A\,2\delta^3(16+B\delta)}{96B^4a^4}. \]

Using the indicated method, an approximate expression of a complex transcendental function by a fractional-rational one is achieved, which corresponds to approximating the curve of the transient process by the integral curve of a certain homogeneous differential equation with constant coefficients.

Fig. 2.

The principal result of applying the method [4] to the problem of analytical approximation of a transcendental transfer function is that the coefficients of the simplified transfer function could be expressed through the design parameters of the furnace. The theoretical results presented in the work agree well with the experimental data.

To confirm the theoretical result, the transfer function of an M-24 enamel furnace was determined. The expression for the transient function, found from the experimentally recorded acceleration curve by Simoiu’s method [3], gives the following dependence:

\[ \bar g[t] = 1 - 1.05 e^{-\frac{t}{2.49}} + 0.06 e^{-\frac{t}{0.133}} . \]

The result of approximating the analytical expression (11) by means of Padé interpolation fractions shows that

\[ \bar g(t) = 1 - 1.47 e^{-\frac{t}{1.76}} + 0.466 e^{-\frac{t}{0.46}} . \]

A comparison of the results of the theoretical analysis with the experimental data is presented in Fig. 2.

It should be noted that the described method makes it possible successfully to calculate the dynamic characteristics of furnaces with good thermal insulation, or when heat losses to the surrounding medium can be readily taken into account.

REFERENCES

1. I. A. Feldman, G. K. Rubin, Z. L. Treizon. Investigation of the thermal inertia of electric resistance furnaces. Transactions of VNIIETO, issue 1, 1965.

2. A. V. Lykov. Theory of Heat Conduction. State Publishing House of Technical-Theoretical Literature, 1952.

3. M. P. Simoiu. Determination of the coefficients of transfer functions of linearized elements and automatic-control systems. Avtomatika i telemekhanika, 18, No. 6, 1957.

4. B. N. Devyatov. Theory of Transient Processes in Technological Apparatuses from the Standpoint of Control Problems. Siberian Branch of the USSR Academy of Sciences, Novosibirsk, 1964.