# A boundary value problem of heat transfer within DBD-based plasma jet setups

Panayiotis Vafeas, Anastassios Skarlatos, Policarpos K. Papadopoulos, Panayiotis Svarnas, Nikolaos Sarmas

► To cite this version:

Panayiotis Vafeas, Anastassios Skarlatos, Policarpos K. Papadopoulos, Panayiotis Svarnas, Nikolaos Sarmas. A boundary value problem of heat transfer within DBD-based plasma jet setups. Mathematical Biosciences and Engineering, 2023, 20 (10), pp.18345-18367. 10.3934/mbe.2023815. cea-04239500

HAL Id: cea-04239500

https://cea.hal.science/cea-04239500

Submitted on 12 Oct 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Research article

# A boundary value problem of heat transfer within DBD-based plasma jet setups

P. Vafeas<fnref n="1" />*, A. Skarlatos<fnref n="2" />, P. K. Papadopoulos<fnref n="3" />, P. Svarnas<fnref n="4" /> and N. Sarmas<fnref n="3" />

**Abstract:** We claim an analytical solution for the thermal boundary value problem that arises in DBD-based plasma jet systems as a preliminary and consistent approach to a simplified geometry. This approach involves the outline of a coaxial plasma jet reactor and the consideration of the heat transfer to the reactor solids, namely, the dielectric barrier and the grounded electrode. The non-homogeneous initial and boundary value thermal problem is solved analytically, while a simple cut-off technique is applied to deal with the appearance of infinite series relationships, being the outcome of merging dual expressions. The results are also implemented numerically, supporting the analytical solution, while a Finite Integration Technique (FIT) is used for the validation. Both the analytical and numerical data reveal the temperature pattern at the cross-section of the solids in perfect agreement. This analytical approach could be of importance for the optimization of plasma jet systems employed in tailored applications where temperature-sensitive materials are involved, like in plasma biomedicine.

**Keywords:** dielectric barrier discharge; cold plasma jet systems; thermal non-homogeneous boundary value problem; analytical methodology; finite integration technique

## 1. Introduction

Cold atmospheric pressure plasmas (CAPPs) refer to the ionized gaseous phase, being out of

<footnote n="1">¹ Department of Chemical Engineering, University of Patras, 26504 Patras, Greece</footnote>

<footnote n="2">² Université Paris-Saclay, CEA, List, 91120 Palaiseau, France</footnote>

<footnote n="3">³ Department of Mechanical Engineering & Aeronautics, University of Patras, 26504 Patras, Greece</footnote>

<footnote n="4">⁴ High Voltage Laboratory, Department of Electrical & Computer Engineering, University of Patras, 26504 Patras, Greece</footnote>

<footnote>* Correspondence: Email: vafeas@chemeng.upatras.gr.</footnote>
thermodynamic equilibrium and thus rich in reactive species, which is sustained close to the ambient pressure. On the other hand, dielectric barrier discharges (DBDs) prevent transition to arcs and thus to electrical discharge regimes associated with raised temperatures. Hence, unique physicochemical features may be combined by means of DBD-based CAPPs, and this concept has led to cutting-edge technologies in the fields of material processing, aeronautics [1], environmental remediation [2], etc. Of special importance is the plasma biomedicine field, where DBD-based CAPPs are explored as tools for wound healing, sterilization [3], cancer treatment [4], tooth bleaching, etc.

Towards such biomedical applications, DBD-based CAPPs are usually engineered in the form of the so-called “plasma jets” [5], where a gas is channeled into the atmospheric air through a dielectric capillary tube, while it is simultaneously subjected to a time-varying electric field. The latter is developed by driving a pair of electrodes of appropriate geometry (e.g., coaxial) with high voltages of various waveforms (e.g., ac-continuous or ns-pulsed in the kHz range). Thus, a complicated interplay between charged species, excited neutrals, chemical radicals, photons, electric fields [6], fluid fields [7–10], etc. and the bio-specimen itself, takes place. However, the benefits of the interaction between this plasma reactive phase and the specimen may be compensated by the rise of the gas temperature with respect to the specimen’s thermal tolerance. In other words, dominant design and operation parameters that are considered during the optimization of the plasma jet setups must be correlated with the thermal fields of various heating mechanisms [11].

Several works have been devoted to the thermal properties of CAPPs, whereas they mostly refer to experimental and numerical techniques [12–20]. Herein, the interest in investigating the thermal properties of plasma jets, as a special case of CAPPs, is combined with the challenge of producing an initial, boundary value problem of heat conduction within such setups. The ultimate motivation for this research is to obtain analytical formulae for the temperature field in the system that can be employed in plasma simulations, instead of time-consuming numerical solutions. These solutions will provide an easy means of estimating the spatial distribution of the temperature for different geometric and operation characteristics, which are important factors for the properties of the plasma. Furthermore, novel physical results can be derived from analytical techniques and they can be utilized for the verification of numerical results. Under these aspects, the present article has endeavored to show that an analytical approach to problems involving CAPP jet setups is conceivable.

Namely, the theoretical analysis considers a two-region problem, which corresponds to the dielectric barrier tube of the DBD and the external grounded electrode, wherein our goal is to calculate the temperature distribution in the axisymmetric configuration of the setup. In order to construct the related initial and boundary value problem and proceed to the solution, our analytical methodology is primarily based on the selection of a suitably located cylindrical coordinate system [21] for modeling purposes with respect to the current investigation. Therein, the domain of field activity is divided into the two subsectors of the dielectric and the electrode-metal, in which classical heat diffusion partial differential equations [22] are considered. These are supplemented by the appropriate conditions, i.e., the Dirichlet and Neumann continuity conditions on the interface of the two solids, an imposed given heat conduction rate assumed on the boundary of the dielectric with the working gas, and the Robin-type heat convection condition on the boundary of the electrode with the outer environment. Otherwise, an initial room temperature condition is imposed for both the areas of thermal activity, completing the well-pose of the problem. Henceforth, the well-known method of separation of variables is applied to both the heat equations, while a standard technique, based on the
method of asymptotic kernels [23], is necessary for handling the non-homogeneous boundary conditions. Following these steps, the time-dependent temperature fields in each region are expanded in terms of cylindrical eigenfunctions [23] and decaying time exponents, in which rotational symmetry is implied, featured by the conditions. Then, the aforementioned problem is solved straightforwardly, resulting in either analytical or semi-analytical expressions for the incorporated unknown constant coefficients, while a pair of dual relations of Fourier-Bessel series [24] for a single unknown also appear due to the different boundaries. The latter are merged to one and only relationship, which is handled appropriately so as to be reduced to an infinite algebraic linear system of equations for the implicated set of constants, by taking its projection to Bessel and Neumann functions. Therein, the obtained handy system is manipulated via well-known cut-off techniques, while the temperature fields are provided in a compact fashion via closed-type convergent series expansions and they are implemented so as to show their response graphically. The results reveal the thermal behavior both in the bulk and on the surfaces of the solid materials, which is of interest for future elaboration either using alternative approaches, e.g., by solving in the Laplace domain and transforming back to the time domain [25,26] or proceeding to possible experimental processing. In order to demonstrate the efficiency of the present method, the analytical solution is juxtaposed with numerical results, being obtained using the advanced Finite Integration Technique (FIT) [27–29], wherein an excellent agreement is obtained.

## 2. Heat Transfer Configuration and Analytical Solution

Let us consider two very long coaxial cylinders, the dielectric tube (alumina) of thermal diffusivity $a_d$ with inner radius $R_d$ and the electrode (brass) of thermal diffusivity $a_e$ with inner radius $R_e$. Furthermore, $k_d$ and $k_e$ are the respective thermal conductivities (see Figure 1). The system is initially set at the same temperature with the free-environment, whose convective heat transfer coefficient is $h$, contains air at room temperature $T_∞$, while plasma kinetics produces a known conduction heat $q$. Hence, in terms of a properly adjusted cylindrical coordinate system $(\rho, \varphi, z)$, in which the $z$-variable is eliminated as a fair approximation of the big length of the plasma jet setup, we consider two different temperature distributions, the first one within the dielectric $T_d(\rho, \varphi, t)$ and the second one within the electrode $T_e(\rho, \varphi, t)$. Then the boundary and initial value problem that we have to solve admits two separate heat partial differential equations in cylindrical geometry, which assume

$$ \alpha_d \left( \frac{\partial^2}{\partial \rho^2} + \frac{1}{\rho} \frac{\partial}{\partial \rho} + \frac{1}{\rho^2} \frac{\partial^2}{\partial \varphi^2} \right) T_d(\rho, \varphi, t) = \frac{\partial T_d(\rho, \varphi, t)}{\partial t} \quad \text{at} \quad \rho \in (R_d, R_e) \quad (1) $$

and

$$ \alpha_e \left( \frac{\partial^2}{\partial \rho^2} + \frac{1}{\rho} \frac{\partial}{\partial \rho} + \frac{1}{\rho^2} \frac{\partial^2}{\partial \varphi^2} \right) T_e(\rho, \varphi, t) = \frac{\partial T_e(\rho, \varphi, t)}{\partial t} \quad \text{at} \quad \rho \in (R_e, R_s), \quad (2) $$
both defined for every $\varphi \in [0,2\pi)$ and $t > 0$, accompanied by the appropriate boundary conditions

$$-k_d \frac{\partial T_d(R_d, \varphi, t)}{\partial \rho} = q \quad \text{for } \varphi \in [0,2\pi) \text{ and } t \ge 0, \qquad (3)$$

$$T_d(R_e, \varphi, t) = T_e(R_e, \varphi, t) \quad \text{for } \varphi \in [0,2\pi) \text{ and } t \ge 0, \qquad (4)$$

