Numerical methods are used to determine the Electron Energy Distribution Function (EEDF) from *I*(*V*) probe characteristics, which are measured using a cylindrical Langmuir probe in the case of weakly ionized plasmas. This task becomes difficult when measurement is complicated by the presence of an external magnetic field or in high pressure plasma because of collision between electrons and heavy particles within the sheath formed around the probe tip. In this case, the electron current must be calculated using the Swift law instead of the Langmuir law. The numerical methods consist of determining the derivative functions of the *I*(*V*) probe characteristics in the case of a noisy signal and correcting the EEDF taking into account the electron diffusion coefficient within the sheath formed around the probe collector. Algorithms are given to detail the methods step by step, which can be used to write homemade codes. The methods are tested in the case of different plasma reactors described in the literature, such as microwave plasma and rf (radio-frequency) and dc (direct current) plasma reactors working at different pressures with or without magnetic field. The results show the effect of pressure or magnetic field on the *I*(*V*) probe characteristics because of the change in the electron diffusion coefficient.

## I. INTRODUCTION

Langmuir probe diagnostics provide useful information about plasma parameters (electron energy or density, ion density, plasma potential or floating potential, etc.) in both basic research and plasma processing.^{1–5} This simple method is widely used to investigate the plasma under different experimental conditions (low and high pressure plasma, thermal plasma, or magnetized plasma). The probe is biased and immersed into the plasma. It is used to measure the current–voltage characteristics of the surrounding plasma in the vicinity of the probe collector. Unfortunately the theoretical analysis of the *I*(*V*) probe characteristics is generally complicated and depends on the experimental conditions and on the type of probe used during the experiment. The *I*(*V*) probe characteristic can be drastically impacted by the probe collector geometry or by the sheath formed around the collector because of the diffusion of charged particles, collision processes within this sheath, external effects such as a magnetic field used to confine the plasma. Theories have been developed taking into account different plasma configurations to determine reliable results,^{6–10} depending on models and hypotheses, but they are available in specific cases only.

The determination of electron energy distributions in plasma is a crucial aspect in plasma physics. This distribution influences the characteristics and the behavior of the plasma, which in turn impacts its applications in material processing, plasma-based propulsion, and semiconductor manufacturing.^{11–13} However, the interpretation of the data acquired from the probe cannot be challenged without numerical methods. In this article, we will discuss how numerical methods can be used to accurately determine the electron energy distribution from Langmuir probe measurements, in different plasmas.

In order to simplify the study, we consider the case of a cylindrical probe collector, which is the most frequent geometry used in weakly ionized plasma. The case of totally ionized plasma is not considered in this article since other phenomena such as the effect of the ion viscosity on the collected current, which complicate the interpretation, should be taken into account.^{14}

## II. THEORY

The basic theory, developed by Druyvesteyn,^{15} shows the dependence of the I(V) probe characteristics on the Electron Energy Distribution Function (EEDF) within the plasma surrounding the probe. In the following part, details of the theories used to determine the EEDF are given, with the different hypotheses used to elaborate the laws, and consequently, these laws can be used under specific conditions only.

Considering a cylindrical probe negatively biased at a voltage (−*V*_{P}) (in reference to the plasma potential) and looking at the electron current impinging on the probe collector, the electrons with velocity **v** contained in the elementary volume of plasma *dV* and impinging on the elementary collector surface δ*A*_{p} in 1 s are contained in the hemispherical volume centered on δ*A*_{p} and of radius *v* (see Fig. 1).

*dV*in cylindrical coordinates is given by

Considering the Electron Velocity Distribution Function *f*(*v*), the number of electrons contained in *dV* with a velocity ranging from *v* and *v* + *dv* is *dn*_{e} = *f*(*v*)*dV*. Electrons of velocity **v** impinging on δ*A*_{p} are contained in the solid angle $d\Omega =\delta Ap\u2061cos\theta v2$, and their number is given by $dne=f(v)dVd\Omega 4\pi $. It corresponds to the elementary current intensity *di*_{e} = *evdn*_{e}. Hence, the total current collected on δ*A*_{p} is $ie=\u222bvmin\u221e\u222b\u0278=02\pi \u222b\theta =0\theta maxdie\theta ,\u0278,vd\theta d\u0278dv$.

In the following part, values of *i*_{e}, *e*, and *V*_{p} are arbitrarily considered positive in the different equations.

The electrons impinging on δ*A*_{p} are collected if the kinetic energy $Ee=12mevcos\theta 2\u2265eVP$ and $\theta max=2eVPmev2$.

The minimum velocity value *v*_{min} is obtained when *θ* = 0 and is given by $vmin=2eVPmev21/2$.

*A*

_{p}is

^{15}It is

*g*(

*ε*

_{e}), $g\epsilon e=12me\epsilon e1/2fv$.

