Modelling of Argon Cold Atmospheric Plasmas for Biomedical Applications

Plasmas for biomedical applications are one of the newest fields of plasma utilization. Especially high is the interest toward plasma usage in medicine. Promising results are achieved in blood coagulation, wound healing, treatment of some forms of cancer, diabetic complications, etc. However, the investigations of the biomedical applications from biological and medical viewpoint are much more advanced than the studies on the dynamics of the plasma. In this work we aim to address some specific challenges in the field of plasma modelling, arising from biomedical applications – what are the plasma reactive species’ and electrical fields’ spatial distributions as well as their production mechanisms; what are the fluxes and energies of the various components of the plasma delivers to the treated surfaces; what is the gas flow pattern? The focus is on two devices, namely the capacitive coupled plasma jet and the microwave surface wave sustained discharge. The devices are representatives of the so called cold atmospheric plasmas (CAPs). These are discharges characterized by low gas temperature – less than 40°C at the point of application – and non-equilibrium chemistry.


Introduction
Plasmas for biomedical applications are one of the newest fields of plasma utilization. Particularly plasma use in health care dates from about the beginning of this century. Promising results are achieved already in blood coagulation, wound healing and treatment of several forms of cancer [1]. There are indications that plasma medicine could offer solutions for lesions that cannot be treated otherwise, like diabetic complications and other severe ulcers [2].
However, the investigations of the biomedical applications from biological and medical viewpoint are much more advanced than the studies on the dynamics of the plasma -the configuration effect, i.e. the practically appropriate design and its effect on the plasma properties, the chemical composition, the transport of species, etc. In this work we aim to address some specific challenges in the field of plasma modelling, which arise from biomedical applications. Three main problems stand out in this research area: what are the plasma reactive species' and electrical fields' spatial distributions and what are the production mechanisms; what are the fluxes and energies of the various species that the plasma delivers to the treated surfaces (cells, tissues, etc.); what is the gas flow pattern?
The recent and most widely used, especially for medical applications, devices are the so called cold atmospheric plasmas (CAPs). These are discharges characterized by low gas temperature -less than 40C at the point of application -and non-equilibrium chemistry. Typically they are RF apparatuses [3]. Microwave devices working at higher gas flow (jets) can also be employed [4]. So far, an extensive theoretical and experimental study has been done only on the plasma needle [5]. Various other RF atmospheric plasma torches with different design are tested for multiple applications but not well studied and characterized (see for instance [6] and the references therein). In this study we focus on two of these devices, namely the capacitive coupled plasma (CCP) jet shortly denoted by plasma shower [7] and on a plasma device, which has barely been investigated for biomedical applications -the surface wave sustained discharge (SWD) [8]. Both of the considered devices operate in argon at atmospheric pressure (p = 1 bar) as it is shown in figure 1. However, the difference in the frequency of the power coupling mechanism induces a big difference in plasma properties.

The SWD Model
The plasma in a surface wave discharge is created and sustained by an electromagnetic wave. The wave is generated by a microwave source. It is with a frequency of 2.45 GHz. The wave propagates along the interface between the plasma and the dielectric tube (with relative permittivity ε d ) keeping the maximum of the field intensity at the plasma surface; hence the name of the wave. The surface electromagnetic wave creates the plasma and at the same time the plasma is a part of the waveguide structure for the wave propagation. The wave gradually transfers its power to the plasma, more precisely to the electrons. The power decreases with the distance from the launcher and so does the plasma density. The energy the electrons absorb from the wave is expended on different elementary processes which sustain the plasma. When the absorbed power is not enough to compensate the losses associated with radiation, diffusion and recombination the discharge ends.
This type of plasma can be sustained in a glass tube or directly in open air as a torch. Figure 2 shows schematic representation of the discharge and the cross sections of the two geometries, which have axial symmetry. A physically sound theoretical model for a discharge created by a surface wave require account in an interconnected manner of both aspects of the problem, namely, the electromagnetic and the plasma i.e. kinetic aspect. The former accounts for the change of wave characteristicsdispersion, spatial distribution of the field components, energy dissipation, while the later takes into account all the processes that ensure the balance between the creation and destruction of particles in the plasma as well as the energy balance of the particles. The full model provides a complete self-consistent description of plasma and wave features. The development of such a model, however, is a rather complex task. Hence usually at some stage of the work it is necessary to adopt some simplifying assumptions. In our model it is assumed that the plasma is in a steady state and radially averaged plasma density is used, building one-dimensional axial model.