$$-k_d \frac{\partial T_d(R_e, \varphi, t)}{\partial \rho} = -k_e \frac{\partial T_e(R_e, \varphi, t)}{\partial \rho} \quad \text{for } \varphi \in [0,2\pi) \text{ and } t \ge 0, \qquad (5)$$

$$-k_e \frac{\partial T_e(R_s, \varphi, t)}{\partial \rho} = h[T_e(R_s, \varphi, t) - T_{\infty}] \quad \text{for } \varphi \in [0,2\pi) \text{ and } t \ge 0 \qquad (6)$$

and the initial condition of common temperature for the dielectric tube and the electrode

$$T_d(\rho, 0) = T_{\infty} \quad \text{for } \rho \in [R_d, R_e] \text{ and } T_e(\rho, 0) = T_{\infty} \quad \text{for } \rho \in [R_e, R_s], \qquad (7)$$

completing thus the set of a well-posed boundary value problem that describes the heat transfer within our system.

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Simplified plasma reactor geometry and the areas of thermal activity, i.e., the dielectric tube with inner and outer radius $R_d$ and $R_e$, respectively and the grounded electrode with inner and outer radius $R_e$ and $R_s$, respectively.</figcaption></figure>

Our first task in order to solve (1)–(7), is to transform (3) and (6) into homogeneous conditions. To this end, we use the method of asymptotic kernels, according to which the solution under consideration is rewritten as the summation of a general solution and a particular time-independent solution (resembling the asymptotic behavior as $t \to +\infty$) via

$$T_d(\rho, \varphi, t) = \tilde{T}_d(\rho, \varphi, t) + K_d(\rho, \varphi) \quad \text{for every } \rho \in (R_d, R_e), \varphi \in [0, 2\pi) \text{ and } t > 0, \qquad (8)$$

and

$$T_e(\rho, \varphi, t) = \tilde{T}_e(\rho, \varphi, t) + K_e(\rho, \varphi) \quad \text{for every } \rho \in (R_e, R_s), \varphi \in [0, 2\pi) \text{ and } t > 0. \qquad (9)$$

According to our argumentation, functions $K_d(\rho, \varphi)$ and $K_e(\rho, \varphi)$ are chosen, such as to solve the problems

$$\left( \frac{\partial^2}{\partial \rho^2} + \frac{1}{\rho} \frac{\partial}{\partial \rho} + \frac{1}{\rho^2} \frac{\partial^2}{\partial \varphi^2} \right) K_d(\rho, \varphi) = 0 \quad \text{for } \rho \in (R_d, R_e) \text{ and } \varphi \in [0, 2\pi) \qquad (10)$$
with boundary condition

$$-k_d \frac{\partial K_d(R_d, \varphi)}{\partial \rho} = q \quad \text{for } \varphi \in [0, 2\pi], \qquad (11)$$

while

$$\left( \frac{\partial^2}{\partial \rho^2} + \frac{1}{\rho} \frac{\partial}{\partial \rho} + \frac{1}{\rho^2} \frac{\partial^2}{\partial \varphi^2} \right) K_e(\rho, \varphi) = 0 \quad \text{for } \rho \in (R_e, R_s) \text{ and } \varphi \in [0, 2\pi] \qquad (12)$$

with boundary condition

$$-k_e \frac{\partial K_e(R_s, \varphi)}{\partial \rho} = h[K_e(R_s, \varphi) - T_{\infty}] \quad \text{for } \varphi \in [0, 2\pi]. \qquad (13)$$

On the other hand, the main fields $\tilde{T}_d(\rho, \varphi, t)$ and $\tilde{T}_e(\rho, \varphi, t)$, implying (10) and (12) into relations (1)–(7) and according to decompositions (8) and (9), satisfy the boundary value problems for every $\varphi \in [0, 2\pi)$ and $t > 0$

$$\alpha_d \left( \frac{\partial^2}{\partial \rho^2} + \frac{1}{\rho} \frac{\partial}{\partial \rho} + \frac{1}{\rho^2} \frac{\partial^2}{\partial \varphi^2} \right) \tilde{T}_d(\rho, \varphi, t) = \frac{\partial \tilde{T}_d(\rho, \varphi, t)}{\partial t} \quad \text{at } \rho \in (R_d, R_e) \qquad (14)$$

and

$$\alpha_e \left( \frac{\partial^2}{\partial \rho^2} + \frac{1}{\rho} \frac{\partial}{\partial \rho} + \frac{1}{\rho^2} \frac{\partial^2}{\partial \varphi^2} \right) \tilde{T}_e(\rho, \varphi, t) = \frac{\partial \tilde{T}_e(\rho, \varphi, t)}{\partial t} \quad \text{at } \rho \in (R_e, R_s) \qquad (15)$$

with boundary conditions

$$\frac{\partial \tilde{T}_d(R_d, \varphi, t)}{\partial \rho} = 0 \quad \text{for } \varphi \in [0, 2\pi) \text{ and } t \ge 0, \qquad (16)$$

$$\tilde{T}_d(R_e, \varphi, t) + K_d(R_e, \varphi) = \tilde{T}_e(R_e, \varphi, t) + K_e(R_e, \varphi) \quad \text{for } \varphi \in [0, 2\pi) \text{ and } t \ge 0, \qquad (17)$$

$$-k_d \frac{\partial[\tilde{T}_d(R_e, \varphi, t) + K_d(R_e, \varphi)]}{\partial \rho} = -k_e \frac{\partial[\tilde{T}_e(R_e, \varphi, t) + K_e(R_e, \varphi)]}{\partial \rho}, \quad \varphi \in [0, 2\pi), \ t \ge 0, \qquad (18)$$

$$k_e \frac{\partial \tilde{T}_e(R_s, \varphi, t)}{\partial \rho} + h \tilde{T}_e(R_s, \varphi, t) = 0, \quad \varphi \in [0, 2\pi), \ t \ge 0 \qquad (19)$$

and initial condition of the system

$$\tilde{T}_d(\rho, \varphi, 0) + K_d(\rho, \varphi) = T_{\infty} \quad \text{for } \rho \in [R_d, R_e] \text{ and } \tilde{T}_e(\rho, \varphi, 0) + K_e(\rho, \varphi) = T_{\infty} \quad \text{for } \rho \in [R_e, R_s], \qquad (20)$$

whereas the common temperature is set by the final solution.

In the sequel, we work as follows. Due to conditions (11) and (13), functions $K_d(\rho, \varphi)$ and $K_e(\rho, \varphi)$ depend only on the $\rho$-variable, i.e., $K_d(\rho, \varphi) \equiv K_d(\rho)$ and $K_e(\rho, \varphi) \equiv K_e(\rho)$, respectively, henceforth (10) and (12) become eventually ordinary differential equations, since $\partial K_d(\rho)/\partial \varphi = 0$
and $\partial K_e(\rho) / \partial \varphi = 0$, whereas applying (11) and (13) independently of $\varphi$-variable, it is readily obtained

$$K_d(\rho) = -q \frac{R_d}{k_d} \ln \rho + c_d \quad \text{for } \rho \in [R_d, R_e] \qquad (21)$$

and

$$K_e(\rho) = T_\infty + c_e \left( \ln \frac{\rho}{R_s} - \frac{k_e}{hR_s} \right) \quad \text{for } \rho \in (R_e, R_s], \qquad (22)$$

where $c_d, c_e \in \mathbb{R}$ are arbitrary constants, appropriately chosen within the forthcoming steps. On the other hand, in order to manipulate the fields $\tilde{T}_d(\rho, \varphi, t)$ and $\tilde{T}_e(\rho, \varphi, t)$ into the partial differential equations (14) and (15), respectively, we apply the well-known method of separation of variables, according to which we suppose a solution of the general form

$$\tilde{T}_j(\rho, \varphi, t) = P(\rho)\Phi(\varphi)T(t) \quad \text{with } j=d \text{ if } \rho \in (R_d, R_e) \text{ and } j=e \text{ if } \rho \in (R_e, R_s) \qquad (23)$$

for $\varphi \in [0, 2\pi)$ and $t > 0$, wherein, since (14) and (15) retain the same structure, it is doable to perform the procedure once. Hence incorporating two separation constants $\kappa$ and $\lambda$, we substitute (23) into either (14) or (15) to get

$$\frac{T'(t)}{\alpha_j T(t)} = \frac{P''(\rho) + \frac{1}{r} P'(\rho)}{P(\rho)} + \frac{1}{\rho^2} \frac{\Phi''(\varphi)}{\Phi(\varphi)} = \kappa \quad \text{with } j=d,e, \qquad (24)$$

accordingly to the domain of interest and

