- Email: [email protected]
Improved fictitious charge method for calculations of electric potential and field generated by point-to-plane electrodes
Journal of Electrostatics 76 (2015) 24e30
Contents lists available at ScienceDirect
Journal of Electrostatics journal homepage: www.elsevier.com/locate/elstat
Improved fictitious charge method for calculations of electric potential and field generated by point-to-plane electrodes Djilali Benyoucef a, *, Mohammed Yousfi b a b
Laboratoire G enie Electrique et Energie Renouvelables, Chlef University, Algeria Universit e de Toulouse, LAPLACE, UMR CNRS 5213, 118 Route de Narbonne, Toulouse 31062, France
a r t i c l e i n f o
a b s t r a c t
Article history: Received 30 November 2014 Received in revised form 8 April 2015 Accepted 5 May 2015 Available online 16 May 2015
In the framework of standard tip-to-plane electrode geometry favorable to corona streamer discharge development at atmospheric pressure, this work is devoted to the improvement of fictitious charge method for calculations of electric potential and field repartition when the tip is powered by a DC voltage. It is in fact dedicated to implement the image charge method (generally used in plane-to-plane electrodes) in the case of a point-to-plane geometry. The numerical method is based on the solution an open system of n equations with m unknowns (n >> m) where m is the number of fictitious charges and n the number of contours at the surface of the tip electrode defining the boundary conditions. This numerical technique can accurately interpolate the shape of the electrode tip whatever its geometry and hence allows us to accurately calculate the electric potential and field even at a position very close to the electrode. It is noteworthy that the solution of such open system of equations cannot be obtained from conventional techniques (Cramer, Gauss, matrix inversion, etc.). We used the method of least squares which enables us to close the equation systems and to find the optimal solution fulfilling all the required boundary conditions. The present method is therefore based on the coupling between the conventional method of fictitious charges using image charge method and the optimization by the Least Squares Method. The results of simulation show that the punctual fictitious charges have given the most accurate results when the electrode has symmetry of revolution like the present geometry of a pen shape anode cylinder ended by a sharp tip set in front of cathode plane. © 2015 Elsevier B.V. All rights reserved.
Keywords: Electrostatic potential and fields Least square method Fictitious charges Image charges Electric discharges
Introduction The calculation of the electric field created by incurved electrodes and/or tip with complex geometry has a great interest in the cases of for instance the electrical networks for sparking point [1], the electric discharge (corona discharges [2e4]) which can be used in multiple applications as waste gas treatment [5,6,7], or ozone production [8], the bio-decontamination and sterilization [9], the technology of Micro-Electro-Mechanical Systems (MEMS) and also in the sensors of the Atomic Force Microscope (AFM) to characterize insulating materials or semiconductors with a resolution at the atomic scale [10]. The algorithm of fictitious charges is one of the fastest algorithms for calculating the electric field as already successfully used
* Corresponding author. E-mail addresses: (D. Benyoucef).
[email protected],
http://dx.doi.org/10.1016/j.elstat.2015.05.005 0304-3886/© 2015 Elsevier B.V. All rights reserved.
[email protected]
in the case of plane-to-plane electrode geometry (see e.g. Refs [11,12]). It consists to replace the electrodes by separate fictitious charges within these electrodes. Their values and their locations will be determined so that their collective effect must fulfill the boundary conditions of the known potential on the electrode surface. Generally, these fictitious charges can be considered as punctual charges uniformly distributed on finite or infinite lines in the case of Cartesian electrode configuration. It can also be considered as uniform charges on circles in the case where the geometry has revolution symmetry. The classical method of fictitious charges is to consider n charges (punctual, linear or circular) with n boundary conditions of potential on the electrode surface [13]. Each boundary condition represents a linear equation of n unknown where each unknown represents one fictitious charge. This leads to a system of n equations with n unknown and the problem can be solved easily (by using the method of Gauss, method of Cramer, matrix inversion, etc.). Because of the discrete nature of the fictitious charges, this technique requires the use of a large number of charges to achieve
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
an acceptable accuracy. Unfortunately, due to some limitation of the charge number, the solution becomes unstable due to the small distance between these charges and cannot reach the required accuracy in the vicinity of the tip which is considered in the present work. This method developed in the present framework is devoted to the improvement of this technique based on the finding of an optimal solution of a system of n equations with m unknown (n >> m) where m is the number of fictitious charges and n the number of boundary conditions on the incurved electrode. The optimization of the solution of such a system is based on the method of Least Squares. It provides the optimal values of fictitious charges which meet the imposed conditions of the potential value on the electrode surface for a very large number of points on the surface of the electrode. Section Description of the method following this introduction is devoted to the description of the method of calculation. Then, results based on parametric study are shown analyzed in Section Simulation results terminated by a short conclusion.
z-axis of revolution symmetry system at distances respectively (þzj) and (zj) with respect to the upper plane is given by equation (1), where ε0 is the vacuum permittivity:
2 up;ij ¼
3
qj 6 1 1 7 4qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 5 ¼ aij qj 4pε0 2 2 r i þ zj zi r i þ zj þ zi
We considered the case of tip-plane electrode geometry with revolution symmetry around the z-axis. The tip is powered by direct current (DC) voltage and the grounded plane is the electrostatic virtual mirror for the application of the image charge method. Fig. 1 displays the chosen electrode configuration with the locations of fictitious charges and their images which are either punctual or circular fictitious charges. The red circles show the shape of the boundary conditions that are considered in this work. The electric potential up,ij created at position Mp(zi,ri) by a punctual fictitious (qj) charge and its image (qj) positioned on the
(1)
In the case of a charge qj uniformly distributed over a circle of radius Rj, the electric potential uc,ij created at the position Mc(zi,ri) by this circular charge distribution and its image is given by the following integral over the polar angle q:
2 uc;ij ¼
qj 4p2 ε0
6 4
Zp 0
dq ffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 zi zj þ Rj þ r 2i 2Rj ri cos q Zp
0
Description of the method
25
3
(2)
dq 7 ffi5 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 2 zi þ zj þ Rj þ r i 2Rj ri cos q
By substituting the variableq byq ¼ p24 and after a mathematical arrangement, we obtain [14]:
2
uc;ij
K kfic 6 ¼ 4qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 4p2 ε0 zi zj þ Rj þ ri qj
3
(3)
Kðkim Þ 7 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 5 ¼ aij qj zi þ zj þ Rj þ ri where K(k) is the complete elliptic integral of the first kind and kfic rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Rj ri and kim are defined as: kfic ¼ 2 and ðzi zj Þ2 þðRj þri Þ2 rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Rj ri kim ¼ 2 2 2 ðzi þzj Þ þðRj þri Þ
Zp=2 with KðkÞ ¼ 0
d4 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ RF 0; 1 k2 ; 1 1 k2 sin2 4
(4)
RF is the symmetric integral of first kind
RF ðx; y; zÞ ¼
1 2
Z∞ 0
dt pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðt þ xÞðt þ yÞðt þ zÞ
(5)
The reader can find the algorithm of this kind of integrals in Ref. [15]. The m fictitious charges qi (with j ¼ 1, m) and their images must verify the boundary conditions in n circles on the surface of the powered electrode. This therefore leads to the following system of equations:
8 m P > > > a1j qj ¼ V0 > > > j¼1 > > >P m < a2j qj ¼ V0 > > j¼1//// > > > m > >P > anj qj ¼ V0 > :
(6)
j¼1
Fig. 1. Configuration of tip-to-plane electrodes with a gap distance D and two examples of image charges of punctual and circular charges (see grey point and circle in the bottom side representing two images symmetric with respect to the grounded plane while a/2 being the tip angle).
where n >> m and V0 is the potential of the powered tip electrode and aji are the potential coefficients. The solution of this system of equations cannot be obtained from a conventional method because
26
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
the number n of equations is greater than the number of m unknowns. Therefore, we have to find an optimal solution that can minimizes the sum of square errors J that is defined as follows:
2
3 3 2 3 2 ε1 V0 a11 … a1m 6 ε2 7 7 6 7 6 6 7 ¼ 6 V0 7 6 a21 … a2m 7 4…5 4…5 4 … … … 5 εn V0 an1 … anm |fflfflffl{zfflfflffl} |fflfflffl{zfflfflffl} |fflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} ε
U
A
J ¼ εT ε ¼ ðU A$Q ÞT ðU A$Q Þ
2
3 q1 6 q2 7 6 7 4…5 qm |fflfflfflffl{zfflfflfflffl}
(7)
Q
(8)
where ε is the vector of errors between the potential of the tip and the calculated potential by using the fictitious charges at the contours of the boundary conditions and εT is the transpose of ε. The b which minimizes this error can be vector of fictitious charges Q calculated as follows:
vJ ¼0 vQ Q ¼ Q b
(9)
After a mathematical transformation, the optimal values of the fictitious charges can be written as:
i1 h b ¼ AT $A AT $U Q
(10)
For optimizing the positions of fictitious charges, we can use one of the methods of optimization either the Genetic Algorithm [13] or the Simplex method [16].
Simulation results The first simulation results are shown an example with interelectrode distance D ¼ 5 mm, L ¼ 45 mm, a¼p/3 and tip curvature radius r ¼ 3 mm in the case of an applied DC voltage V0 ¼ 5 kV. These are followed by another results using with different tip curvature radius and also negative polarity for applied voltage V0. Fig. 2 displays the relative error between the calculated tip potential and the fixed potential at the electrode surface versus the distance z at the plane and the ratio h between the ring radius of fictitious charges and the radius of the contours of boundary conditions (radius of the electrode). h ¼ 0 means that the ring of fictitious charge becomes a punctual fictitious charge. Fig. 2 shows near the tip (z ¼ 5 mm) the clear decrease of the relative error when h decreases. This means that the optimal distribution of the fictitious charges which gives the most accurate results when the
Fig. 2. Relative error on potential for different values of h for a gap distance D ¼ 5 mm.
electrode has revolution symmetry corresponds to the case of the punctual charge Fig. 3 displays the relative error on the potential between the calculated and fixed values at the electrode tip. As this relative error is less than 0.003, our method allows us to calculate with a very good accuracy the maximum electric field and hence the Schwaiger factor, the capacity or the electrostatic force, etc. Furthermore, the following results are calculated using 50 punctual fictitious charges with 1500 contours of boundary conditions uniformly distributed according to the z-axis. Fig. 4 displays the spatial variation of electric potential. As expected in the case of a tip electrode, the electric potential rapidly falls down and reaches about two fifth of applied voltage at a distance of 1 mm far from the tip. Obviously, along the z-axis, this variation depends of the curvature radius. Fig. 5 displays the absolute spatial variation of the axial components of the electric field modulus while Fig. 6 shows the radial one and Fig. 7 the variation of the total electric field modulus. As expected, the maximum of the axial component is positioned close to the tip on the z-axis (i.e. when radial position r ¼ 0). For the present configuration, its absolute value is equal to 21.54 kV/cm and Schwaiger factor to 0.46. The maximum of the radial component (13.28 kV/cm) is located at the vicinity of the powered electrode at the end of curvature. The maximum and the position of the total electric field are those of its axial component. The total electric field and the potential in the gap along the z-axis are shown in Fig. 8. Fig. 9 displays the location and the variation of the punctual fictitious charges. Due to the image charge method, the charge magnitudes are alternately distributed (i.e. with negative and positive values) with a decreasing magnitude for higher z values close to the plane. Fig. 10 displays the variation of the surface charge density son the plane surface along r-axis. It was calculated from the classical equation:
E¼
s n ε0
(11)
E is the electric field vector perpendicular to the plane electrode. The surface charge density on the plane is proportional to the axial component of the electric field calculated close to the plane. The negative sign of the charge indicates that the electric field is directed from the tip to the plane. Capacity C of the electrode configuration versus the inter-electrode distance is shown in Fig. 11.It corresponds to the ratio between the absolute value of the accumulated charges on the plane electrode overV0the applied DC voltage. C obviously decreases with the increasing of inter-
Fig. 3. Relative error of potential in the case of punctual fictitious charges for a gap distance D ¼ 5 mm.
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
27
Fig. 4. Two representations of the spatial variation of electric potential for a gap distance D ¼ 5 mm and a tip voltage V0 ¼ 5 kV: (a) mesh, (b) maps.
Fig. 5. Two representations of the spatial variation of axial electric field modulus for a gap distance D ¼ 5 mm and a tip voltage V0 ¼ 5 kV: (a) mesh, (b) maps.
Fig. 6. Two representations of the spatial variation of radial electric field modulus for a gap distance D ¼ 5 mm and a tip voltage V0 ¼ 5 kV: (a) mesh, (b) maps.
electrode distance [17]. This is due to the decrease of the electric field in the vicinity of the plane and therefore the concomitant decrease of the surface charge density on the plane. In this work, the capacity is calculated from the classical equation:
Z C¼
0
R
sðrÞdr V0
(12)
Fig. 12 shows the absolute value of the total electric field along the surface of the tip for different radii of curvature (r). The electric field maxima Emax, as expected, decrease when r increases. These maxima are Emax ¼ 33.80 kV cm1, Emax ¼ 27.80 kV cm1,
Emax ¼ 21.54 kV.cm-1and Emax ¼ 18.92 kV cm1 respectively for r ¼ 1 mm, r ¼ 2 mm, r ¼ 3mmand r ¼ 4 mm. This decrease of the electric field is obviously due to the expansion of the equipotential lines near the tip which is due to the increase of the radius tip. In the same time, this causes the increase of electric field near the opposite plane electrode. This in turn increases the plane surface charge and therefore the increase of the capacitance of the tip-plane gap. This capacitance is respectively equal to C¼0.85 pF, C¼0.90 pF, C¼0.93 pF and C¼0.95 pF in the case of the considered radius (r ¼ 1 mm, r ¼ 2 mm, r ¼ 3 mm and r ¼ 4 mm). Last, for sake of a validation, we show a comparison between the results of the present improved method and some analytical results of the electric potential and field obtained in a specific case versus
28
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
Fig. 7. Two representations of the spatial variation of the total electric field modulus for a gap distance D ¼ 5 mm and a tip voltage V0 ¼ 5 kV: (a) mesh, (b) maps.
Fig. 8. Spatial variation of electric potential (green line) and electric field (blue dashed line) in gap along z-axis for DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 3 mm.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
z-axis. The results of improved fictitious charge method are obtained for a powered electrode when the tip is terminated by a half sphere of radius r ¼ 4 mm and for a gap distance D between the tip of the semi-sphere and the plane. The analytical expression of the electric potential U(z)generated along the z-axis in this semispherical/plane configuration can be written as follows:
Fig. 9. Locations and magnitudes of punctual fictitious charges on z-axis in the powered electrode axis for DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 3 mm.
Fig. 10. Variation of the surface charge density on the plane for DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 3 mm.
UðzÞz
V0 ð2D þ rÞ 1 1 2D jD þ r zj jD þ r þ zj
(13)
Figs. 13 (a) and 14 (a) display the comparison between the analytical and the present numerical solution using fictitious image charges of the electric potential for positive and negative DC voltage while Figs. 13 (b) and 14 (b) display the same comparison in the case of the electric field along the z-axis (i.e. for r ¼ 0). These
Fig. 11. Variation of the configuration capacity as a function of gap distance for curvature radius r ¼ 3 mm.
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
29
any change on the relative error between analytical and numerical results. Last, if the tip radius is chosen much smaller than the interelectrode distance (<1 mm) for conditions close to tip radius currently used in corona discharge setups at atmospheric pressure, the present image charge numerical method is very easily adaptable to such extreme operating parameters. Fig. 16a and b displays an illustration of the sharp electric field variation at the close vicinity of the tip. The magnitude of electric field reaches for instance 165 kV/cm for r ¼ 0.1 mm against 33.8 kV/cm for r ¼ 1 mm. Conclusion
Fig. 12. Variation of the total Electric field modulus along the surface of the tip axis for DC voltage V0 ¼ 5 kV and a gap distance D ¼ 5 mm.
comparisons clearly show a very good agreement in both cases (electric field and potential) since the relative deviation between analytical and numerical results is less than 0.2% for electric potential and less than 1.5% for electric field. The latter is a little bit higher because it accumulates an additional numerical error due to the potential derivative required for the electric field. Therefore, these comparisons can be considered as a reliable validation test for the fictitious image charge method developed in the present work When the applied DC voltage is increased from for instance 5 kV (Fig. 13) to 10 kV (Fig. 15), only the magnitude of potential and electric field are increased as expected by a factor twice without
This paper was devoted to the implementation the image charge method in the case of tip-to-plane electrode configuration and to the improvement and the validation of this interesting method in the case of revolution symmetry geometry. It is shown that this method simulates with a satisfactory manner the electric potential and field. The major advantage of this technique is that it requires no boundary conditions for solving the Laplace equation. The calculation time is very short in the configuration studied in the present work (less than 1 s in a standard laptop with core Intel I7 process). Furthermore, this image charge method allows us very easily to introduce the effect of space charges, very useful for simulation of the streamer corona discharges, provided their spatial distribution is known. This image method is validated from comparison with analytical results known in the case of a semi-sphereto-plane configuration. Finally, it was shown that the punctual fictitious charges give, in comparison to for instance planar or circular distribution, the most accurate results when the electrode has symmetry of revolution.
Fig. 13. Comparison between the results of the improved fictitious charge method and the analytical results in the gap along the z-axis (a) electric potential, (b) electric field for positive DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 4 mm.
Fig. 14. Comparison between the results of the improved fictitious charge method and the analytical results in the gap along the z-axis (a) electric potential, (b) electric field for negative DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 4 mm.
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
V (kV)
8
-0.2
-0.4
4
-0.5
0 0
Analytic method Improved method Error on electric field (%)
-0.3
6
2
b
-10
-20
2
3
4
-0.7 5
z (mm)
1.5 1 0.5
-30
0
-0.6
1
2
error on electric field(%)
Analytic method Improved method Error on potential (%)
E (kV/cm)
a
10
error on potential(%)
30
-0.5 -40 0
1
2
3
4
-1 5
z (mm)
Fig. 15. Comparison between the results of the improved fictitious charge method and the analytical results in the gap along the z-axis (a) electric potential, (b) electric field for positive DC voltage V0 ¼ 10 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 4 mm.
Fig. 16. Variation of the total Electric field modulus in the case of 3 small tip radius for DC voltage V0 ¼ 5 kV and gap distance D ¼ 5 mm: (a) along the z-axis (b) along the surface of the tip axis.
References [1] A. Younes, M. Younes, H. Sayah, M. Bilici, A. Samuila, L. Dascalescu, Effect of spark discharges on the trajectories of insulating particles in roll-type coronaelectrostatic separators. Experimental and numerical study, J. Electrost. 71 (1) (2013) 84e91. [2] Rakshit Tirumala, Go David B, Comparative study of corona discharge simulation techniques for electrode configurations inducing non-uniform electric fields, J. Electrost. 72 (2) (2014) 99e106. [3] X. Wang, N. Mingjiang, X. Gang, Z. Jiapeng, G. Xiang, C. Kefa, An analytical method for DC negative corona discharge in a wire-cylinder device at high temperatures, J. Electrost. 72 (4) (2014) 270e284. [4] M. Kachl, L. Dascalescu, Corona discharges in asymmetric electrode configurations, J. Electrost. 72 (1) (2014) 6e12. [5] B.M. Penetrante, S.E. Schultheis, “Non-thermal Plasma Techniques for Pollution Control,” Part a and B, Editors, Springer-Verlag, New York, 1993. [6] Y.S. Mok, S.W. Ham, I.S. Nam, Evaluation of energy-utilization efficiencies for SO2 and NO removal by pulsed corona discharge process, Plasma Chem. Plasma Process. 18 (4) (1998) 535e550. rimentales des de charges couronne pour la [7] A. Abahazem, Etudes expe pollution des gaz, The se de doctorat de l’, Universite Paul Sabatier, 2009. de [8] J. Chen, P. Wang, Effect of relative humidity on electron distribution and ozone production by DC coronas in air, Plasma Sci. IEEE Trans. 33 (2) (2005) 808e812.
[9] I.V. Timoshkin, M. Maclean, M.P. Wilson, M.J. Given, S.J. Mac Gregor, Tao Wang, J.G. Anderson, Bactericidal effect of corona discharges in atmospheric air, Plasma Sci. IEEE Trans. 2012 40 (10) (2012) 2322e2333. [10] A. Boularas, F. Baudoin, C. Villeneuve-Faure, S. Clain, G. Teyssedre, Multidimensional modelling of electrostatic force distance curve over dielectric surface: influence of tip geometry and correlation with experiment, J. Appl. Phys. 116 (8) (2014) 084106. [11] A.J. Davies, C.S. Davies, C.J. Evans, Computer simulation of rapidly developing gaseous discharges, Proc. IEE 118 (1971) 816e823. [12] M. Yousfi, A. Poinsignon, A. Hamani, Finite element method for conservation equations in electric discharge areas, J. Comp. Phys. 113 (2) (1994) 268e278. [13] M. Talaat, Electrostatic field calculation in air gaps with a transverse layer of dielectric barrier, J. Electrost. 72 (5) (2014) 422e427. [14] D. Benyoucef, M. Yousfi, effects of increasing magnetic field and decreasing pressure on asymmetric magnetron radio frequency discharges, Plasma Sci. IEEE Trans. 41 (4) (2013) 829e838. [15] B.C. Carlson, Numerical computation of real and complex elliptic integrals, Numer. Algorithms 10 (1995) 13e26. [16] M. Yousfi, M.D. Benabdessadok, Boltzmann equation analysis of electronmolecule collision cross sections in water vapor and ammonia, J. Appl. Phys. 80 (12) (1996) 6619e6631. [17] D. Raouti, S. Flazi, D. Benyoucef, electrical modelling of a positive point to plane corona discharge at atmospheric pressure, Contributions Plasma Phys. 54 (10) (2014) 851e858.
Contents lists available at ScienceDirect
Journal of Electrostatics journal homepage: www.elsevier.com/locate/elstat
Improved fictitious charge method for calculations of electric potential and field generated by point-to-plane electrodes Djilali Benyoucef a, *, Mohammed Yousfi b a b
Laboratoire G enie Electrique et Energie Renouvelables, Chlef University, Algeria Universit e de Toulouse, LAPLACE, UMR CNRS 5213, 118 Route de Narbonne, Toulouse 31062, France
a r t i c l e i n f o
a b s t r a c t
Article history: Received 30 November 2014 Received in revised form 8 April 2015 Accepted 5 May 2015 Available online 16 May 2015
In the framework of standard tip-to-plane electrode geometry favorable to corona streamer discharge development at atmospheric pressure, this work is devoted to the improvement of fictitious charge method for calculations of electric potential and field repartition when the tip is powered by a DC voltage. It is in fact dedicated to implement the image charge method (generally used in plane-to-plane electrodes) in the case of a point-to-plane geometry. The numerical method is based on the solution an open system of n equations with m unknowns (n >> m) where m is the number of fictitious charges and n the number of contours at the surface of the tip electrode defining the boundary conditions. This numerical technique can accurately interpolate the shape of the electrode tip whatever its geometry and hence allows us to accurately calculate the electric potential and field even at a position very close to the electrode. It is noteworthy that the solution of such open system of equations cannot be obtained from conventional techniques (Cramer, Gauss, matrix inversion, etc.). We used the method of least squares which enables us to close the equation systems and to find the optimal solution fulfilling all the required boundary conditions. The present method is therefore based on the coupling between the conventional method of fictitious charges using image charge method and the optimization by the Least Squares Method. The results of simulation show that the punctual fictitious charges have given the most accurate results when the electrode has symmetry of revolution like the present geometry of a pen shape anode cylinder ended by a sharp tip set in front of cathode plane. © 2015 Elsevier B.V. All rights reserved.
Keywords: Electrostatic potential and fields Least square method Fictitious charges Image charges Electric discharges
Introduction The calculation of the electric field created by incurved electrodes and/or tip with complex geometry has a great interest in the cases of for instance the electrical networks for sparking point [1], the electric discharge (corona discharges [2e4]) which can be used in multiple applications as waste gas treatment [5,6,7], or ozone production [8], the bio-decontamination and sterilization [9], the technology of Micro-Electro-Mechanical Systems (MEMS) and also in the sensors of the Atomic Force Microscope (AFM) to characterize insulating materials or semiconductors with a resolution at the atomic scale [10]. The algorithm of fictitious charges is one of the fastest algorithms for calculating the electric field as already successfully used
* Corresponding author. E-mail addresses: (D. Benyoucef).
[email protected],
http://dx.doi.org/10.1016/j.elstat.2015.05.005 0304-3886/© 2015 Elsevier B.V. All rights reserved.
[email protected]
in the case of plane-to-plane electrode geometry (see e.g. Refs [11,12]). It consists to replace the electrodes by separate fictitious charges within these electrodes. Their values and their locations will be determined so that their collective effect must fulfill the boundary conditions of the known potential on the electrode surface. Generally, these fictitious charges can be considered as punctual charges uniformly distributed on finite or infinite lines in the case of Cartesian electrode configuration. It can also be considered as uniform charges on circles in the case where the geometry has revolution symmetry. The classical method of fictitious charges is to consider n charges (punctual, linear or circular) with n boundary conditions of potential on the electrode surface [13]. Each boundary condition represents a linear equation of n unknown where each unknown represents one fictitious charge. This leads to a system of n equations with n unknown and the problem can be solved easily (by using the method of Gauss, method of Cramer, matrix inversion, etc.). Because of the discrete nature of the fictitious charges, this technique requires the use of a large number of charges to achieve
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
an acceptable accuracy. Unfortunately, due to some limitation of the charge number, the solution becomes unstable due to the small distance between these charges and cannot reach the required accuracy in the vicinity of the tip which is considered in the present work. This method developed in the present framework is devoted to the improvement of this technique based on the finding of an optimal solution of a system of n equations with m unknown (n >> m) where m is the number of fictitious charges and n the number of boundary conditions on the incurved electrode. The optimization of the solution of such a system is based on the method of Least Squares. It provides the optimal values of fictitious charges which meet the imposed conditions of the potential value on the electrode surface for a very large number of points on the surface of the electrode. Section Description of the method following this introduction is devoted to the description of the method of calculation. Then, results based on parametric study are shown analyzed in Section Simulation results terminated by a short conclusion.
z-axis of revolution symmetry system at distances respectively (þzj) and (zj) with respect to the upper plane is given by equation (1), where ε0 is the vacuum permittivity:
2 up;ij ¼
3
qj 6 1 1 7 4qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 5 ¼ aij qj 4pε0 2 2 r i þ zj zi r i þ zj þ zi
We considered the case of tip-plane electrode geometry with revolution symmetry around the z-axis. The tip is powered by direct current (DC) voltage and the grounded plane is the electrostatic virtual mirror for the application of the image charge method. Fig. 1 displays the chosen electrode configuration with the locations of fictitious charges and their images which are either punctual or circular fictitious charges. The red circles show the shape of the boundary conditions that are considered in this work. The electric potential up,ij created at position Mp(zi,ri) by a punctual fictitious (qj) charge and its image (qj) positioned on the
(1)
In the case of a charge qj uniformly distributed over a circle of radius Rj, the electric potential uc,ij created at the position Mc(zi,ri) by this circular charge distribution and its image is given by the following integral over the polar angle q:
2 uc;ij ¼
qj 4p2 ε0
6 4
Zp 0
dq ffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 zi zj þ Rj þ r 2i 2Rj ri cos q Zp
0
Description of the method
25
3
(2)
dq 7 ffi5 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 2 zi þ zj þ Rj þ r i 2Rj ri cos q
By substituting the variableq byq ¼ p24 and after a mathematical arrangement, we obtain [14]:
2
uc;ij
K kfic 6 ¼ 4qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 4p2 ε0 zi zj þ Rj þ ri qj
3
(3)
Kðkim Þ 7 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 5 ¼ aij qj zi þ zj þ Rj þ ri where K(k) is the complete elliptic integral of the first kind and kfic rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Rj ri and kim are defined as: kfic ¼ 2 and ðzi zj Þ2 þðRj þri Þ2 rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Rj ri kim ¼ 2 2 2 ðzi þzj Þ þðRj þri Þ
Zp=2 with KðkÞ ¼ 0
d4 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ RF 0; 1 k2 ; 1 1 k2 sin2 4
(4)
RF is the symmetric integral of first kind
RF ðx; y; zÞ ¼
1 2
Z∞ 0
dt pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðt þ xÞðt þ yÞðt þ zÞ
(5)
The reader can find the algorithm of this kind of integrals in Ref. [15]. The m fictitious charges qi (with j ¼ 1, m) and their images must verify the boundary conditions in n circles on the surface of the powered electrode. This therefore leads to the following system of equations:
8 m P > > > a1j qj ¼ V0 > > > j¼1 > > >P m < a2j qj ¼ V0 > > j¼1//// > > > m > >P > anj qj ¼ V0 > :
(6)
j¼1
Fig. 1. Configuration of tip-to-plane electrodes with a gap distance D and two examples of image charges of punctual and circular charges (see grey point and circle in the bottom side representing two images symmetric with respect to the grounded plane while a/2 being the tip angle).
where n >> m and V0 is the potential of the powered tip electrode and aji are the potential coefficients. The solution of this system of equations cannot be obtained from a conventional method because
26
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
the number n of equations is greater than the number of m unknowns. Therefore, we have to find an optimal solution that can minimizes the sum of square errors J that is defined as follows:
2
3 3 2 3 2 ε1 V0 a11 … a1m 6 ε2 7 7 6 7 6 6 7 ¼ 6 V0 7 6 a21 … a2m 7 4…5 4…5 4 … … … 5 εn V0 an1 … anm |fflfflffl{zfflfflffl} |fflfflffl{zfflfflffl} |fflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} ε
U
A
J ¼ εT ε ¼ ðU A$Q ÞT ðU A$Q Þ
2
3 q1 6 q2 7 6 7 4…5 qm |fflfflfflffl{zfflfflfflffl}
(7)
Q
(8)
where ε is the vector of errors between the potential of the tip and the calculated potential by using the fictitious charges at the contours of the boundary conditions and εT is the transpose of ε. The b which minimizes this error can be vector of fictitious charges Q calculated as follows:
vJ ¼0 vQ Q ¼ Q b
(9)
After a mathematical transformation, the optimal values of the fictitious charges can be written as:
i1 h b ¼ AT $A AT $U Q
(10)
For optimizing the positions of fictitious charges, we can use one of the methods of optimization either the Genetic Algorithm [13] or the Simplex method [16].
Simulation results The first simulation results are shown an example with interelectrode distance D ¼ 5 mm, L ¼ 45 mm, a¼p/3 and tip curvature radius r ¼ 3 mm in the case of an applied DC voltage V0 ¼ 5 kV. These are followed by another results using with different tip curvature radius and also negative polarity for applied voltage V0. Fig. 2 displays the relative error between the calculated tip potential and the fixed potential at the electrode surface versus the distance z at the plane and the ratio h between the ring radius of fictitious charges and the radius of the contours of boundary conditions (radius of the electrode). h ¼ 0 means that the ring of fictitious charge becomes a punctual fictitious charge. Fig. 2 shows near the tip (z ¼ 5 mm) the clear decrease of the relative error when h decreases. This means that the optimal distribution of the fictitious charges which gives the most accurate results when the
Fig. 2. Relative error on potential for different values of h for a gap distance D ¼ 5 mm.
electrode has revolution symmetry corresponds to the case of the punctual charge Fig. 3 displays the relative error on the potential between the calculated and fixed values at the electrode tip. As this relative error is less than 0.003, our method allows us to calculate with a very good accuracy the maximum electric field and hence the Schwaiger factor, the capacity or the electrostatic force, etc. Furthermore, the following results are calculated using 50 punctual fictitious charges with 1500 contours of boundary conditions uniformly distributed according to the z-axis. Fig. 4 displays the spatial variation of electric potential. As expected in the case of a tip electrode, the electric potential rapidly falls down and reaches about two fifth of applied voltage at a distance of 1 mm far from the tip. Obviously, along the z-axis, this variation depends of the curvature radius. Fig. 5 displays the absolute spatial variation of the axial components of the electric field modulus while Fig. 6 shows the radial one and Fig. 7 the variation of the total electric field modulus. As expected, the maximum of the axial component is positioned close to the tip on the z-axis (i.e. when radial position r ¼ 0). For the present configuration, its absolute value is equal to 21.54 kV/cm and Schwaiger factor to 0.46. The maximum of the radial component (13.28 kV/cm) is located at the vicinity of the powered electrode at the end of curvature. The maximum and the position of the total electric field are those of its axial component. The total electric field and the potential in the gap along the z-axis are shown in Fig. 8. Fig. 9 displays the location and the variation of the punctual fictitious charges. Due to the image charge method, the charge magnitudes are alternately distributed (i.e. with negative and positive values) with a decreasing magnitude for higher z values close to the plane. Fig. 10 displays the variation of the surface charge density son the plane surface along r-axis. It was calculated from the classical equation:
E¼
s n ε0
(11)
E is the electric field vector perpendicular to the plane electrode. The surface charge density on the plane is proportional to the axial component of the electric field calculated close to the plane. The negative sign of the charge indicates that the electric field is directed from the tip to the plane. Capacity C of the electrode configuration versus the inter-electrode distance is shown in Fig. 11.It corresponds to the ratio between the absolute value of the accumulated charges on the plane electrode overV0the applied DC voltage. C obviously decreases with the increasing of inter-
Fig. 3. Relative error of potential in the case of punctual fictitious charges for a gap distance D ¼ 5 mm.
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
27
Fig. 4. Two representations of the spatial variation of electric potential for a gap distance D ¼ 5 mm and a tip voltage V0 ¼ 5 kV: (a) mesh, (b) maps.
Fig. 5. Two representations of the spatial variation of axial electric field modulus for a gap distance D ¼ 5 mm and a tip voltage V0 ¼ 5 kV: (a) mesh, (b) maps.
Fig. 6. Two representations of the spatial variation of radial electric field modulus for a gap distance D ¼ 5 mm and a tip voltage V0 ¼ 5 kV: (a) mesh, (b) maps.
electrode distance [17]. This is due to the decrease of the electric field in the vicinity of the plane and therefore the concomitant decrease of the surface charge density on the plane. In this work, the capacity is calculated from the classical equation:
Z C¼
0
R
sðrÞdr V0
(12)
Fig. 12 shows the absolute value of the total electric field along the surface of the tip for different radii of curvature (r). The electric field maxima Emax, as expected, decrease when r increases. These maxima are Emax ¼ 33.80 kV cm1, Emax ¼ 27.80 kV cm1,
Emax ¼ 21.54 kV.cm-1and Emax ¼ 18.92 kV cm1 respectively for r ¼ 1 mm, r ¼ 2 mm, r ¼ 3mmand r ¼ 4 mm. This decrease of the electric field is obviously due to the expansion of the equipotential lines near the tip which is due to the increase of the radius tip. In the same time, this causes the increase of electric field near the opposite plane electrode. This in turn increases the plane surface charge and therefore the increase of the capacitance of the tip-plane gap. This capacitance is respectively equal to C¼0.85 pF, C¼0.90 pF, C¼0.93 pF and C¼0.95 pF in the case of the considered radius (r ¼ 1 mm, r ¼ 2 mm, r ¼ 3 mm and r ¼ 4 mm). Last, for sake of a validation, we show a comparison between the results of the present improved method and some analytical results of the electric potential and field obtained in a specific case versus
28
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
Fig. 7. Two representations of the spatial variation of the total electric field modulus for a gap distance D ¼ 5 mm and a tip voltage V0 ¼ 5 kV: (a) mesh, (b) maps.
Fig. 8. Spatial variation of electric potential (green line) and electric field (blue dashed line) in gap along z-axis for DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 3 mm.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
z-axis. The results of improved fictitious charge method are obtained for a powered electrode when the tip is terminated by a half sphere of radius r ¼ 4 mm and for a gap distance D between the tip of the semi-sphere and the plane. The analytical expression of the electric potential U(z)generated along the z-axis in this semispherical/plane configuration can be written as follows:
Fig. 9. Locations and magnitudes of punctual fictitious charges on z-axis in the powered electrode axis for DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 3 mm.
Fig. 10. Variation of the surface charge density on the plane for DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 3 mm.
UðzÞz
V0 ð2D þ rÞ 1 1 2D jD þ r zj jD þ r þ zj
(13)
Figs. 13 (a) and 14 (a) display the comparison between the analytical and the present numerical solution using fictitious image charges of the electric potential for positive and negative DC voltage while Figs. 13 (b) and 14 (b) display the same comparison in the case of the electric field along the z-axis (i.e. for r ¼ 0). These
Fig. 11. Variation of the configuration capacity as a function of gap distance for curvature radius r ¼ 3 mm.
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
29
any change on the relative error between analytical and numerical results. Last, if the tip radius is chosen much smaller than the interelectrode distance (<1 mm) for conditions close to tip radius currently used in corona discharge setups at atmospheric pressure, the present image charge numerical method is very easily adaptable to such extreme operating parameters. Fig. 16a and b displays an illustration of the sharp electric field variation at the close vicinity of the tip. The magnitude of electric field reaches for instance 165 kV/cm for r ¼ 0.1 mm against 33.8 kV/cm for r ¼ 1 mm. Conclusion
Fig. 12. Variation of the total Electric field modulus along the surface of the tip axis for DC voltage V0 ¼ 5 kV and a gap distance D ¼ 5 mm.
comparisons clearly show a very good agreement in both cases (electric field and potential) since the relative deviation between analytical and numerical results is less than 0.2% for electric potential and less than 1.5% for electric field. The latter is a little bit higher because it accumulates an additional numerical error due to the potential derivative required for the electric field. Therefore, these comparisons can be considered as a reliable validation test for the fictitious image charge method developed in the present work When the applied DC voltage is increased from for instance 5 kV (Fig. 13) to 10 kV (Fig. 15), only the magnitude of potential and electric field are increased as expected by a factor twice without
This paper was devoted to the implementation the image charge method in the case of tip-to-plane electrode configuration and to the improvement and the validation of this interesting method in the case of revolution symmetry geometry. It is shown that this method simulates with a satisfactory manner the electric potential and field. The major advantage of this technique is that it requires no boundary conditions for solving the Laplace equation. The calculation time is very short in the configuration studied in the present work (less than 1 s in a standard laptop with core Intel I7 process). Furthermore, this image charge method allows us very easily to introduce the effect of space charges, very useful for simulation of the streamer corona discharges, provided their spatial distribution is known. This image method is validated from comparison with analytical results known in the case of a semi-sphereto-plane configuration. Finally, it was shown that the punctual fictitious charges give, in comparison to for instance planar or circular distribution, the most accurate results when the electrode has symmetry of revolution.
Fig. 13. Comparison between the results of the improved fictitious charge method and the analytical results in the gap along the z-axis (a) electric potential, (b) electric field for positive DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 4 mm.
Fig. 14. Comparison between the results of the improved fictitious charge method and the analytical results in the gap along the z-axis (a) electric potential, (b) electric field for negative DC voltage V0 ¼ 5 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 4 mm.
D. Benyoucef, M. Yousfi / Journal of Electrostatics 76 (2015) 24e30
V (kV)
8
-0.2
-0.4
4
-0.5
0 0
Analytic method Improved method Error on electric field (%)
-0.3
6
2
b
-10
-20
2
3
4
-0.7 5
z (mm)
1.5 1 0.5
-30
0
-0.6
1
2
error on electric field(%)
Analytic method Improved method Error on potential (%)
E (kV/cm)
a
10
error on potential(%)
30
-0.5 -40 0
1
2
3
4
-1 5
z (mm)
Fig. 15. Comparison between the results of the improved fictitious charge method and the analytical results in the gap along the z-axis (a) electric potential, (b) electric field for positive DC voltage V0 ¼ 10 kV, gap distance D ¼ 5 mm and curvature radius r ¼ 4 mm.
Fig. 16. Variation of the total Electric field modulus in the case of 3 small tip radius for DC voltage V0 ¼ 5 kV and gap distance D ¼ 5 mm: (a) along the z-axis (b) along the surface of the tip axis.
References [1] A. Younes, M. Younes, H. Sayah, M. Bilici, A. Samuila, L. Dascalescu, Effect of spark discharges on the trajectories of insulating particles in roll-type coronaelectrostatic separators. Experimental and numerical study, J. Electrost. 71 (1) (2013) 84e91. [2] Rakshit Tirumala, Go David B, Comparative study of corona discharge simulation techniques for electrode configurations inducing non-uniform electric fields, J. Electrost. 72 (2) (2014) 99e106. [3] X. Wang, N. Mingjiang, X. Gang, Z. Jiapeng, G. Xiang, C. Kefa, An analytical method for DC negative corona discharge in a wire-cylinder device at high temperatures, J. Electrost. 72 (4) (2014) 270e284. [4] M. Kachl, L. Dascalescu, Corona discharges in asymmetric electrode configurations, J. Electrost. 72 (1) (2014) 6e12. [5] B.M. Penetrante, S.E. Schultheis, “Non-thermal Plasma Techniques for Pollution Control,” Part a and B, Editors, Springer-Verlag, New York, 1993. [6] Y.S. Mok, S.W. Ham, I.S. Nam, Evaluation of energy-utilization efficiencies for SO2 and NO removal by pulsed corona discharge process, Plasma Chem. Plasma Process. 18 (4) (1998) 535e550. rimentales des de charges couronne pour la [7] A. Abahazem, Etudes expe pollution des gaz, The se de doctorat de l’, Universite Paul Sabatier, 2009. de [8] J. Chen, P. Wang, Effect of relative humidity on electron distribution and ozone production by DC coronas in air, Plasma Sci. IEEE Trans. 33 (2) (2005) 808e812.
[9] I.V. Timoshkin, M. Maclean, M.P. Wilson, M.J. Given, S.J. Mac Gregor, Tao Wang, J.G. Anderson, Bactericidal effect of corona discharges in atmospheric air, Plasma Sci. IEEE Trans. 2012 40 (10) (2012) 2322e2333. [10] A. Boularas, F. Baudoin, C. Villeneuve-Faure, S. Clain, G. Teyssedre, Multidimensional modelling of electrostatic force distance curve over dielectric surface: influence of tip geometry and correlation with experiment, J. Appl. Phys. 116 (8) (2014) 084106. [11] A.J. Davies, C.S. Davies, C.J. Evans, Computer simulation of rapidly developing gaseous discharges, Proc. IEE 118 (1971) 816e823. [12] M. Yousfi, A. Poinsignon, A. Hamani, Finite element method for conservation equations in electric discharge areas, J. Comp. Phys. 113 (2) (1994) 268e278. [13] M. Talaat, Electrostatic field calculation in air gaps with a transverse layer of dielectric barrier, J. Electrost. 72 (5) (2014) 422e427. [14] D. Benyoucef, M. Yousfi, effects of increasing magnetic field and decreasing pressure on asymmetric magnetron radio frequency discharges, Plasma Sci. IEEE Trans. 41 (4) (2013) 829e838. [15] B.C. Carlson, Numerical computation of real and complex elliptic integrals, Numer. Algorithms 10 (1995) 13e26. [16] M. Yousfi, M.D. Benabdessadok, Boltzmann equation analysis of electronmolecule collision cross sections in water vapor and ammonia, J. Appl. Phys. 80 (12) (1996) 6619e6631. [17] D. Raouti, S. Flazi, D. Benyoucef, electrical modelling of a positive point to plane corona discharge at atmospheric pressure, Contributions Plasma Phys. 54 (10) (2014) 851e858.