The CCP Model
A model of an RF capacitively coupled discharge with the form of a shower has been developed and used to obtain key plasma parameters as well as to gain insight into the physical processes that determine the plasma behaviour. The configuration we are aiming at should mimic the plasma shower as well as possible. However this shower is intrinsically 3D while our model is 2D. To describe reality in a reasonable way we replaced the shower plates by two parallel plates equipped with slits. These slits are shifted with respect to each other in order to represent the fact that the shower holes are not in line with each other. The width of these slits is the same as the diameter of the holes. The flow rate is selected such that the gas velocity through these slits is the same as that through the shower holes.
In the modelling the device is considered as a periodic 2D Cartesian structure. The structural unitcell, shown in the figure, includes half of a hole and half of the metal surface between two holes from the powered electrode on the left, the plasma bulk in the middle and a similar electrode pattern but displaced in vertical direction for the grounded electrode on the right. On the top and the bottom of the unit-cell mirror boundary conditions are applied. In the figure, l p = 2 mm, denotes the distance  between the electrodes, w h /2 = 0.5 mm is the half width of the hole, w e /2 = 1.0 mm is the half width of the electrode between two holes, d = 1 mm is the depth of the holes. Three different sub-studies were carried out in order to explore step-by-step the complex cooperation between the different mechanisms: A) a flowless classical parallel plates' geometry, B) a flowless perforated parallel plates' configuration and C) the plasma shower. With this three-stage approach we were able to investigate the role of the sheaths and space charge, the influence of the shower holes on the field distribution and the effect of the flow on the distribution of species [4].

Modelling approaches comparison
We shell start with a comparison of the main features of the shower and the SWD and a discussion how that affects the corresponding modelling approaches.
Capacitively coupled plasmas (CCP), like the plasma shower, work on electrostatic principles; the functioning of a surface-wave discharge, on the other hand, is based on electrodynamics, i.e. the propagation of EM waves. Thus for a CCP Poisson's law is essential, while the description of the electromagnetic behaviour of SWD is based on the Maxwell equation.
The shower is an electrode-driven discharge; the surface-wave plasma on the other hand is electrodeless. In the RF discharges the sheath next to the electrodes plays an important role in accelerating ions and creating secondary electrons upon collisions with the cathode. This implies that space charges play an important role in the shower plasma, while the surface-wave plasma can be treated as being quasi neutral. The existence of electrodes in the plasma shower implies that surface processes such as secondary electron emission and/or sputtering have to be taken into account. This is not needed for the surface-wave plasma.
In a CCP the electron density and the sheath voltage follow the oscillations of the externally applied electric field and change both in time and space. In a SWD the wave frequency ω is smaller than the electron plasma frequency ω p in the bulk plasma, ω < ω p . Since the time scale for the sheath formation is ω p -1 , it follows that the sheath in a SWD remains largely unchanged during one time period of the oscillating EM field. For modelling this means that in a SWD one can consider the electron density as being constant in time [3] so that a stationary model can be used. Moreover we can work in a SWD with space averaged plasma quantities. On the other hand the CCP requires a timedependent model and for the plasma shower we need at least a two-dimensional geometry to get an adequate description of the plasma in and around the shower holes.
At the industrial driving frequencies of the SWD (typically 2.45 GHz) and the required electron densities (typically 10 20 m -3 ) the wavelengths and skin-depths are comparable or shorter than the typical equipment size. This is not so in the classical CCP that is usually operated at 13.56 MHz. Thus, a uniform amplitude approximation of the oscillating electromagnetic (EM) field of the SWD becomes impossible while in RF discharge it is practically the case.
In the shower the convective flow is essential to bring the active species created in between the plate to the region downstream. In the SWD modeling the gas flow rate can be neglected in many cases.
In both model structures we follow the subdivision of plasma features into three main blocks: configuration, transport and gas discharge chemistry. Configuration deals with the impact of the environment on the plasma and thus, among others, with the shape and sizes of the plasma, the boundary conditions and energy coupling modus. Transport describes the transport of species, momentum and energy in the plasma (and afterglow). Chemistry deals with the creation and destruction of plasma species; hence it has an important impact on the sources and the coefficients that drive and facilitate transport.
A schematic representation, giving the similarities and the differences between the model blocks is given in table 1. The numerical approach is also built by means of three main units. They are denoted as chemical kinetics, plasma dynamics and fluid dynamics modules. A schematic representation of the numerical procedure is shown in the figure. The chemical kinetics module, based of Boltzmann solver combined with CRM, provides the necessary data of reaction and transport coefficients. In the case of the plasma shower we work with 4 active species handled by the Boltzmann solver Bolsig+, developed at Toulouse. For the description of the SWD a solver, called Sofie, developed in the group of Sofia University is used. It consists of an EEDF solver linked to a CRM.
The plasma dynamics of the shower plasma and SWD is also treated with different software packages. MD2D, a time-dependent two-dimensional fluid model developed at the Technical University of Eindhoven and part of the PLASIMO modelling platform, deals with the configurationtransport aspects of the RF plasma. The SWD dynamics is handled by electrodynamic description of the wave propagation that is self-consistently linked with the kinetic model. In order to account for the effect of the collisions on the wave propagation the full expression for the plasma dielectric permittivity, resulting in complex algebra calculations, is considered.
Regarding the fluid dynamics we calculated the convective velocity entering the convection term in the drift-diffusion-convection equation of MD2D with the Navier-Stokes solver of the PLASIMO platform. As stated before there is no convection modelling applied to the SWD.