$$\frac{\rho^2 P''(\rho) + \rho P'(\rho)}{P(\rho)} = -\frac{\Phi''(\varphi)}{\Phi(\varphi)} = \lambda, \qquad (25)$$

whereas either $\rho \in (R_d, R_e)$ ($j=d$) or $\rho \in (R_e, R_s)$ ($j=e$). Solving the second equality of (25) and demanding the proper periodic azimuthal behavior, we are led to

$$\Phi(\varphi) \equiv \Phi_n(\varphi) = \begin{cases} \cos(n\varphi) & , \text{where } \lambda \equiv \lambda_n = n^2 \text{ for } n \ge 0 \text{ and } \varphi \in [0, 2\pi), \\ \sin(n\varphi) & , \end{cases} \qquad (26)$$

whilst by combination of (26) with (24) for $\kappa = -\mu_j^2$ with $j=d,e$, we obtain the Bessel-type ordinary differential equation

$$P''(\rho) + \frac{1}{\rho} P'(\rho) + \left( \mu_j^2 - \frac{n^2}{\rho^2} \right) P(\rho) = 0, \qquad (27)$$

whereas, again, either $\rho \in (R_d, R_e)$ ($j=d$) or $\rho \in (R_e, R_s)$ ($j=e$), whose solutions are the Bessel functions $J_n$ and the Neumann functions $N_n$, that is
$$P(\rho) \equiv P_n(\rho) = \begin{cases} J_n(\mu_j \rho) \\ N_n(\mu_j \rho) \end{cases} \quad \text{with } j=d \text{ if } \rho \in (R_d, R_e) \text{ and } j=e \text{ if } \rho \in (R_e, R_s). \qquad (28)$$

Finally, the time-dependence for $\kappa = -\mu_j^2$ with $j=d,e$ within (24) yields

$$T'(t) + \alpha_j \mu_j^2 T(t) = 0 \Rightarrow T(t) \equiv T_{\mu_j}(t) = e^{-\alpha_j \mu_j^2 t} \quad \text{with } j=d,e. \qquad (29)$$

Gathering all the formulae obtained above through (26), (28) and (29), the general solutions are provided via

$$\tilde{T}_j(\rho, \varphi, t) = \sum_{n=0}^{\infty} \int \sum_{\mu_j} \left[ a_n^{\mu_j} J_n(\mu_j \rho) + b_n^{\mu_j} N_n(\mu_j \rho) \right] \left[ c_n \cos(n\varphi) + d_n \sin(n\varphi) \right] e^{-\alpha_j \mu_j^2 t}, \qquad (30)$$

for every $\varphi \in [0, 2\pi)$ and $t > 0$, where $j=d$ if $\rho \in (R_d, R_e)$ and $j=e$ if $\rho \in (R_e, R_s)$. The symbolic feature “$\int \sum_{\mu_j} \dots$” corresponds either to series definition in the case where $\mu_j$ take discrete values or to integral notation if $\mu_j$ are continuous parameters for $j=d,e$. Due to the nature of conditions (17) and (18), the azimuthal dependence must be omitted by imposing $c_n = d_n = 0$ for $n \ge 1$ and keeping the terms for $n=0$. Consequently, by the new definition of $A_0^{\mu_j} \equiv a_0^{\mu_j} c_0$ and $B_0^{\mu_j} \equiv b_0^{\mu_j} c_0$, the general solution (30), which corresponds to the solution of (14) and (15), reduces to

$$\tilde{T}_d(\rho, t) = \int \sum_{\mu_d} \left[ A_0^{\mu_d} J_0(\mu_d \rho) + B_0^{\mu_d} N_0(\mu_d \rho) \right] e^{-\alpha_d \mu_d^2 t} \quad \text{for } \rho \in (R_d, R_e) \text{ and } t > 0 \qquad (31)$$

and

$$\tilde{T}_e(\rho, t) = \int \sum_{\mu_e} \left[ A_0^{\mu_e} J_0(\mu_e \rho) + B_0^{\mu_e} N_0(\mu_e \rho) \right] e^{-\alpha_e \mu_e^2 t} \quad \text{for } \rho \in (R_e, R_s) \text{ and } t > 0, \qquad (32)$$

whereas the set of $A_0^{\mu_d}$, $B_0^{\mu_d}$ with parameter $\mu_d$ and $A_0^{\mu_e}$, $B_0^{\mu_e}$ with parameter $\mu_e$ must be evaluated by (16)–(20) with respect to (21) and (22).

Bearing in mind that $\partial/\partial\rho = \mu_j\partial/\partial(\mu_j\rho)$ and the derivative property $J'_0(\mu_j\rho) = -J_1(\mu_j\rho)$ and $N'_0(\mu_j\rho) = -N_1(\mu_j\rho)$ over the argument for $j=d,e$, putting (16) into (31), we reach

$$B_0^{\mu_d} = -A_0^{\mu_d} \frac{J_1(\mu_d R_d)}{N_1(\mu_d R_d)}, \qquad (33)$$

while applying (19) into (32), we get

$$B_0^{\mu_e} = -A_0^{\mu_e} \frac{[hJ_0(\mu_e R_s) - k_e \mu_e J_1(\mu_e R_s)]}{[hN_0(\mu_e R_s) - k_e \mu_e N_1(\mu_e R_s)]} \qquad (34)$$
for every value of the parameters $\mu_d$ and $\mu_e$, respectively. Substituting (33) into (31) and (34) into (32), we obtain the homogeneous part of the solutions via

$$ \tilde{T}_d(\rho, t) = \int \sum_{\mu_d} C_0^{\mu_d} \left[ \frac{J_0(\mu_d \rho)}{J_1(\mu_d R_d)} - \frac{N_0(\mu_d \rho)}{N_1(\mu_d R_d)} \right] e^{-\alpha_d \mu_d^2 t} \quad (35) $$

for $\rho \in (R_d, R_e)$ and $t > 0$, while

$$ \tilde{T}_e(\rho, t) = \int \sum_{\mu_e} C_0^{\mu_e} \left[ \frac{J_0(\mu_e \rho)}{[hJ_0(\mu_e R_s) - k_e \mu_e J_1(\mu_e R_s)]} - \frac{N_0(\mu_e \rho)}{[hN_0(\mu_e R_s) - k_e \mu_e N_1(\mu_e R_s)]} \right] e^{-\alpha_e \mu_e^2 t} \quad (36) $$

for $\rho \in (R_e, R_s)$ and $t > 0$, wherein $C_0^{\mu_d} = A_0^{\mu_d} J_1(\mu_d R_d)$ and $C_0^{\mu_e} = A_0^{\mu_e} [hJ_0(\mu_e R_s) - k_e \mu_e J_1(\mu_e R_s)]$ are the newly defined constants, respectively. We precede with the interface conditions (17) and (18), which, by virtue of (35) and (36), as well as (21) and (22), yield

$$ \begin{aligned} & \int \sum_{\mu_d} C_0^{\mu_d} \left[ \frac{J_0(\mu_d R_e)}{J_1(\mu_d R_d)} - \frac{N_0(\mu_d R_e)}{N_1(\mu_d R_d)} \right] e^{-\alpha_d \mu_d^2 t} - \left[ q \frac{R_d}{k_d} \ln R_e - c_d + T_\infty + c_e \left( \ln \frac{R_e}{R_s} - \frac{k_e}{hR_s} \right) \right] \\ &= \int \sum_{\mu_e} C_0^{\mu_e} \left[ \frac{J_0(\mu_e R_e)}{[hJ_0(\mu_e R_s) - k_e \mu_e J_1(\mu_e R_s)]} - \frac{N_0(\mu_e R_e)}{[hN_0(\mu_e R_s) - k_e \mu_e N_1(\mu_e R_s)]} \right] e^{-\alpha_e \mu_e^2 t} \end{aligned} \quad (37) $$

and

$$ \begin{aligned} & k_d \int \sum_{\mu_d} \mu_d C_0^{\mu_d} \left[ \frac{J'_0(\mu_d R_e)}{J_1(\mu_d R_d)} - \frac{N'_0(\mu_d R_e)}{N_1(\mu_d R_d)} \right] e^{-\alpha_d \mu_d^2 t} - \left[ q \frac{R_d}{R_e} + c_e \frac{k_e}{R_e} \right] \\ &= k_e \int \sum_{\mu_e} \mu_e C_0^{\mu_e} \left[ \frac{J'_0(\mu_e R_e)}{[hJ_0(\mu_e R_s) - k_e \mu_e J_1(\mu_e R_s)]} - \frac{N'_0(\mu_e R_e)}{[hN_0(\mu_e R_s) - k_e \mu_e N_1(\mu_e R_s)]} \right] e^{-\alpha_e \mu_e^2 t} \end{aligned} \quad (38) $$

with $t \ge 0$, respectively. The utility of relationships (37) and (38) is twofold, since the proper matching for any $t \ge 0$ acquires two things. Firstly, we are forced to consider that

$$ q \frac{R_d}{k_d} \ln R_e - c_d + T_\infty + c_e \left( \ln \frac{R_e}{R_s} - \frac{k_e}{hR_s} \right) = 0 \quad \text{and} \quad q \frac{R_d}{R_e} + c_e \frac{k_e}{R_e} = 0, \quad (39) $$

which immediately offer us

$$ c_e = -q \frac{R_d}{k_e} \quad \text{and} \quad c_d = T_\infty + \frac{qR_d}{k_d} \left( \ln R_e - \frac{k_d}{k_e} \ln \frac{R_e}{R_s} + \frac{k_d}{hR_s} \right), \quad (40) $$

providing the two arbitrary constants in the fields (21) and (22). Once done, implying (39) into (37) and (38), we proceed to the second step, i.e., in order to retain the same exponents in the time variable, we imply

$$ \alpha_d \mu_d^2 = \alpha_e \mu_e^2 \Rightarrow \frac{\mu_d^2}{\mu_e^2} = \frac{\alpha_e}{\alpha_d} \quad \text{or} \quad \mu_d = \mu_e \sqrt{\frac{\alpha_e}{\alpha_d}} \quad (41) $$

and therein, for notational convenience, we define as $\mu \equiv \mu_e$, therefore $\mu_d = \mu \sqrt{\alpha_e / \alpha_d}$. Hence,
following the new parameter notation, we interchange $F_0^\mu \equiv C_0^{\mu_d}$ and $G_0^\mu \equiv C_0^{\mu_e}$ for the same reason and we readily define convenient functions that enter into the fields (35) and (36), as well as into conditions (37) and (38), those being

$$f(\mu\sqrt{\alpha_e / \alpha_d} \rho) = \frac{J_0(\mu\sqrt{\alpha_e / \alpha_d} \rho)}{J_1(\mu\sqrt{\alpha_e / \alpha_d} R_d)} - \frac{N_0(\mu\sqrt{\alpha_e / \alpha_d} \rho)}{N_1(\mu\sqrt{\alpha_e / \alpha_d} R_d)} \quad (42)$$

and

$$g(\mu\rho) = \frac{J_0(\mu\rho)}{[hJ_0(\mu R_s) - k_e\mu J_1(\mu R_s)]} - \frac{N_0(\mu\rho)}{[hN_0(\mu R_s) - k_e\mu N_1(\mu R_s)]} \quad (43)$$

Within the same aspect, we introduce the ratio $C_0^\mu$ of the newly-defined defined constants $C_0^{\mu_d}$ and $C_0^{\mu_e}$ as

$$C_0^\mu \equiv \frac{C_0^{\mu_d}}{C_0^{\mu_e}} = \frac{F_0^\mu}{G_0^\mu} \quad \text{or} \quad F_0^\mu = G_0^\mu C_0^\mu, \quad (44)$$

so from (37) and (38) we readily obtain

$$C_0^\mu f(\mu\sqrt{\alpha_e / \alpha_d} R_e) = g(\mu R_e) \quad \text{and} \quad C_0^\mu f'(\mu\sqrt{\alpha_e / \alpha_d} R_e) = \frac{k_e}{k_d} \sqrt{\frac{\alpha_d}{\alpha_e}} g'(\mu R_e), \quad (45)$$

where the prime denotes differentiation with respect to the argument for each case. The latter are two transcendental conditions, being solved numerically for the evaluation of the ratio $C_0^\mu$ from (44) and the parameter $\mu$. Once done, parameter $\mu_d$ is readily recovered via (41), while we remain with only one coefficient to calculate, either $F_0^\mu$ or $G_0^\mu$, since the ratio (44) is a known quantity. Incorporating the above relationships (37)–(45) to the fields (35) and (36), we may directly rewrite them as

$$\tilde{T}_d(\rho, t) = \int \sum_{\mu} G_0^{\mu} C_0^{\mu} f(\mu\sqrt{\alpha_e / \alpha_d} \rho) \exp(-\alpha_e \mu^2 t) \text{ for } \rho \in (R_d, R_e) \text{ and } t > 0 \quad (46)$$

and

$$\tilde{T}_e(\rho, t) = \int \sum_{\mu} G_0^{\mu} g(\mu\rho) \exp(-\alpha_e \mu^2 t) \text{ for } \rho \in (R_e, R_s) \text{ and } t > 0. \quad (47)$$

Substituting (46) and (47) into (8) and (9) (taking into account the azimuthal independence), while attaching (21) and (22) with (40), we get

$$T_d(\rho, t) = T_\infty - \frac{qR_d}{k_d} \left( \ln \frac{\rho}{R_e} - \frac{k_d}{k_e} \ln \frac{R_s}{R_e} - \frac{k_d}{hR_s} \right) + \int \sum_{\mu} G_0^{\mu} C_0^{\mu} f(\mu\sqrt{\alpha_e / \alpha_d} \rho) \exp(-\alpha_e \mu^2 t) \quad (48)$$

for every $\rho \in [R_d, R_e]$ and $t > 0$, while
$$ T_e(\rho, t) = T_\infty - \frac{qR_d}{k_e} \left( \ln \frac{\rho}{R_s} - \frac{k_e}{hR_s} \right) + \int \sum_{\mu} G_0^{\mu} g(\mu\rho) \exp(-\alpha_e \mu^2 t) \quad (49) $$

for every $\rho \in [R_e, R_s]$ and $t > 0$.

Thus, the initial conditions (7) or equivalently (20), by virtue of (48) and (49), render the expressions

$$ \int \sum_{\mu} G_0^{\mu} C_0^{\mu} f(\mu \sqrt{\alpha_e / \alpha_d} \rho) = \frac{qR_d}{k_d} \ln \frac{\rho}{R_e} - qR_d \left( \frac{1}{k_e} \ln \frac{R_s}{R_e} + \frac{1}{hR_s} \right) \quad \text{for every } \rho \in [R_d, R_e] \quad (50) $$

and

$$ \int \sum_{\mu} G_0^{\mu} g(\mu\rho) = \frac{qR_d}{k_e} \ln \frac{\rho}{R_s} - \frac{qR_d}{hR_s} \quad \text{for every } \rho \in [R_e, R_s], \quad (51) $$

which comprise two relations that, beyond the fact that they contain the $\rho$-variable for $\rho \in [R_d, R_s]$, they also depend upon the physical and geometrical characteristics of the system under consideration. The last tricky step is to calculate the remaining coefficient $G_0^\mu$ from the dual series conditions (50) and (51), where despite their simplicity, special treatment is required according to the following. Since we attain two relationships with the same unknown, our first task includes the merging of (50) and (51) into a single expression, which means that they can be written as

$$ \int \sum_{\mu} G_0^{\mu} P(\mu\rho) = p(\rho) \quad \text{for every } \rho \in [R_d, R_e] \cup [R_e, R_s] \quad \text{or} \quad \rho \in [R_d, R_s], \quad (52) $$

where

$$ P(\mu\rho) = \begin{cases} C_0^\mu f(\mu\sqrt{\alpha_e / \alpha_d}\rho), & \rho \in [R_d, R_e] \\ g(\mu\rho), & \rho \in [R_e, R_s] \end{cases} \quad (53) $$

and

$$ p(\rho) = \begin{cases} \frac{qR_d}{k_d} \ln \frac{\rho}{R_e} - qR_d \left( \frac{1}{k_e} \ln \frac{R_s}{R_e} + \frac{1}{hR_s} \right), & \rho \in [R_d, R_e] \\ \frac{qR_d}{k_e} \ln \frac{\rho}{R_s} - \frac{qR_d}{hR_s}, & \rho \in [R_e, R_s] \end{cases} \quad (54) $$

that conveniently reduces our problem to that of solving one condition (52) with the single unknown coefficient $G_0^\mu$. Our second task is now to eliminate the $\rho$-variable from (52) and towards this direction we multiply relationship (52) by $J_0(\mu'\rho)d\rho$ and we integrate overall from $R_d$ to $R_s$ in an appropriate manner, based on integration either from $R_d$ to $R_e$ or from $R_e$ to $R_s$, according to (53) and (54). During this procedure, we observe that functions (42) and (43), which appear on the left-hand side of (52), due to definition (53), involve the zeroth-order Bessel functions $J_0$ and
Neumann functions $N_0$, which are defined both for the same and for different arguments. Hence, it is apparent that we need to analytically calculate the definite integrals

$$ 
\begin{aligned}
P_{\mu,\mu'} &\equiv \int_{R_d}^{R_e} J_0(\mu'\rho) J_0(\mu\sqrt{\alpha_e / \alpha_d}\rho) \rho d\rho \\
&= R_e \frac{\mu\sqrt{\alpha_e / \alpha_d} J_0(\mu' R_e) J_1(\mu\sqrt{\alpha_e / \alpha_d} R_e) - \mu' J_1(\mu' R_e) J_0(\mu\sqrt{\alpha_e / \alpha_d} R_e)}{(\alpha_e / \alpha_d)\mu^2 - \mu'^2} \\
&\quad - R_d \frac{\mu\sqrt{\alpha_e / \alpha_d} J_0(\mu' R_d) J_1(\mu\sqrt{\alpha_e / \alpha_d} R_d) - \mu' J_1(\mu' R_d) J_0(\mu\sqrt{\alpha_e / \alpha_d} R_d)}{(\alpha_e / \alpha_d)\mu^2 - \mu'^2},
\end{aligned} 
\quad (55) $$

$$ 
\begin{aligned}
Q_{\mu,\mu'} &\equiv \int_{R_d}^{R_e} J_0(\mu'\rho) N_0(\mu\sqrt{\alpha_e / \alpha_d}\rho) \rho d\rho \\
&= R_e \frac{\mu\sqrt{\alpha_e / \alpha_d} J_0(\mu' R_e) N_1(\mu\sqrt{\alpha_e / \alpha_d} R_e) - \mu' J_1(\mu' R_e) N_0(\mu\sqrt{\alpha_e / \alpha_d} R_e)}{(\alpha_e / \alpha_d)\mu^2 - \mu'^2} \\
&\quad - R_d \frac{\mu\sqrt{\alpha_e / \alpha_d} J_0(\mu' R_d) N_1(\mu\sqrt{\alpha_e / \alpha_d} R_d) - \mu' J_1(\mu' R_d) N_0(\mu\sqrt{\alpha_e / \alpha_d} R_d)}{(\alpha_e / \alpha_d)\mu^2 - \mu'^2},
\end{aligned} 
\quad (56) $$

for every value of $\mu, \mu'$ and

$$ 
\begin{aligned}
U_{\mu,\mu'} &\equiv \int_{R_e}^{R_s} J_0(\mu'\rho) J_0(\mu\rho) \rho d\rho = R_s \frac{\mu J_0(\mu' R_s) J_1(\mu R_s) - \mu' J_1(\mu' R_s) J_0(\mu R_s)}{\mu^2 - \mu'^2} \\
&\quad - R_e \frac{\mu J_0(\mu' R_e) J_1(\mu R_e) - \mu' J_1(\mu' R_e) J_0(\mu R_e)}{\mu^2 - \mu'^2},
\end{aligned} 
\quad (57) $$

$$ 
\begin{aligned}
V_{\mu,\mu'} &\equiv \int_{R_e}^{R_s} J_0(\mu'\rho) N_0(\mu\rho) \rho d\rho = R_s \frac{\mu J_0(\mu' R_s) N_1(\mu R_s) - \mu' J_1(\mu' R_s) N_0(\mu R_s)}{\mu^2 - \mu'^2} \\
&\quad - R_e \frac{\mu J_0(\mu' R_e) N_1(\mu R_e) - \mu' J_1(\mu' R_e) N_0(\mu R_e)}{\mu^2 - \mu'^2},
\end{aligned} 
\quad (58) $$

for $\mu \neq \mu'$, while

$$ U_{\mu,\mu} \equiv \int_{R_e}^{R_s} J_0^2(\mu\rho) \rho d\rho = \frac{R_s^2}{2} [J_0^2(\mu R_s) + J_1^2(\mu R_s)] - \frac{R_e^2}{2} [J_0^2(\mu R_e) + J_1^2(\mu R_e)], \quad (59) $$

$$ 
\begin{aligned}
V_{\mu,\mu} &\equiv \int_{R_e}^{R_s} J_0(\mu\rho) N_0(\mu\rho) \rho d\rho = \frac{R_s^2}{2} [J_0(\mu R_s) N_0(\mu R_s) + J_1(\mu R_s) N_1(\mu R_s)] \\
&\quad - \frac{R_e^2}{2} [J_0(\mu R_e) N_0(\mu R_e) + J_1(\mu R_e) N_1(\mu R_e)]
\end{aligned} 
\quad (60) $$

when $\mu = \mu'$, being given in a convenient fashion. Otherwise, the right-hand side of (52) includes
simple logarithmic and constant functions, due to definition (54), consequently, we need to evaluate
the trivial formulae

$$
\begin{align}
X_{\mu'} &\equiv \int_{R_d}^{R} \ln \frac{\rho}{R_e} J_0(\mu' \rho) \rho d\rho = \int_{R_d}^{R} \ln \frac{\rho}{R_e} \frac{d}{d\rho} \left[ \frac{\rho J_1(\mu' \rho)}{\mu'} \right] d\rho = \frac{1}{\mu'} \left[ \rho J_1(\mu' \rho) \ln \frac{\rho}{R_e} \right]_{R_d}^{R_e} - \frac{1}{\mu'} \int_{R_d}^{R} J_1(\mu' \rho) d\rho \nonumber \\
&= \frac{1}{\mu'} \left[ R_e J_1(\mu' R_e) \ln \frac{R_e}{R_e} - R_d J_1(\mu' R_d) \ln \frac{R_d}{R_e} \right] - \frac{1}{\mu'^2} \int_{\mu' R_d}^{\mu' R_e} J_1(\mu' \rho) d(\mu' \rho) \nonumber \\
&= \frac{1}{\mu'^2} \left\{ [J_0(\mu' R_e) - J_0(\mu' R_d)] + \mu' R_d J_1(\mu' R_d) \ln \frac{R_e}{R_d} \right\} \tag{61}
\end{align}
$$

and, similarly,

$$
Y_{\mu'} \equiv \int_{R_e}^{R} \ln \frac{\rho}{R_s} J_0(\mu' \rho) \rho d\rho = \frac{1}{\mu'^2} \left\{ [J_0(\mu' R_s) - J_0(\mu' R_e)] + \mu' R_e J_1(\mu' R_e) \ln \frac{R_s}{R_e} \right\}, \quad (62)
$$

as well as

$$
Z_{\mu'}(\rho) \equiv \int J_0(\mu'\rho)\rho d\rho = \frac{1}{\mu'^2} \int (\mu'\rho)J_0(\mu'\rho)d(\mu'\rho) = \frac{1}{\mu'^2}(\mu'\rho)J_1(\mu'\rho) = \frac{\rho}{\mu'}J_1(\mu'\rho), \quad (63)
$$

where the latter is defined for $\rho \in [R_d, R_s]$. In the sequel, the outcomes (55)–(63), whose evaluation
is mainly based on the use of standard recurrence relations of Bessel functions [23], are readily
implemented into (52) (after being multiplied by $J_0(\mu'\rho)\rho d\rho$ and integrated with respect to the
interval $[R_d, R_s]$), taking into account (53) and (54) as well, hence, we obtain straightforwardly

$$
\int \sum_{\mu} (u_{0}^{\mu, \mu'} + v_{0}^{\mu, \mu'}) G_{0}^{\mu} = x_{0}^{\mu'} + y_{0}^{\mu'} \quad \text{for every } \mu' \in \mathbb{R}, \qquad (64)
$$

wherein

$$
u_0^{\mu, \mu'} \equiv C_0^{\mu} \left[ \frac{P_{\mu, \mu'}}{J_1(\mu \sqrt{\alpha_e / \alpha_d} R_d)} - \frac{Q_{\mu, \mu'}}{N_1(\mu \sqrt{\alpha_e / \alpha_d} R_d)} \right], \quad (65)
$$

$$
x_0^{\mu'} \equiv \frac{q R_d}{k_d} X_{\mu'} - q R_d \left( \frac{1}{k_e} \ln \frac{R_s}{R_e} + \frac{1}{h R_s} \right) [Z_{\mu'}(R_e) - Z_{\mu'}(R_d)] \quad (66)
$$

and

$$
v_0^{\mu, \mu'} \equiv \frac{U_{\mu, \mu'}}{[hJ_0(\mu R_s) - k_e \mu J_1(\mu R_s)]} - \frac{V_{\mu, \mu'}}{[hN_0(\mu R_s) - k_e \mu N_1(\mu R_s)]}, \quad (67)
$$

$$
y_0^{\mu'} \equiv \frac{q R_d}{k_e} Y_{\mu'} - \frac{q R_d}{h R_s} [Z_{\mu'}(R_s) - Z_{\mu'}(R_e)] \quad (68)
$$

with $\mu, \mu' \in \mathbb{R}$. Here, we observe that in order to manipulate condition (64), we are obliged to set the
parameter $\mu$ so as to obtain discrete values, that is $\mu \equiv \mu_{0,m}$ for every $m \ge 1$ (similarly $\mu' \equiv \mu_{0,m'}$ with $m' \ge 1$), whereas the zeroth order of the Bessel functions is readily implied to the notation. Then, the aforementioned series-integral symbolism becomes now a series, due to the eigenvalue property of parameter $\mu$, thus $\int \sum_{\mu} \cdots \to \sum_{m=1}^{+\infty} \cdots$, and therein our analysis is restricted to a classical infinite series semi-analytical solution.

Under this aim, we recapitulate for convenience our findings in the present analysis, hence the temperature distributions within the dielectric tube and the electrode that come from expressions (48) and (49), admit

$$T_d(\rho, t) = T_\infty - \frac{qR_d}{k_d} \left( \ln \frac{\rho}{R_e} - \frac{k_d}{k_e} \ln \frac{R_s}{R_e} - \frac{k_d}{hR_s} \right) + \sum_{m=1}^{+\infty} G_0^m C_0^m f(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} \rho) \exp(-\alpha_e \mu_{0,m}^2 t) \quad (69)$$

for every $\rho \in [R_d, R_e]$ and $t \ge 0$, while

$$T_e(\rho, t) = T_\infty - \frac{qR_d}{k_e} \left( \ln \frac{\rho}{R_s} - \frac{k_e}{hR_s} \right) + \sum_{m=1}^{+\infty} G_0^m g(\mu_{0,m}\rho) \exp(-\alpha_e \mu_{0,m}^2 t) \quad (70)$$

for every $\rho \in [R_d, R_s]$ and $t \ge 0$, where formulae (42) and (43) become

$$f(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} \rho) = \frac{J_0(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} \rho)}{J_1(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_d)} - \frac{N_0(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} \rho)}{N_1(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_d)} \quad (71)$$

and

$$g(\mu_{0,m}\rho) = \frac{J_0(\mu_{0,m}\rho)}{[hJ_0(\mu_{0,m}R_s) - k_e\mu_{0,m}J_1(\mu_{0,m}R_s)]} - \frac{N_0(\mu_{0,m}\rho)}{[hN_0(\mu_{0,m}R_s) - k_e\mu_{0,m}N_1(\mu_{0,m}R_s)]} \quad (72)$$

for $m \ge 1$, while $C_0^m$ and $\mu_{0,m}$ are the solutions of the two transcendental relationships (45), which provide us with

$$C_0^m = \frac{g(\mu_{0,m} R_e)}{f(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_e)} \quad \text{and} \quad \frac{f'(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_e)}{f(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_e)} = \frac{k_e}{k_d} \sqrt{\frac{\alpha_d}{\alpha_e}} \frac{g'(\mu_{0,m} R_e)}{g(\mu_{0,m} R_e)} \quad \text{for } m \ge 1, \quad (73)$$

wherein $C_0^m$ for $m \ge 1$ is obtained analytically from the left-hand side of (73), while the evaluation of $\mu_{0,m}$ for $m \ge 1$ is subject to further numerical handling (see in the sequel for a more detailed and accurate explanation on that matter), as it is revealed from the right-hand side of the set of relationships (73). Otherwise, $G_0^m$ for $m \ge 1$ can be evaluated as the outcome of the infinite system of linear algebraic equations

$$\sum_{m=1}^{+\infty} (u_0^{m,m'} + v_0^{m,m'}) G_0^m = x_0^{m'} + y_0^{m'} \quad \text{for every } m' \ge 1, \quad (74)$$
resulting from relation (64) with formulae (65)–(68) and (55)–(63), in which the involved
mathematical quantities assume the following complicated, yet handy, expressions in terms of the
geometrical and physical parameters of the system, i.e.,

$$
\begin{align}
x_0^{m'} &= \frac{qR_d}{k_d} X_{m'} - qR_d \left( \frac{1}{k_e} \ln \frac{R_s}{R_e} + \frac{1}{hR_s} \right) [Z_{m'}(R_e) - Z_{m'}(R_d)] \nonumber \\
&= \frac{qR_d}{k_d} \frac{1}{\mu_{0,m'}^2} \left\{ \left[ J_0(\mu_{0,m'} R_e) - J_0(\mu_{0,m'} R_d) \right] + \mu_{0,m'} R_d J_1(\mu_{0,m'} R_d) \ln \frac{R_e}{R_d} \right\} \nonumber \\
&\quad - \frac{qR_d}{\mu_{0,m'}} \left( \frac{1}{k_e} \ln \frac{R_s}{R_e} + \frac{1}{hR_s} \right) [R_e J_1(\mu_{0,m'} R_e) - R_d J_1(\mu_{0,m'} R_d)] \quad \text{for } m' \ge 1 \tag{75}
\end{align}
$$

and

$$
\begin{align}
y_0^{m'} &= \frac{qR_d}{k_e} Y_{m'} - \frac{qR_d}{hR_s} [Z_{m'}(R_s) - Z_{m'}(R_e)] \nonumber \\
&= \frac{qR_d}{k_e} \frac{1}{\mu_{0,m'}^2} \left\{ \left[ J_0(\mu_{0,m'} R_s) - J_0(\mu_{0,m'} R_e) \right] + \mu_{0,m'} R_e J_1(\mu_{0,m'} R_e) \ln \frac{R_s}{R_e} \right\} \nonumber \\
&\quad - \frac{qR_d}{\mu_{0,m'} hR_s} [R_s J_1(\mu_{0,m'} R_s) - R_e J_1(\mu_{0,m'} R_e)] \quad \text{for } m' \ge 1, \tag{76}
\end{align}
$$

while

$$
u_0^{m,m'} = C_0^m \left[ \frac{P_{m,m'}}{J_1(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_d)} - \frac{Q_{m,m'}}{N_1(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_d)} \right], \quad (77)
$$

$$
v_0^{m,m'} = \frac{U_{m,m'}}{[hJ_0(\mu_{0,m}R_s) - k_e\mu_{0,m}J_1(\mu_{0,m}R_s)]} - \frac{V_{m,m'}}{[hN_0(\mu_{0,m}R_s) - k_e\mu_{0,m}N_1(\mu_{0,m}R_s)]} \quad (78)
$$

with

$$
\begin{equation}
\begin{split}
P_{m,m'} = R_e \frac{\mu_{0,m} \sqrt{\alpha_e / \alpha_d} J_0(\mu_{0,m'} R_e) J_1(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_e) - \mu_{0,m'} J_1(\mu_{0,m'} R_e) J_0(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_e)}{(\alpha_e / \alpha_d) \mu_{0,m}^2 - \mu_{0,m'}^2} \\
- R_d \frac{\mu_{0,m} \sqrt{\alpha_e / \alpha_d} J_0(\mu_{0,m'} R_d) J_1(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_d) - \mu_{0,m'} J_1(\mu_{0,m'} R_d) J_0(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_d)}{(\alpha_e / \alpha_d) \mu_{0,m}^2 - \mu_{0,m'}^2},
\end{split}
\tag{79}
\end{equation}
$$

$$
\begin{equation}
\begin{split}
Q_{m,m'} = R_e \frac{\mu_{0,m} \sqrt{\alpha_e / \alpha_d} J_0(\mu_{0,m'} R_e) N_1(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_e) - \mu_{0,m'} J_1(\mu_{0,m'} R_e) N_0(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_e)}{(\alpha_e / \alpha_d) \mu_{0,m}^2 - \mu_{0,m'}^2} \\
- R_d \frac{\mu_{0,m} \sqrt{\alpha_e / \alpha_d} J_0(\mu_{0,m'} R_d) N_1(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_d) - \mu_{0,m'} J_1(\mu_{0,m'} R_d) N_0(\mu_{0,m} \sqrt{\alpha_e / \alpha_d} R_d)}{(\alpha_e / \alpha_d) \mu_{0,m}^2 - \mu_{0,m'}^2}
\end{split}
\tag{80}
\end{equation}
$$

and
$$ 
\begin{aligned}
U_{m,m'} &= R_s \frac{\mu_{0,m} J_0(\mu_{0,m'} R_s) J_1(\mu_{0,m} R_s) - \mu_{0,m'} J_1(\mu_{0,m'} R_s) J_0(\mu_{0,m} R_s)}{\mu_{0,m}^2 - \mu_{0,m'}^2} \\
&\quad - R_e \frac{\mu_{0,m} J_0(\mu_{0,m'} R_e) J_1(\mu_{0,m} R_e) - \mu_{0,m'} J_1(\mu_{0,m'} R_e) J_0(\mu_{0,m} R_e)}{\mu_{0,m}^2 - \mu_{0,m'}^2},
\end{aligned}
\quad (81) $$

$$ 
\begin{aligned}
V_{m,m'} &= R_s \frac{\mu_{0,m} J_0(\mu_{0,m'} R_s) N_1(\mu_{0,m} R_s) - \mu_{0,m'} J_1(\mu_{0,m'} R_s) N_0(\mu_{0,m} R_s)}{\mu_{0,m}^2 - \mu_{0,m'}^2} \\
&\quad - R_e \frac{\mu_{0,m} J_0(\mu_{0,m'} R_e) N_1(\mu_{0,m} R_e) - \mu_{0,m'} J_1(\mu_{0,m'} R_e) N_0(\mu_{0,m} R_e)}{\mu_{0,m}^2 - \mu_{0,m'}^2},
\end{aligned}
\quad (82) $$

for $m \neq m'$, whilst

$$ U_{m,m} = \frac{R_s^2}{2} [J_0^2(\mu_{0,m} R_s) + J_1^2(\mu_{0,m} R_s)] - \frac{R_e^2}{2} [J_0^2(\mu_{0,m} R_e) + J_1^2(\mu_{0,m} R_e)], \quad (83) $$

$$ 
\begin{aligned}
V_{m,m} &= \frac{R_s^2}{2} [J_0(\mu_{0,m} R_s) N_0(\mu_{0,m} R_s) + J_1(\mu_{0,m} R_s) N_1(\mu_{0,m} R_s)] \\
&\quad - \frac{R_e^2}{2} [J_0(\mu_{0,m} R_e) N_0(\mu_{0,m} R_e) + J_1(\mu_{0,m} R_e) N_1(\mu_{0,m} R_e)]
\end{aligned}
\quad (84) $$

for $m=m'$, all above provided for every $m, m' \ge 1$.

This sequence of systems of equations (74) with (75)–(84) appears often in problems of mathematical physics and our effort is then limited in developing a solution, derived from a partial differential equation problem with mixed boundaries. To this end, we introduce an idea, based on usual cut-off techniques in order to solve infinite linear systems. According to this methodology, we truncate the series (74) by the imposition of an indispensable upper limit $M \in \mathbb{N}^*$, so as the series become finite, due to the fact that we may write $\sum_{m=1}^{+\infty} \dots \to \sum_{m=1}^{M} \dots$, and facilitate the procedure.

Thereafter, (74) reduces to a finite algebraic linear system of equations if we force the number of the incorporated equations to be finite as well, that is $m'=1,2,...,M'$, wherein $M' \in \mathbb{N}^*$. This system becomes quadrature only if $M=M'=L \in \mathbb{N}^*$ and once done, it is rewritten in matrix form as

$$ A\mathbf{x} = \mathbf{b}, \text{ where } A = \begin{bmatrix} \cdots & \vdots & \cdots \\ \cdots & u_0^{m,m'} + v_0^{m,m'} & \cdots \\ \vdots & \vdots & \vdots \end{bmatrix}, \quad \mathbf{x} = \begin{bmatrix} \vdots \\ G_0^m \\ \vdots \end{bmatrix}, \quad \text{and} \quad \mathbf{b} = \begin{bmatrix} \vdots \\ x_0^{m'} + y_0^{m'} \\ \vdots \end{bmatrix}, \quad (85) $$

in which $m,m'=1,2,...,L$, while $A$ is the squared-type invertible matrix of the coefficients of the unknowns, $\mathbf{x}$ is the vector of the unknown coefficients and the $\mathbf{b}$ is vector of the known constants. Thus, each time we repetitively solve quadrature systems numerically for increasing $L$ and the standard cut-off is applied in order to decide the truncation amplitude $L_{\max}$, for which the number of the evaluated constant coefficients are sufficient to obtain the expected accuracy in order for the series to converge. The solution of the system (85) yields the sought coefficients $\mathbf{x}$, i.e., the unknown constant coefficients $G_0^m$ for every $m=1,2,...,L_{\max}$.
The final solution of the problem is eventually achieved via the temperature functions (69) and (70), followed by the conveniently chosen functions (71) and (72), in which the constant coefficients $C_0^m$ and $\mu_{0,m}$ for $m=1,2,...,L_{max}$ are the solutions of relationships (73), while the set of constants $G_0^m$ for $m=1,2,...,L_{max}$ that satisfies the easily amenable infinite set of linear relationships (74) with (75)–(84), are given by (85), i.e., $\mathbf{x}=\mathbf{A}^{-1}\mathbf{b}$, when the incorporated inverse matrix $\mathbf{A}^{-1}$ is given. Under the aim of presenting a compact expression for the temperature field, we utilize the definition of the Heaviside step function $H(x)=1$ for $x \ge 0$, otherwise zero, and we incorporate fields (69) and (70) into

$$T(\rho,t) = H(R_e - \rho)T_d(\rho,t) + [H(R_s - \rho) - H(R_e - \rho)]T_e(\rho,t) \quad \text{for } \rho \in [R_d, R_s] \text{ and } t \ge 0, \quad (86)$$

which provides the temperature in both the dielectric tube and the external electrode in a handy analytical fashion. Therefore, the value of the temperature on the outer cover of the electrode for $\rho = R_s$, being of our main concern, is provided theoretically according to (86) and (70) (or from (70) directly) as

$$T_s(t) = T_{\infty} + \frac{qR_d}{hR_s} + \sum_{m=1}^{+\infty} G_0^m g(\mu_{0,m} R_s) \exp(-\alpha_e \mu_{0,m}^2 t), \quad (87)$$

where $T_s(t) \equiv T_e(R_s,t)$ for any $t \ge 0$, while in practice the upper infinite limit of the series (87) must be truncated to $L_{max}$ for engineering applications.

In order to demonstrate the above analysis, we devote the rest of the analysis to the numerical implementation and validation of the obtained analytical expressions for the temperature distribution in the domains of interest, i.e., within the dielectric tube and the external electrode. Towards this direction, we provide specific indicative values of the involved physical parameters and properties, appearing during our analysis [22]. To this end, we consider a representative value for the heat conduction $q=18624 \text{ Wm}^{-2}$ in plasma jets and $h=13.2 \text{ Wm}^{-2}\text{K}^{-1}$ for such systems, while we have $k_d=30 \text{ Wm}^{-1}\text{K}^{-1}$ and $k_e=110 \text{ Wm}^{-1}\text{K}^{-1}$, wherein we directl calculate $a_d=k_d/\rho_d c_{p,d}=9.9 \times 10^{-6} \text{ m}^2\text{s}^{-1}$ and $a_e=k_e/\rho_e c_{p,e}=3.3 \times 10^{-5} \text{ m}^2\text{s}^{-1}$, since it is $\rho_d=3800 \text{ kgm}^{-3}$, $\rho_e=8600 \text{ kgm}^{-3}$ and $c_{p,d}=800 \text{ Jkg}^{-1}\text{K}^{-1}$, $c_{p,e}=388 \text{ Jkg}^{-1}\text{K}^{-1}$, respectively for each case and the geometrical characteristics of the system yield $R_d=5.7 \times 10^{-4} \text{ m}$, $R_e=1.25 \times 10^{-3} \text{ m}$ and $R_s=5 \times 10^{-3} \text{ m}$. Note that these values assume units in K for the temperature, however our plots use the corresponding units in °C. Given these representative values, we initially have to numerically solve the right-hand-side of (73) so as to find the proper values of $\mu_{0,m}$ for $m \ge 1$ that force the new-defined function
$$h(x) = \frac{f'(\sqrt{\alpha_e / \alpha_d} x)}{f(\sqrt{\alpha_e / \alpha_d} x)} - \frac{k_e}{k_d} \sqrt{\frac{\alpha_d}{\alpha_e}} \frac{g'(x)}{g(x)} \quad \text{with } x = \mu_{0,m} R_e \quad \text{for } m \ge 1 \qquad (88)$$

to become zero, i.e., $h(x)=0$. To ensure that all the roots were computed, the sign of the equation was checked at every interval between two consecutive roots. It was found that the roots do not follow a specific pattern, so it was necessary to use very small steps in the control procedure. These roots are given in Figure 2, hence we have solved the equivalent relationship of the second equation in (73) for the evaluation of $\mu_{0,m}$, $m \ge 1$.

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Visualization of the transcendental function $h(x)$ with $x = \mu_{0,m} R_e$ for $m \ge 1$ behavior for the lowest values of the wavenumber. The function roots are marked with blue circles.</figcaption></figure>

Having found the roots $\mu_{0,m}$ for $m \ge 1$, it is necessary to calculate the terms (75)–(84), which are necessary for computing the linear system of algebraic equations (85). For different values of $m$ and $m'$ the procedure was straightforward in (75)–(82), while for the diagonal elements, where $m = m'$, it is obvious that we have to refer to relations (83) and (84). The resulting coefficients for the system (85) are obtained via the application of the cut-off method. Thereafter, the temperature distribution $T$ in the dielectric tube and the external electrode follows from relationship (86), wherein its variation with respect to time $t$ at the middle thickness of the dielectric tube $(R_d + R_e)/2$ and the external electrode $(R_e + R_s)/2$, as well as with respect to the radius $\rho$ of the system at 3000s, which is approximately steady state, is implemented and depicted in Figure 3 and
Figure 4, respectively. Obviously, the derivation of the temperature refers to the transient heat transfer problem, as shown in Figure 3, even though the verification in Figure 4 was for steady temperature. The results of the finalized semi-analytical solution have been compared against numerical simulations carried out using the Finite Integration Technique (FIT). FIT is a grid-based method very similar to the Finite Difference in Time Domain (FDTD), originally introduced for the solution of the Maxwell equations, using a system of staggered orthogonal grids [27] and which has been extended to acoustics, elastodynamics (the extension is known as Elastodynamics Finite Integration Technique (EFIT)) [28], as well as to heat equation problems [29]. We indicate that for numerical convenience purposes the temperature $T$ has been re-gauged so as the initial temperature of the system is 0°C instead of room temperature.

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Temperature transients $T$ as function of time $t$ at the middle thickness of the dielectric tube (left) and the external electrode (right), using the analytical solution (solid line) and the Finite Integration Technique (FIT) (dash line).</figcaption></figure>

The analytical solution computations are performed in two steps. The first (and most tedious one) concerns the calculation of the eigenvalues. For the (non-optimized) root finding algorithm used in the present implementation, this time was 17s. The second step comprises the inversion of the system matrix for the estimation of the development coefficients and the evaluation of the series. The computation time for this second part has been counted to 80ms. These computational times have to be compared with the respective ones required for the numerical solution of the problem. Since we are dealing with a diffusion problem in the time domain, the considered time window has to be sampled and the discrete linear system of equation obtained by the application of the method's grid (either using FIT or FEM) must be inverted at each time instance (implicit formulation). For the specific example examined here, the total computational time reached is 4s. Clearly, there is a computational overhead when we need to evaluate the eigenvalues of the analytical solution. This calculation however has to be carried out only once, independently from the excitation signal. Every other evaluation is of the order of ms, whereas the numerical solution has to be re-run for each new
excitation. The gain in performance increases dramatically when moving at higher dimensions, the cost, however, that one has to pay is that the analytical solution becomes cumbersome.

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Temperature distribution $T$ as function of the radius $\rho$ at the steady state, using the analytical solution (solid line) and the Finite Integration Technique (FIT) (dash line).</figcaption></figure>

At this stage, we provide the necessary discussion about the relative errors during our methodology. There are three sources of error that enter in the computation of the analytical solution, i.e., the error of the eigenvalue calculation, the numerical evaluation of the special functions and the truncation of the development series. The first one is controlled by the numerical solution of the transcendental equation (88), which is determined by the method threshold and the round-off error. The evaluation error of the special (Bessel) functions is determined by the specific numerical implementation of the computational platform (in our case Matlab) to which one usually has little control, unless dedicated implementations are used. With higher orders, overflows occur; therefore, in our implementation, the Bessel functions are always evaluated in form of ratios. The function evaluation error can be assumed having negligible impact to the results. The third one (number of series term taken into account) is the one with the strongest impact to the solution precision. There is a trade-off between precision and efficiency and one has to resort to numerical experimentation in order to pick-up the optimal number. In addition, care must be taken with higher-order terms, since the root finding algorithm becomes less accurate due to overflows (despite the normalization of the Bessel functions). Taking all these factors into consideration and after numerical experimentation, we have concluded that a number of 36 terms gives the most satisfactory results in terms of the computational time and precision. Concerning the issue of error variation as function of the position raised, one can expect that it may be higher at the interfaces, where the continuity is imposed by matching the series terms.

The analytical solution is compared with results provided by the FIT method for a given test
configuration, whose parameters have been carefully chosen in order to represent a realistic case. Since the two calculations (analytical and numerical one) are independent and given that the FIT computations have been carried out using a generic code, validated for a number of different test problems (in fact the FIT solver is part of the CIVA commercial simulation platform), we accept that the agreement of the two data sets validates the correctness of the presented solution. It can be seen that the calculated initial condition displays a great variation, probably due to the over defined nature of the algebraic system. According to the problem formulation and expression (87), the temperature on the surface of the electrode as $t \to +\infty$ (or large time scaling for applications) converges to the temperature 160.84 °C, which is the plateau in Figure 3. On the other hand, for $t=0$ (initial condition), the room temperature is readily recovered from (87), as it is shown in Figure 3.

The analytical technique presented herein, deals with the 1D three-layer problem of the plasma cell, described in Figure 1. It can be easily extended to address finite multilayer configurations with a consequent increase, however, of the solution complexity. A basic limitation of the method is that it remains applicable to piecewise homogeneous geometries, i.e., material gradients cannot be tackled without drastic modifications, such as the introduction of hybrid development bases mixing the herein applied Fourier basis with a spatial mesh.

### 3. Conclusions and discussion

The present work was devoted to the presentation of an analytical technique for the determination of the temperature in a representative plasma jet reactor. We studied the heat transfer process from the plasma area to the outer environment through the solid components of the reactor, i.e., the dielectric tube and the grounded electrode.

The physical problem itself was mathematically formulated with respect to the cylindrical geometry, which is employed, since it fits the plasma jet reactor. We separated the area of thermal activity into two distinct domains, matching the one on the dielectric tube and the other of the external electrode, considering them as coaxial cylinders. Under this aim, we constructed appropriately the corresponding initial and boundary value problems in each domain with either Dirichlet or Neumann conditions and the solution technique was based on the method of separation of variables in the cylindrical coordinate system and the imposition of a particular method of asymptotic kernels for the manipulation of the non-homogeneous boundary conditions. The time-dependent temperature fields within the tube and the electrode were constructed in terms of cylindrical harmonic eigenfunctions and the analytical solution was obtained in view of the restriction of rotational symmetry, due to the nature of the boundary conditions. During the analysis, all the implicated unknown constant coefficients were evaluated straightforwardly, except one. For this last unknown, we were obliged to solve a pair of dual Fourier-Bessel series relations, which were handled simultaneously as a single relationship and led to infinite linear systems for the unknowns via integration techniques, involving mostly Bessel and Neumann eigenfunctions among other simple functions. These systems were solved with the aim of a repetitive cut-off method until the expected accuracy was achieved.

The efficiency of the analysis was demonstrated by implementing numerically the expressions of the final formula, which satisfy the initial and boundary value problem under consideration, while a Finite Integration Technique (FIT) was applied for the validation of the methodology. It was shown
that the distribution of the temperature behaves as expected, as far as the spatial dependence is concerned, while as the time evolves, it converges.

## Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

## Conflict of interest

The authors declare there is no conflict of interest.

## References

1. V. Papadimas, C. Doudesis, P. Svarnas, P. K. Papadopoulos, G. P. Vafakos, P. Vafeas, SDBD flexible plasma actuator with Ag-Ink electrodes: experimental assessment, *Appl. Sci.*, **11** (2021), 11930. https://doi.org/10.3390/app112411930

2. P. Svarnas, E. Giannakopoulos, I. Kalavrouziotis, C. Krontiras, S. Georga, R. S. Pasolari, et al., Sanitary effect of FE-DBD cold plasma in ambient air on sewage biosolids, *Sci. Total Environ.*, **705** (2020), 135940. https://doi.org/10.1016/j.scitotenv.2019.135940

3. P. Svarnas, A. Spiliopoulou, P. G. Koutsoukos, K. Gazeli, E. D. Anastassiou, Acinetobacter baumannii deactivation by means of DBD-Based helium plasma jet, *Plasma*, **2** (2019), 77-90. https://doi.org/10.3390/plasma2020008

4. K. Pefani-Antimisiari, D. K. Athanasopoulos, A. Marazioti, K. Sklias, M. Rodi, A. L. de Lastic, et al., Synergistic effect of cold atmospheric pressure plasma and free or liposomal doxorubicin on melanoma cells, *Sci. Rep.*, **11** (2021), 14788. https://doi.org/10.1038/s41598-021-94130-7

5. K. Gazeli, P. Svarnas, P. Vafeas, P. K. Papadopoulos, A. Gkelios, F. Clément, Investigation on streamers propagating into a helium jet in air at atmospheric pressure: Electrical and optical emission analysis, *J. Appl. Phys.*, **114** (2013), 103304. https://doi.org/10.1063/1.4820570

6. P. Vafeas, P. K. Papadopoulos, G. P. Vafakos, P. Svarnas, M. Doschoris, Modelling the electric field in reactors yielding cold atmospheric-pressure plasma jets, *Sci. Rep.*, **10** (2020), 5694. https://doi.org/10.1038/s41598-020-61939-7

7. P. K. Papadopoulos, P. Vafeas, P. Svarnas, K. Gazeli, P. M. Hatzikonstantinou, A. Gkelios, et al., Interpretation of the gas flow field modification induced by guided streamer ('plasma bullet') propagation, *J. Phys. D: Appl. Phys.*, **47** (2014), 425203. https://doi.org/10.1088/0022-3727/47/42/425203

8. P. Svarnas, P. K. Papadopoulos, P. Vafeas, A. Gkelios, F. Clément, A. Mavon, Influence of atmospheric pressure guided streamers (plasma bullets) on the working gas pattern in air, *IEEE Trans. Plasma Sci.*, **42** (2014), 2430-2431. https://doi.org/10.1109/TPS.2014.2322098

9. D. K. Logothetis, P. K. Papadopoulos, P. Svarnas, P. Vafeas, Numerical simulation of the interaction between helium jet flow and an atmospheric-pressure "plasma jet", *Comput. Fluids*, **140** (2016), 11-18. https://doi.org/10.1016/j.compfluid.2016.09.006
10. P. K. Papadopoulos, D. K. Athanasopoulos, K. Sklias, P. Svarnas, N. Mourousias, K. Vratsinis, et al., Generic residual charge based model for the interpretation of the electro-hydrodynamic effects in cold atmospheric pressure plasmas, *Plasma Sources Sci. Technol.*, **28** (2019), 065005.
https://doi.org/10.1088/1361-6595/ab0a3c

11. P. Svarnas, P. K. Papadopoulos, D. Athanasopoulos, K. Sklias, K. Gazeli, P. Vafeas, Parametric study of thermal effects in a capillary dielectric-barrier discharge related to plasma jet production: Experiments and numerical modelling, *J. Appl. Phys.*, **124** (2018), 064902.
https://doi.org/10.1063/1.5037141

12. T. Nozaki, Y. Miyazaki, Y. Unno, K. Okazaki, Energy distribution and heat transfer mechanisms in atmospheric pressure non-equilibrium plasmas, *J. Phys. D: Appl. Phys.*, **34** (2001), 3383-3390.
https://doi.org/10.1088/0022-3727/34/23/310

13. S. Y. Moon, W. A. Choe, Comparative study of rotational temperatures using diatomic OH, O2 and N2+ molecular spectra emitted from atmospheric plasmas, *Spectrochim. Acta, Part B*, **58** (2003), 249-257.
https://doi.org/10.1016/S0584-8547(02)00259-8

14. J. H. Kim, Y. H. Kim, Y. H. Choi, W. Choe, J. J. Choi, Y. S. Hwang, Optical measurements of gas temperatures in atmospheric pressure RF cold plasmas, *Surf. Coat. Technol.*, **171** (2003), 211-215.
https://doi.org/10.1016/S0257-8972(03)00273-1

15. C. Yubero, M. S. Dimitrijević, M. C. García, M. D. Calzada, Using the van der Waals broadening of the spectral atomic lines to measure the gas temperature of an argon microwave plasma at atmospheric pressure, *Spectrochim. Acta, Part B*, **62** (2007), 169-176.
https://doi.org/10.1016/j.sab.2007.02.008

16. A. Ionascut-Nedelcescu, C. Carlone, U. Kogelschatz, D. V. Gravelle, M. I. Boulos, Calculation of the gas temperature in a throughflow atmospheric pressure dielectric barrier discharge torch by spectral line shape analysis, *J. Appl. Phys.*, **103** (2008), 063305.
https://doi.org/10.1063/1.2891419

17. S. Hofmann, A. F. H. van Gessel, T. Verreycken, P. Bruggeman, Power dissipation, gas temperatures and electron densities of cold atmospheric pressure helium and argon RF plasma jets, *Plasma Sources Sci. Technol.*, **20** (2011), 065010.
https://doi.org/10.1088/0963-0252/20/6/065010

18. Z. S. Chang, G. J. Zhang, X. J. Shao, Z. H. Zhang, Diagnosis of gas temperature, electron temperature, and electron density in helium atmospheric pressure plasma jet, *Phys. Plasmas*, **19** (2012), 073513.
https://doi.org/10.1063/1.4739060

19. S. J. Doyle, K. G. Xu, Usof thermocouples and argon line broadening for gas temperature measurement in a radio frequency atmospheric microplasma jet, *Rev. Sci. Instrum.*, **88** (2017), 023114.
https://doi.org/10.1063/1.4976683

20. C. Yubero, A. Rodero, M. S. Dimitrijevic, A. Gamero, M. C. García, Gas temperature determination in an argon non-thermal plasma at atmospheric pressure from broadenings of atomic emission lines, *Spectrochim. Acta, Part B*, **129** (2017), 14-20.
https://doi.org/10.1016/j.sab.2017.01.002

21. P. Moon, E. Spencer, *Field Theory Handbook*, Springer-Verlag, Berlin, Heidelberg, 1988.
https://doi.org/10.1007/978-3-642-83243-7

22. G. Nellis, S. Klein, *Heat Transfer*, Cambridge University Press, Cambridge, 2012.
https://doi.org/10.1017/CBO9780511841606
23. E. W. Hobson, *The Theory of Spherical and Ellipsoidal Harmonics*, Chelsea Publishing Company, New York, 1965.

24. I. N. Sneddon, R. P. Srivastav, Dual series relations I – Dual relations involving Fourier-Bessel series, in *Proceedings of the Royal Society of Edinburgh*, **66** (1963), 150–160. https://doi.org/10.1017/S0080454100007809

25. T. Theodoulidis, A. Skarlatos, Efficient calculation of transient eddy current response from multilayer cylindrical conductive media, *Phil. Trans. R. Soc. A*, **378** (2020), 20190588. https://doi.org/10.1098/rsta.2019.0588

26. A. Ratsakou, A. Skarlatos, C. Reboud, D. Lesselier, Shape reconstruction of delamination defects using thermographic infrared signals based on an enhanced Canny approach, *Infrared Phys. Technol.*, **111** (2020), 103527. https://doi.org/10.1016/j.infrared.2020.103527

27. T. Weiland, Time domain electromagnetic field computation with the finite difference methods, *Int. J. Numer. Modell. Electron. Networks Devices Fields*, **9** (1996), 295–319. https://doi.org/10.1002/(SICI)1099-1204(199607)9:4<295::AID-JNM240>3.0.CO;2-8

28. R. Marklein, The finite integration technique as a general tool to compute acoustic, electromagnetic, elastodynamic, and coupled wave fields, in *Review of Radio Science*, (1999), 201–244. Available from: https://www.researchgate.net/publication/228540772.

29. A. Ratsakou, C. Reboud, A. Skarlatos, D. Lesselier, Fast models dedicated to simulation of eddy current thermography, in *Electromagnetic Nondestructive Evaluation (XXI)*, **43** (2018), 175–182. https://doi.org/10.3233/978-1-61499-836-5-175

<figure><img src="image_7.png" /></figure>

AIMS Press

©2023 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)