*F*(

*x*,

*y*) and the integral $u=\u222b\phi (x)f(x)F(x,y)dy$, the first derivative function is given by

*u*=

*i*

_{e},

*x*=

*V*

_{P}, and

*y*=

*v*, the first derivative function vs

*V*

_{p}of the Langmuir law is given by

_{e}=

*eV*

_{P}.

Equations (5) and (6) are the Druyvesteyn Equations used to determine the Electron Velocity Distribution Function (EVDF) or the Electron Energy Distribution Function (EEDF) from the electron current intensity collected by a negatively biased probe of area δ*A*_{P}.^{15}

Details of the theory given above show that both the Langmuir law and the Druyvesteyn law are available only in the case where electrons are not disturbed in the sheath formed around the probe. The energy of electrons collected by the probe depends on the initial electron energy and on the biased voltage only. Hence, these equations can be used when the disturbance of the plasma due to the probe can be neglected, i.e., when the electron diffusion rate is high enough to cope with the drain of electrons from the plasma to the probe surface. If the diffusion coefficient is too low, the measured distribution function will be underestimated compared to the real distribution function of the undisturbed plasma. The disturbance is greater when the gas pressure and the probe area increase and as the electron velocity decreases. This effect is also observed in magnetized plasma because of the electron gyration within the magnetic field. In all these cases, corrections of the Langmuir and Druyvesteyn equations are necessary, taking into account the electron diffusion within the sheath formed around the probe. This has been done by Swift.^{16}

*d*Φ(

*ε*

_{e}) of energy ranging from

*ε*

_{e}and

*ε*

_{e}+

*dε*

_{e}does not only depend on the biased voltage but also on the electron diffusion through the sheath surrounding the probe collector. Assuming a radial diffusion through the sheath, the electron flow is given by

*r*,

*f*(

*ε*

_{e}), and

*D*(

*ε*

_{e}) are the sheath radius, the EEDF, and the electron diffusion coefficient, respectively. Consequently, the EEDF change is due to the diffusion through the sheath and is $f\epsilon e\u2212f0\epsilon ed\epsilon e=di\epsilon e\u222br0rdreD\epsilon eS$, where $di\epsilon e=eSd\Phi \epsilon e$ and

*r*

_{0}is the sheath thickness.

From this equation, it can be seen that the larger the diffusion coefficient, the lower is the change in the EEDF.

*r*

_{p}biased at the potential

*V*

_{p}with a sheath radius

*r*

_{s}is now

^{16}

*ε*

_{e}) is called the diffusion parameter and