The surface-wave plasma
The investigation of the SWD starts with a study of the propagation characteristics of the surface wave maintaining the plasma [8]. We have analysed the wave dispersion relation for azimuthally symmetric wave at arbitrary ratio of the electron-atom collision frequency and the wave frequency. The analysis includes propagation and damping or attenuation diagrams. The formation of the propagation-damping diagrams for fixed ratio of  ea /ω (corresponding to a given atmospheric pressure plasma column) using particular dispersion curves is demonstrated in figure 5 (a) and (b). In the case of an exciter that launches a wave with constant frequency in non-homogeneous plasma, the dispersion relation (weakly) depends on the axial coordinate via plasma density and thus ω p . It must be satisfied at each point along the column. So, in contrast to the classical dispersion, the wave propagation characteristics cannot be described by means of a constant ω p . Instead, the wave excitation frequency is constant (denoted as  * in this case). The relation between the (changing)  p and the corresponding k can be presented by means of a phase (propagation) diagram. Just as in the case of the dispersion diagram it will give  / p versus R k . But the variation in  * / p, is now due to a change in  p , not in .

Collisional-radiative model
In each dispersion and attenuation curve ( figure 5 (a) and (b)) there is one point corresponding to wave frequency 2.45 GHz. These points are shown by arrows and they belong to the phase (propagation) and damping diagrams of inhomogeneous plasma sustained by the surface wave with frequency 2.45 GHz. Connecting all such points one actually obtains the phase (propagation) diagram thick black line in figure 5 (a) and damping diagramthick black line in figure 5 (  ) is found to be much higher than that corresponding to the turning-back point n e ( k max , [9]) considered by these authors as a criterion for the column end at the same conditions. The electron density n e ( k k     ) is also much higher than the resonance n e,res (  √ ), corresponding to the plasma density at the end of low pressure SWD.
The spatial distributions of the components of EM field sustaining a SWD are obtained and the influence of the electron-ion collision frequency  ea on the wave components radial distribution is studied. It is found out that the increase of  ea facilitates the penetration of the wave into the plasma and diminishes the role of the dielectric tube [7].
(a) (b) Figure 6. Dependence of the radial distribution of the axial (a) and radial (b) surface wave field components on the collision frequency for configuration plasma-dielectric-air -the solid lines and plasma-air -the dotted lines. By means of the complete model the axial distributions of all discharge characteristics can be delivered for the two configurations -plasma inside the dielectric tube and a torch. Figure 7 presents the distributions from the launcher toward the end of the discharge (z = 0 m) along the plasma column of the electron density and temperature, the electron-atom collision frequency and the wave power sustaining the discharge. The red points in figure 7 correspond to the actual end of plasma where the damping coefficient is equal to the propagation coefficient and becomes higher after that points. The gas temperature is considered a fixed parameter and is taken to be T g = 300 K. The internal radius of the tube (or the plasma radius) equals R = 1 mm, the quartz wall thickness is d = 1 mm and the relative high frequency dielectric permittivity is ε d = 3.8.
In contrast to low and intermediate pressure, we find for atmospheric pressure SWDs that the wave energy is not always completely dissipated in the plasma ( figure 7(c)). The reason is that the plasma sustaining by the wave ends when the damping coefficient becomes higher than the propagation coefficient. At that region the remaining wave power is not enough for sustaining the plasma but the wave can still propagate along the waveguide structure, at least theoretically. The plasma density n e (figure 7(a)) decreases from the wave launcher to the column end remaining quite high at the end (red points) which is also experimentally observed [10]. Figure 7. Axial distribution of electron density (a), electron temperature (b), electron-atom collision frequency (c) and the wave power sustaining the plasma (d). The red points mark the calculated end of the column corresponding to It is also found that T e ( figure 7 (b)) increases approaching the end of the plasma column. The reason is that towards the plasma end the ionization ratio  = n e /N is lower and at lower  the tail of the EEDF is more depleted. A depleted tail demands for an enhanced bulk of the EEDF in order to keep the ionization ongoing, which explains the increase of the mean energy of electrons. Similar is the explanation why in the tube plasma where n e is higher electrons acquire lower T e and  ea (compare figure 7(a) with (b) and (c)).
Finally the axial variation of the atomic states axial distribution is given in figure 8. The predominant species are the molecular ions Ar 2 + which are an order higher than the most populated Ar state Ar(4s) followed by Ar 3 + and Ar + . The density of argon dimmer is also substantial, n(Ar 2 * ) ~ 10 19 m -3 .