Equation (7) is the modified form of the Langmuir law (2) proposed by Swift,^{16} taking into account the diffusion coefficient of the electron through the sheath. When the electron diffusion through the sheath can be neglected, ψ(*ε*_{e}) = 0, and consequently, *i*_{e}(*V*_{P}) is given by Eq. (2). The Swift law [Eq. (7)] has been used by different authors^{6–10} to study the effect of magnetic fields or moderate-collisional plasma on I(V) probe characteristics. These authors have developed models assuming hypotheses and considering the product of the electron velocity and the diffusion coefficient, *v*_{e}(*ε*_{e})*D*(*ε*_{e}), constant within the sheath thickness.^{8} Consequently, the diffusion parameter Ψ(*ε*_{e}, *V*) = Ψ(*ε*_{e}). It does not depend on the probe bias voltage. For a probe tip of radius *r*_{p} and length l, the Druyvesteyn method is reasonably applicable if $\lambda e>0.75rpln\pi l4rp$, and if $\lambda e<rp7ln\pi l4rp$, the EEDF can be calculated using the first derivative function of the electron current *di*_{e}/*dV*.^{6,8} Hence, it is necessary to adjust the probe radius to satisfy the corresponding validity conditions. By this way, authors consider a mean value for the electron mean free path all over the EEDF. To prevent error due to these assumptions and approximations, the following part gives a numerical method that does not need any additional hypotheses.

*ε*

_{e}) can be calculated using the ratio A/B,

^{17,18}

It is worth noting that in Eq. (10), because $diedVp=\u2212diedVapp\u22640$ and $d2iedVp2=d2iedVapp2\u22650$, the ratio $die/dVpd2ie/dVp2\u22640$.

The diffusion parameter is expected to be positive. Otherwise, the electron flow should increase when it reaches the collector. The condition $\Psi \epsilon e\u22650$ is fulfilled only if $die/dVpd2ie/dVp2\u2264\u2212\epsilon e\u2212eVpe$.

These results show that the diffusion parameter can be calculated directly using the ratio of the first to the second derivative functions of the experimental *I*_{e}(*V*_{P}) probe characteristic without any additional model. However, the value of Ψ(*ε*_{e}) depends on the noise measured with the experimental signal during the data acquisition, and it is necessary to determine accurate values of the derivative functions. The larger the signal/noise ratio, the better the determination of Ψ(*ε*_{e}).

The next part of this article presents numerical methods to determine the EEDF from the *I*(*V*) probe characteristics considering the general case of the electron current given by Swift Eq. (7).

## III. DETERMINATION OF THE SECOND DERIVATIVE FUNCTION

The EEDF can be calculated using the values of the first and second derivatives of the electron current as a function of the applied potential (in reference to the plasma potential). Different methods can be used to determine the derivative functions. The easiest should be to measure the change in *i*_{e} vs applied voltage _{Vapp}. This method was used by Medicus,^{19} but it is generally not efficient because of the low signal/noise ratio. “Analog differentiation” methods can also be used.^{20,21} However, additional equipment such as active differentiators composed of cascading operational amplifiers is needed, which generally decreases the apparatus response time necessary for a high energy resolution.^{22} The second derivative function of the probe measurement is not a direct value, but it results in a convolution product between the “true” second derivative and the instrumental (or apparatus) functions. It can have an effect on the results as a shift of the plasma potential position or a broadening of the derivative function shape. A correction is necessary, which can be done with the knowledge of the instrumental function only.^{8,23} Numerical methods can also be used to smooth data by suppressing the noise of the recorded signals. These methods do not need additional electronic equipment. The data smoothing method proposed by Hayden^{24} or the optimal Wiener filtering method^{25} considers the convolution product between the signal and the apparatus function of the probe. The peculiarity of this convolution product is that the noise spectrum is statistically uncorrelated with the instrumental function. However, these methods need the knowledge of the apparatus function of the setup, which is also generally unknown. In the case of a Langmuir probe, it depends on the sheath formed around the probe, which changes with the biased voltage. Hence, these numerical methods that involve the knowledge of the apparatus function are also not easy to use.

Other numerical methods based on polynomial interpolation, such as the Lagrange polynomial interpolation, the cubic spline interpolation method,^{25} or the Savitzky–Golay algorithm,^{26} are more suitable to smooth the *I*(*V*) probe characteristic. These methods do not depend on apparatus functions. Moreover, the results can be improved by gradually changing the polynomial parameters, and any change in the curve (shift of the plasma potential value or change in the derivative function shape) can be detected. Hence, any error due to the smoothing effect can be prevented.

The following part presents two efficient numerical methods based on polynomial interpolations and Fourier transforms. We compare these methods in the case of different experiments proposed in the literature.

### A. The smoothing differentiation method (SDM)

The method was first used to determine the second derivative of the probe characteristic by Fujita *et al.*.^{27} This method is derived from the Savitzky–Golay algorithm and is based on the least squares principle to fit the experimental data with a polynomial curve and on the differential operation by convolution with weighting functions. The step $h=xi\u2212xi\u22121$ between two successive points $xi,yi$ must be constant. It is assumed that the smoothed value *y*(*j*) in the vicinity of a point i can be estimated by a quadratic polynomial $yj=a2j\u2212i2+a1j\u2212i+a0$, where *j* = −*m* + *i*, …, −1 + *i*, *i*, *i* + 1, ….*i* + *m*. *m* is a suitable integer selected as the fitting width.

*a*

_{0},

_{1,2}are determined considering that the squared error

*I*between the estimated value

*y*(

*j*) and the measured value

*y*

_{0}(

*j*) is minimized,

Consequently, the fitted value $yi=a0$, the first derivative $y\u2032i=\u2202yi\u2202i=a1h$, and the second derivative $y\u2033i=\u22022yi\u2202i2=2a2h2$. This method consists of considering Taylor development of the second order for each experimental point *f*(*x*):

$fx=fx0+hf\u2032x+h22f\u2033x++R(h)$, where *h* is the step between two points (*j* − *i*), which must be constant, and *R*(*h*) is the rest obtained after minimization of squared error I. The smaller the step between the two points, the more precise are the *a*_{0}, *a*_{1}, and *a*_{2} values.

Figure 2 shows the SDM algorithm.

This method has been tested to determine the second derivative function of the *I*(*V*) probe characteristics under different experimental conditions. Figure 3 shows the second derivative function of the *I*(*V*) probe characteristic measured in argon microwave plasma working at a frequency of 2.45 GHz, a pressure of 2 Torr, and an incident power of 100 W. The results are obtained using different m values corresponding to a fitting width of 2*m* + 1 values. It can be seen that the noise decreases with increasing m value, and good results are obtained for m larger than 10. However, there appears a broadening of the curve, with increasing m producing a slight shift of the potential plasma value determined at the zero-crossing point with the *V* axis when the second derivative becomes negative (indicated by the vertical arrow). The plasma potential value shifts from 13.7 V for *m* = 5 to 14.4 V for *m* = 20. This is the main problem observed with this method. This shift decreases with increasing density of points in the data array. Hence, the m value should be gradually increased to control any shifting or broadening effect and to check the results of integrating the second derivative and compare the results to the initial experimental *I*(*V*) curve. Thus, any shift in the plasma potential or broadening of the second derivative will be detected.

### B. The simulation of harmonic components (SHCs)

This second numerical method can be compared to the “analog method” of Kortshagen and Schlüter,^{20} but it does not involve additional apparatus. Hence, there is no need of the knowledge of the apparatus function of new equipment. The principle of the method consists in using a Fourier transformation to isolate harmonic components, which are obtained based on the original experimental *I*(*V*) probe characteristic by the simulation of the effect of a sinusoidal superimposed perturbation on the measured signal. This method is efficient in the case of noisy signals.^{28–31}

Let us consider the part of the experimental *I*(*V*) probe characteristic corresponding to the negatively biased probe with respect to the plasma potential. In this case, electrons are repulsed and only electrons of kinetic energy larger than *eV*_{p} (*V*_{p} is the biased voltage) are collected.

Suppose that a sinusoidal component $ut=u0sin(\omega t)$ is added to the applied voltage, *V*(*t*) = *V*_{P} + *u*(*t*).

*I*(

*V*

_{P}) curve as infinitely derivable, the Taylor expansion gives

*I*

^{i}[

*V*

_{p}+

*u*(

*t*)] is the ith derivative of

*I*[

*V*

_{p}+

*u*(

*t*)] with respect to

*t*and

*u*

_{0}

^{n}

*ε*[

*u*(

*t*)] is the rest, $limn\u2192\u221eu0n\epsilon ut\u21920$.

This equation shows that the second derivative appears in the term related to the second harmonic 2ω. Isolating these components, it is possible to obtain a signal whose amplitude is proportional to the second derivative of the electron current *I*(*V*_{p}).

*I*(

*V*

_{p}) probe characteristics are a data array [(

*x*

_{1},

*y*

_{1}), …, (

*x*

_{n},

*y*

_{n})], corresponding to n acquisition points. We simulate for each of these points (

*x*

_{i},

*y*

_{i}) using a sinusoidal signal $ut=u0sin(\omega t)$ with

*w*= 2

*πf*, which is added to the applied voltage

*x*

_{i}=

*V*

_{P}. We obtain for each point a new array containing

*n*

_{s}points [(

*x*

_{i}+

*u*(

*t*)], [

*I*

_{p}(

*x*

_{i}+

*u*(

*t*)]. The value of the current intensity

*I*

_{p}[

*x*

_{i}+

*u*(

*t*)] can be deduced from the experimental

*I*(

*V*

_{p}) curve using a Lagrange polynomial interpolation of degree (

*m*− 1) through m consecutive points (

*x*

_{k},

*y*

_{k}) of the initial experimental data array for

*x*=

*x*

_{i}+

*u*

_{0}sin(

*ωt*) ranging from

*x*

_{k=1}to

*x*

_{k=m}, i.e., $xi+uo<xm$ and $xi\u2212uo>x1$, with

*x*

_{n}>

*x*

_{n−1},

^{25}which gives

*m*= 4.

Hence, at each point (*x*_{i}, *y*_{i}) of the initial data array, a new data array containing *n*_{s} points corresponding to the simulated probe characteristic is determined, which would be obtained when a sinusoidal superimposed component is added to the initial applied potential (*V*_{p}). This new signal consists of multiple harmonic components, and one of them is the second harmonic, which is proportional to cos(2*wt*) = cos(4*πft*). Using the direct Fourier transform, this component can be isolated in the frequency domain. Using the inverse Fourier transform and dividing it by $14u02\u2061cos4\pi ft$, we obtain *n*_{s} points corresponding to $xi=VP,yi\u2033=\u22022IP\u2202VP2$. This procedure must be done for any (*x*_{i}, *y*_{i}) points of the initial experimental data array, and the output is $xi=VP,yi\u2033=1ns\u2211ns\u22022IP\u2202VP2$. The second derivative values calculated for the n_{s} points are close to each other, and the final value used is the mean value of these second derivative values calculated over the *n*_{s} points.

It is worth noting that in this calculation, the frequency f must be between $f\u2208\u2212fNy,fNy$ (where *f*_{NY} is the Nyquist frequency) to prevent any “aliasing effect.”^{25} The Nyquist frequency is given by $fNy=12\Delta $ with Δ = *t*_{ech}, where t_{ech} is the arbitrary value of the time between two consecutive points [(*x*_{i} + u(*t*)], [*I*_{p}(*x*_{i} + *u*(*t*)] of the simulated new array.^{25}

Contrary to the previous method, this one does not need a constant interval between two successive points. Figure 4 shows the SHC algorithm.

Figure 5 shows the second derivative function calculated using this method for the same case previously seen for the SDM method on the same *I*(*V*) probe characteristic obtained in argon microwave plasma working at a frequency of 2.45 GHz, a pressure of 2 Torr, and an incident power of 100 W. Measurements are performed using the following parameters: sampling time tech = 10^{−6} s, *n*_{s} = 250 points, and frequency *f* = 242187 Hz (the Nyquist frequency is 0.5 MHz). Different *u*_{o} values are used and compared to the results previously obtained using the SDM method for *m* = 20. No broadening of the curve and shift in the plasma potential are observed with increasing *u*_{o}. It is worth noting that by increasing *u*_{o}, we can filter the noise contained in the data array. *u*_{o} has a filtering effect on the noise.

Concerning the rest of the Taylor expansion depending on (*u*_{o})^{n}, this one can always converge to 0 by changing the coordinate system to obtain *u*_{o} < 1. With the numerical method, it is always possible to define a new coordinate system for the *I*(*V*) probe characteristic after the data acquisition, so the second derivative can be calculated using a *u*_{o} value lower than 1. This means that the same results are obtained using *u*_{o} > 1 and *u*_{o}/*α* < 1 for the second derivative values. In the second case, the X-axis coordinates (*V* values) are divided by *α* and the Y-axis coordinates (second derivative values) are multiplied by *α*^{2}. Thus, the condition on the rest of the Taylor expansion $limn\u2192\u221euo/\alpha n\epsilon ut\u21920$ is fulfilled.

The filtering effect on the noise obtained with different *u*_{o} values can be seen in Fig. 6. The figure shows the spectral density measured in the frequency domain in the case of *I*(*V*) probe characteristic measured in the plasma “bubble.”^{32} The rate used to apply the bias voltage is 165.74 V/s (6.033 ms/V). This probe characteristic will be considered later in the text (paragraph V, application case C).

It can be seen that the large values of the frequency oscillations disappear with increasing *u*_{o}. For *u*_{o} = 0.01, oscillations are observed up to about 24 kHz; for *u*_{o} = 0.02, oscillations are observed up to about 13kHz; and for *u*_{o} = 0.5, only very low frequency remains (<2.5 kHz).

Regarding both methods (SDM or SHC), the best way to determine smoothed data is to gradually increase the parameter values (*m* and *u*_{o} corresponding to SDM and SHC methods, respectively) to prevent any distortion or shift in the second derivative function and then to integrate the second derivative function and compare the results to the experimental *I*(*V*) characteristic. Keeping in mind that the lower the step between two experimental points the more accurate the results are, it is worth noting that results obtained by one of the two previous methods (SDM or SHC) can be improved by fitting these results using the SDM method. In this case, the fitted second derivative value corresponds to parameter a_{o} (see details previously given).

In the following part, the SHC method is used to determine the EEDF from *I*(*V*) probe characteristics measured in different plasma proposed in the literature. The effect of electron diffusion through the sheath formed around the probe collector is considered.

## IV. DETERMINATION OF THE EEDF

As previously shown, the Druyvesteyn equation is available under collisionless conditions only. In this case, calculations using the Langmuir law [Eq. (2)] fit the experimental *I*(*V*) probe characteristic. Otherwise, the Swift law [Eq. (7)] is required. In these conditions, the EEDF is obtained taking into account the diffusion parameter by means of 7 steps:

The first and second derivative functions $\u2202ieV\u2202V$ and $\u22022ieV\u2202V2$ of the experimental

*I*_{e}(*V*) characteristic are calculated using one of the previous numerical methods, and the Druyvesteyn Eq. (6) is used to calculate $f0\epsilon e$, which corresponds to $\Psi \epsilon e=0$ in Eq. (7).$\Psi \epsilon e$ is calculated using Eq. (11) for each

*V*value of the data array.Using Eq. (10) and the previous values of $\Psi \epsilon e$ and $f0\epsilon e$, the new values of $\u22022ieV\u2202V2$ are calculated. Then steps (2) and (3) are repeated for all

*V*values of the experimental data array*I*_{e}(*V*).The values of $\u22022ieV\u2202V2$ calculated are compared to the experimental values (exact values) previously calculated in step (1). The difference $\u22022ie\epsilon e\u2202V2exp\u2212\u22022ie\epsilon e\u2202V2EQ(10)$ is due to the effect of the electron diffusion through the sheath.

- The EEDF $f0\epsilon e$ is incremented using the difference between the two second derivatives,(12)$\epsilon e\u22121/2f\epsilon e=\epsilon e\u22121/2f0\epsilon e+\u22022ie\epsilon e\u2202V2exp\u2212\u22022ie\epsilon e\u2202V2EQ(10).$
The new values of $\u22022ieV\u2202V2$ are calculated using Eq. (10) and the new value of the EEDF previously calculated Eq. (12). Again, $\u22022ie\epsilon e\u2202V2exp$ and $\u22022ie\epsilon e\u2202V2EQ(10)$ values are compared, and the increase in the EEDF using Eq. (12), by changing $f0\epsilon e$ to $f\epsilon e$, is still realized as long as the difference between both second derivative functions remains larger than an arbitrary limit. This happens generally after two or three iterations.

Finally, integration of the last EEDF, corresponding to the best fit $\u22022ie\epsilon e\u2202V2EQ(10)$ versus

*V*using Eq. (7), is done to check the good agreement between the calculated*I*(*V*) probe characteristic and the experimental values.

We have tested the method on different *I*(*V*) probe characteristics obtained in different plasma reactors and different experimental conditions: microwave plasma, without and with confining magnetic field, plasma “bubble,” and rf plasma with confining magnetic field. All the *I*(*V*) probe characteristics that we use have been published in the literature, except the first one, which has been measured in our own reactor. The results are given and discussed in the following part.

## V. APPLICATION TO MEASUREMENTS PERFORMED IN VARIOUS PLASMAS

Different points must be checked before data treatment to obtain reliable *I*(*V*) probe characteristics using Langmuir probes. One of the phenomena becomes important for the correct interpretation of the results: the influence of the collision of charged particles with neutrons in the space charged sheath surrounding the probe.^{33} Collision effects are important in plasma with a gas pressure larger than several 100 Pa. When the gas pressure is too high, a very large probe radius provides distortions of the second derivative, and the position of the plasma potential determined at the zero crossing point of the second derivative function can be changed. According to David *et al.*,^{33} the determination of the plasma potential using the second derivative can be done using a probe of radius *r*_{p} and length *L*_{p} if the condition $S=keke+1ke+lnLprp>1$, where *k*_{e} = *λ*_{e}/*r*_{p}.^{33}

According to Godyak *et al.*,^{22} reliable measurements by means of Langmuir probes must have a good energy resolution. This one is defined by the energy interval Δ*ε*_{e} between the zero crossing point of the second derivative function and the maximum of the EEDF. Information concerning low energy electrons can be lost if the energy resolution is too low (depending on the data acquisition system). The criterion for acceptable energy resolution is Δ*ε*_{e} < T_{e}. In practice, Δ*ε*_{e}/T_{e} < [0.3–0.5] corresponds to a good resolution.

Another important point concerning the probe measurements is the dynamic range, i.e., the ratio between the maximum value of the EEDF and the minimum value. According to Godyak *et al.*,^{22} a good dynamic range corresponds to 3–4 orders of magnitude. This corresponds to ε_{emax} ranging [7–9]T_{e}. In the case of a Maxwellian electron energy distribution function, this corresponds to a ratio $f\epsilon e\u222b0\u221ef\epsilon ed\epsilon e$, ranging from 4 × 10^{−4} to 3 × 10^{−3} assuming kT_{e} = 1 eV. This condition seems excessive, and a dynamic range corresponding to 2–3 orders of magnitude could also be appropriate and more easily accessible.

Small energy resolution and a good dynamic range in measurements with control of probe contamination give confidence in data acquisition.^{22}

### A. Expanding microwave plasma

Figure 7 shows the electron current *i*_{e} vs *V* = *V*_{p} − *V*_{app} [*V*_{p} is the plasma potential (here it is 14 V) and *V*_{app} is the applied voltage], and measurements are performed in the expansion of a microwave argon plasma working at a frequency of 2.45 GHz, an incident power of 100 W, and a pressure of 2 Torr. The Langmuir probe is a cylindrical single probe (tungsten wire of 5 mm length and 0.1 mm diameter) located in the plasma expansion. The ion current collected at saturation is determined at large values of the negative biased voltage. The collected electron current is calculated by removing the ion current from the total current. In this case, collision has no effect on the plasma potential of the probe characteristic because the ratio S ranges from 2.1 to 7.2, considering *λ*_{e} ranges from 2 10^{−4} to 5 10^{−4} m. The energy resolution corresponds to Δ*ε*_{e}/T_{e} = 0.07.

The figure also shows the first and second derivative functions calculated using the SHC method.

Figure 8 shows the change in $\Psi \epsilon e$ versus ε_{e} for different values of *V*_{app}. The values are calculated using Eq. (11). It can be seen that the diffusion parameter decreases with increasing electron energy and converges to 0 at large electron energy values. The low energy electrons are mainly affected by the disturbance due to electron collision in the sheath.

Under these conditions and using the process previously described, the real EEDF corresponding to the Swift law can be calculated. Figure 9 compares the EEDF calculated using the Druyvesteyn equation [Eq. (6)] and the corrected one calculated using [Eq. (10)]. The change between both the EEDFs is mainly observed at low electron energy. Figure 10 compares the experimental *i*_{e}(*V*) probe characteristic obtained first with values calculated using the Swift law [Eq. (7)], taking into account the electron diffusion parameter, and second with values calculated using the Langmuir law [Eq. (2)], which corresponds to $\Psi \epsilon e=0$. The results show that better agreement is obtained when the electron diffusion parameter is taken into account. Nevertheless, this effect is not important in the present case because the diffusion coefficient does not change drastically at low pressure. Under these conditions, the error bar corresponding to the calculated first and second derivative functions of *I*(*V*) can have an important effect on the result.

### B. Microwave plasma confined by a magnetic field

The method has also been tested in the case of magnetized plasma. Figure 11 gives the results obtained in magnetized plasma sustained in hydrogen. The probe characteristic has been measured by Cortazar *et al.*^{34} in a microwave plasma working at 2.45 GHz with a power equal to 1500 W and sustained in hydrogen at 0.19 Pa. The Langmuir probe is a cylindrical single probe (tungsten wire of 6 mm length and 0.5 mm diameter) located in the middle of the reactor along the Z-axis of the cylindrical plasma chamber. It is parallel to the magnetic field B_{z}, whose intensity is equal to 97 mT. A schematic of the experimental setup is shown in Ref. 34.

Figure 11 shows the corrected EEDF calculated taking into account the Swift law; using Eq. (10), it is compared to the EEDF calculated using the Druyvesteyn equation [Eq. (6)]. Figure 12 compares the experimental electron current values with those calculated using the Swift law (7) and using the Langmuir law (2). It can be seen that the best fit with experimental values is obtained after correction using the Swift law. The magnetic field acts on the electron motion and consequently on their radial diffusion coefficient component, which depends on the Larmor radius.^{35} The radial diffusion coefficient is lower when the probe is oriented along the magnetic field than perpendicularly and at low electron energy. Consequently, the electron energy within the sheath around the probe collector is changed because it depends on the applied voltage and also on the electron diffusion coefficient.

### C. Plasma “bubble”

The method has also been tested in a plasma “bubble” using the results obtained by Stenzel and Urrutia.^{32} The *I*_{e}(*V*) probe characteristic is measured in a plasma produced by a filamentary cathode (V_{discharge} = 30 V and I_{discharge} = 100 mA) in argon at 5 × 10^{−4} m Torr. A magnetic field of 2 × 10^{−3}T is applied, and a spherical grid of coarse mesh (0.25 mm) around the plasma bubble is biased at 1 V. A schematic of the experimental setup is shown in Ref. 32. As done previously, Figure 13 we compare the EEDF calculated using Eq. (6) and the one corrected in the case of Eq. (10). Figure 14 compares the experimental *I*_{e}(*V*) probe characteristics with the one calculated first using the Langmuir law and second calculated using the Swift law. The magnetic field and the pressure gas are too low to have a significant effect on the *I*_{e}(*V*) probe characteristic, and because of numerical errors on the values of the first and second derivative due to the signal noise, the best agreement between experiment and theory is obtained using the Langmuir law.

### D. RF plasma confined with a magnetic field

The method has also been tested for measurements performed by means of a Langmuir probe (*r*_{p} = 25 *μ*m and *l*_{p} = 3 mm) by Calderelli *et al.*,^{20} in an RF plasma reactor working at a frequency of 27.12 MHz, at a power of 1 kW, and in argon at pressure ranging from 0.51 × 10^{−3} to 1 × 10^{−3} Torr. In this article, the authors compare different methods to determine the second derivative function of the *I*(*V*) characteristic. We have tested the SHC method on the *I*(*V*) curve used in Ref. 20 and obtained the results at 200 W using a low magnetic field intensity equal to 0.03 T. Calculations have been performed to determine the second derivative function for different values of u_{0} ranging from 0.9 to 4 (*u*_{0} is the adjustable parameter of the SHC method). Results vs (*V*_{p} − *V*_{app}) are compared and shown in Fig. 15. It can be seen that the second derivative is strongly noisy for low values of u_{0}, up to *u*_{0} = 3, and is well defined for larger values. It is defined only on one decade for *u*_{0} = 0.9 and over two decades for *u*_{0} = 4. For *u*_{0} larger than 3, only 1% of the signal intensity cannot be correctly treated. The plasma potential *V*_{p} is equal to 67.0, 66.5, and 65.4 V for *u*_{0} = 0.9, 3, and 4, respectively. Integrating the second derivative function vs *V* and comparing it to the experimental *I*(*V*) curve, it can be seen that the best fit is obtained for *u*_{0} = 4 (see Fig. 16).

Calderelli *et al.*^{20} tested the same method (named the ac superimposed method in Ref. 20) using *u*_{0} < 1. The results shown in Fig. 6(b) in Ref. 20 are strongly noisy because of the low value of *u*_{0}. These are similar to ours obtained for *u*_{0} = 0.9. It is worth noting that the shoulder observed at (*V*_{p} − *V*_{app}) ranging from 20 to 25 eV, which is observed in Ref. 20 using the Savitzky–Golay method, the Gaussian filter, the b spline, and the Blackman window methods and in the present work (see Fig. 15 for *u*_{0} = 4), disappeared in Ref. 20 when using the “analog” method, probably because of a smoothing effect.

In Ref. 20, the authors measured a plasma potential value of *V*_{p} = 66.95 V using the same method, *V*_{p} = 67, 14 V using the Savitzky–Golay method, and *V*_{p} = 68.08 V using the “analog” method.

Figure 17 compares the experimental values of the *I*_{e}(*V*) probe characteristic first with the values calculated using the Langmuir law and second using the Swift law after correction of the EEDF using Eq. (10). It can be seen that in the present case, the two laws give similar results. The change due to the correction of the EEDF taking into account the effect of electron diffusion parameter is not important. This is due to the low pressure and low intensity value of the magnetic field (0.03 T) used to confine the plasma and due to the error bars on the electron current value and consequently on the first and second derivative functions.

## VI. CONCLUSION

This article reports on numerical methods used to determine the derivative functions and the EEDF from the *I*(*V*) probe characteristics measured using a cylindrical Langmuir probe. The effect of electron diffusion within the sheath formed around the probe is taken into account. The electron diffusion coefficient depends on the plasma pressure (e-neutral collision) and also on the magnetic field, which can be used to confine the plasma (electron gyration in the field). In these cases, the electron current collected by the probe is not representative of the EEDF within the plasma bulk, and some corrections are necessary.

After a theoretical review of the electron current collected by a negatively biased probe (retardation conditions), numerical methods are used to determine the second derivative function of the *I*(*V*) probe characteristic and to correct the EEDF, taking into account the electron diffusion within the probe sheath. These methods are applied to various cases reported in the literature or to our own experiments carried out at the laboratory and corresponding to different plasma conditions (microwave, radio frequency, and direct-current plasmas with or without confining magnetic field). Results show that when a low magnetic field intensity (2 or 30 mT) is applied along an axis parallel to the probe, no greater change is observed on the *I*(*V*) probe characteristic. The change is of the order of the measurement errors due to the signal noise. However, for a larger magnetic field value (97 mT), a significant change is obtained compared to the *I*(*V*) probe measured without magnetic field, so the Druyvesteyn and Langmuir equations cannot be used because of the effect of the low electron diffusion through the sheath (Swift law). Results show that the noise observed on the curve when *I*_{e} is calculated using the Swift law depends strongly on the error due to the calculation of the ratio $die/dVpd2ie/dVp2$.

Regardless of the theory used to measure the EEDF (with or without correction due to the electron diffusion parameter), the electrostatic probes can be used to study the plasma only if the disturbed volume by the probe has a characteristic length (*r*_{d}) much smaller than the electron energy relaxation length λ_{ε}, which is the characteristic length of the undisturbed EEDF formation and depends on the electron collision cross sections with any other species present in the plasma.^{6,11,36} Hence, as long as *r*_{d} ≪ λ_{ε}*,* the EEDF measured by the probe corresponds to the undisturbed EEDF within the plasma. For larger *r*_{d} values, the effect of collision processes within the sheath drastically changes the EEDF. Hence, measurements are not representative of the EEDF within the plasma bulk.

The probe disturbed length *r*_{d} is of the order of the probe sheath radius (probe radius + sheath thickness) only in the case of a collisionless regime in the sheath. This condition is generally not fulfilled, especially because of ions of mean free path being much lower than the electron mean free path. Hence, for a cylindrical probe of length *L* and radius *r*_{p}, we use $rd\u2248rPlnLrP$.

## ACKNOWLEDGMENTS

The authors would like to express their deep gratitude to O.D. Cortázar, A. Megía-Macías, R.L. Stenzel, J.M. Urrutia, and F. Filleul for their contributions in providing the I(V) data files recorded in different plasma reactors and for information and discussions concerning these experiments. The accessibility to these data has granted the opportunity to test these numerical methods under different plasma conditions.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**J. L. Jauberteau**: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Software (equal); Supervision (equal); Validation (equal); Writing – original draft (equal). **I. Jauberteau**: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Software (equal); Supervision (equal); Validation (equal); Writing – original draft (equal).

## DATA AVAILABILITY

Data are available if they are requested.

## REFERENCES

*Numerical Recipes: The Art of Scientific Computing*

*Plasma Chemistry*