The CCP plasma
The study of the classical parallel plate configuration CCP reveals that strongly non-Maxwellian kinetics is manifested. It results in lowering the direct ionization and a promoting the stepwise ionization. This leads to a two-loop chemistry: 1) the excitation-radiative decay loop and 2) the metastable ionization-recombination circuit. The consequence of this two-loop chemistry and the underlying reaction rates is that the density of the molecular ion is much higher than that of the atomic ion n(Ar 2 + ) >> n(Ar + ) and also that the excited species are much more abundant than the charged particles n(Ar * ) >> n(Ar 2 + ) as can be seen from figure 9.
(a) (b) Figure 9. Spatial distributions of the densities of (a) electrons and molecular ions, (b) atomic ions and excited atoms at maximum negative potential on the powered electrode Figure 10. Spatial distributions in m -3 of the densities of (a) electrons, (b) atomic ions, (c) molecular ions and (d) excited atoms in a plasma shower without flow. All snapshots are taken for maximal negative potential on the powered electrode. Note the difference in the scales.   Figure 11. Spatial distributions in m -3 of the densities of (a) electrons, (b) atomic ions, (c) molecular ions and (d) excited atoms in a plasma shower with flow. All snapshots are taken for maximal negative potential on the powered electrode. Note the difference in the scales.
Finally we find that the flow plays a crucial role: the comparison between figure 10 and figure 11 reveals this. As the most effective outward transport mechanism it is the basis of the formation of the post-discharge region, the region that is of most interest for applications. However as the time scale of the chemical processes is smaller than that of transport even with a flow we get no Ar + flux in the post-discharge region (see figure 11(b)). In addition to modifying the distribution of the plasma particles in the spatial afterglow the flow recirculation reduces the wall losses and facilitates the ionization and excitation processes in between the plates. This results in higher densities of the plasma constituents there and more homogeneous distribution. Thus it plays a key-role in the stabilization of the discharge.

Conclusions
Two CAP devices, the plasma shower and the surface-wave discharge have been investigated in view of their biomedical applications. The main focus was on the three main problems that stand out in this field: the plasma reactive species' and electrical fields' spatial distributions, the fluxes and energies of the various species which the plasma delivers to the treated surfaces and the gas flow pattern. At the same time a discussion on the differences in the modelling approaches was done.
The SWD is found to maintain plasma with higher density n e ~ 10 20 m -3 while for the CCP it is n e ~ 10 18 m -3 . The predominant species in surface-wave discharge are the charged ones while in the plasma shower they are the excited atoms. Though with lower than the charged particles density the various excited species in the SWD are substantial n(Ar * ) ~ 10 17 ÷ 10 17 m -3 . The spatial distributions of the plasma constituents and their mean energies for both devices are obtained. The mean energy of the electrons in the SWD is found to be about 1-2 eV versus 3.5 eV in the CCP. The two high frequency excitation fields and the different creating mechanisms related to them are investigated and linked to the different functioning and properties of the discharges. The role of the flow in the CCP is extensively studied and found to play a two-fold role -as the transport mechanism creating the discharge application zone and as a stabilizing factor. The flow in the SWD is not investigated at this stage.