Fast computation of complex error function of the real argument
A new efficient technique for evaluation of complex error function of the real argument on the base of the Euler-Maclaurin formula and a non-singular formula for the principal value of the Cauchy integral is given. . It is demonstrated that for the real values of argument, which are most “expensive”...
Збережено в:
| Дата: | 2017 |
|---|---|
| Автори: | , , |
| Формат: | Стаття |
| Мова: | English |
| Опубліковано: |
Національний науковий центр «Харківський фізико-технічний інститут» НАН України
2017
|
| Назва видання: | Вопросы атомной науки и техники |
| Теми: | |
| Онлайн доступ: | http://dspace.nbuv.gov.ua/handle/123456789/122121 |
| Теги: |
Додати тег
Немає тегів, Будьте першим, хто поставить тег для цього запису!
|
| Назва журналу: | Digital Library of Periodicals of National Academy of Sciences of Ukraine |
| Цитувати: | Fast computation of complex error function of the real argument / S.V. Malko, S.S. Pavlov, K.K. Tretiak // Вопросы атомной науки и техники. — 2017. — № 1. — С. 76-79. — Бібліогр.: 13 назв. — англ. |
Репозитарії
Digital Library of Periodicals of National Academy of Sciences of Ukraine| id |
irk-123456789-122121 |
|---|---|
| record_format |
dspace |
| spelling |
irk-123456789-1221212017-06-28T03:02:54Z Fast computation of complex error function of the real argument Malko, S.V. Pavlov, S.S. Tretiak, K.K. Фундаментальная физика плазмы A new efficient technique for evaluation of complex error function of the real argument on the base of the Euler-Maclaurin formula and a non-singular formula for the principal value of the Cauchy integral is given. . It is demonstrated that for the real values of argument, which are most “expensive” in the computational time, the new algorithm is almost two times faster than algorithm 680 with the same accuracy. The code was successfully implemented to ray tracing code that is used for application for investigation of ICR plasma heating in JET-type tokamaks. Предлагается новый эффективный метод вычисления комплексной функции ошибок реального аргумента на основе несингулярной интегральной формы Коши и формулы Эйлера-Маклорена. Показано, что для действительных значений аргумента, которые являются наиболее "дорогими" по времени вычислений, новый алгоритм примерно в два раза быстрее, чем алгоритм 680 с той же точностью. Алгоритм был успешно использован в коде лучевых траекторий, который применялся для исследования ИЦР-нагрева плазмы в токамаке типа JET. Пропонується новий ефективний метод обчислення комплексної функції помилок реального аргументу на основі несінгулярної інтегральної форми Коші та формули Ейлера-Маклорена. Показано, що для дійсних значень аргументу, які є найбільш "дорогими" за часом обчислень, новий алгоритм приблизно в два рази швидший, ніж алгоритм 680 з тією ж точністю. Алгоритм був успішно використаний в коді променевих траєкторій, який застосовувався для дослідження ІЦР-нагріву плазми в токамаці типу JET. 2017 Article Fast computation of complex error function of the real argument / S.V. Malko, S.S. Pavlov, K.K. Tretiak // Вопросы атомной науки и техники. — 2017. — № 1. — С. 76-79. — Бібліогр.: 13 назв. — англ. 1562-6016 PACS: 52.27.Ny http://dspace.nbuv.gov.ua/handle/123456789/122121 en Вопросы атомной науки и техники Національний науковий центр «Харківський фізико-технічний інститут» НАН України |
| institution |
Digital Library of Periodicals of National Academy of Sciences of Ukraine |
| collection |
DSpace DC |
| language |
English |
| topic |
Фундаментальная физика плазмы Фундаментальная физика плазмы |
| spellingShingle |
Фундаментальная физика плазмы Фундаментальная физика плазмы Malko, S.V. Pavlov, S.S. Tretiak, K.K. Fast computation of complex error function of the real argument Вопросы атомной науки и техники |
| description |
A new efficient technique for evaluation of complex error function of the real argument on the base of the Euler-Maclaurin formula and a non-singular formula for the principal value of the Cauchy integral is given. . It is demonstrated that for the real values of argument, which are most “expensive” in the computational time, the new algorithm is almost two times faster than algorithm 680 with the same accuracy. The code was successfully implemented to ray tracing code that is used for application for investigation of ICR plasma heating in JET-type tokamaks. |
| format |
Article |
| author |
Malko, S.V. Pavlov, S.S. Tretiak, K.K. |
| author_facet |
Malko, S.V. Pavlov, S.S. Tretiak, K.K. |
| author_sort |
Malko, S.V. |
| title |
Fast computation of complex error function of the real argument |
| title_short |
Fast computation of complex error function of the real argument |
| title_full |
Fast computation of complex error function of the real argument |
| title_fullStr |
Fast computation of complex error function of the real argument |
| title_full_unstemmed |
Fast computation of complex error function of the real argument |
| title_sort |
fast computation of complex error function of the real argument |
| publisher |
Національний науковий центр «Харківський фізико-технічний інститут» НАН України |
| publishDate |
2017 |
| topic_facet |
Фундаментальная физика плазмы |
| url |
http://dspace.nbuv.gov.ua/handle/123456789/122121 |
| citation_txt |
Fast computation of complex error function of the real argument / S.V. Malko, S.S. Pavlov, K.K. Tretiak // Вопросы атомной науки и техники. — 2017. — № 1. — С. 76-79. — Бібліогр.: 13 назв. — англ. |
| series |
Вопросы атомной науки и техники |
| work_keys_str_mv |
AT malkosv fastcomputationofcomplexerrorfunctionoftherealargument AT pavlovss fastcomputationofcomplexerrorfunctionoftherealargument AT tretiakkk fastcomputationofcomplexerrorfunctionoftherealargument |
| first_indexed |
2025-07-08T21:09:45Z |
| last_indexed |
2025-07-08T21:09:45Z |
| _version_ |
1837114575664709632 |
| fulltext |
ISSN 1562-6016. ВАНТ. 2017. №1(107)
76 PROBLEMS OF ATOMIC SCIENCE AND TECHNOLOGY. 2017, № 1. Series: Plasma Physics (23), p. 76-79.
FAST COMPUTATION OF COMPLEX ERROR FUNCTION OF THE
REAL ARGUMENT
S.V. Malko1, S.S. Pavlov2, K.K. Tretiak2
1Centro de Laseres Pulsados (CLPU), Salamanca, Spain;
2Institute of Plasma Physics of the NSC KIPT, Kharkov, Ukraine
A new efficient technique for evaluation of complex error function of the real argument on the base of the Euler-
Maclaurin formula and a non-singular formula for the principal value of the Cauchy integral is given. . It is
demonstrated that for the real values of argument, which are most “expensive” in the computational time, the new
algorithm is almost two times faster than algorithm 680 with the same accuracy. The code was successfully
implemented to ray tracing code that is used for application for investigation of ICR plasma heating in JET-type
tokamaks.
PACS: 52.27.Ny
INTRODUCTION
The function ])exp(/21)[exp()(
0
22
z
dttizzw
is known as the complex error function or as Fadeeva
function. The name “Fadeeva function” was probably
introduced by Poppe and Wijers [1], since the function
was firstly tabulated by Fadeeva and Terent’ev in 1954
[2]. This function appears in many areas of physics and
mathematics in different forms. For example,
considering such applications as plasma spectroscopy,
nuclear physics, radiative heat transfer and nuclear
magnetic resonance it is necessary efficient and accurate
evaluation of this function. Depending on the problem
there are two main trends of the evaluation techniques:
high accuracy and computational time. We will consider
applications when computational time is more important
than accuracy.
Along with the applications, rather abundant are the
methods evaluating this function, from tables [2, 3] to
modern software [1, 4-7]. Still the main requirements
are high efficiency or somewhat accuracy that should be
investigated in the developing of any successful
algorithm for the computation of this function.
For instance, this function is common in the theory
of plasma waves, and its evaluation is a necessary
ground of ion cyclotron wave analysis of the laboratory
fusion plasmas. The nonrelativistic plasma dispersion
function )(zZ , which describes the absorption and
dispersion properties of plasmas along the magnetic
field is related to the function )(zw as )()( zwizZ .
For this reason the function )(zw is named as well as
the plasma dispersion function, which appears in all
elements of the plasma dielectric tensor. As a rule, in
plasma wave applications the function )(zw is
evaluated massively, therefore the efficiency of
involved numerical algorithm is of primary importance.
So, plasma applications can be related to the one of two
main trends mentioned above: they are connected with
the provision of high efficiency of calculations with a
previously specified accuracy.
At present time, the most successful method of this
trend is the method of Gautschi-Poppe-Wijers [1, 4].
The main idea of this method is to use the continued
fractions. For example, the continued fraction of Jacobi
...
2/312/11
)(
zzzz
zw
has been proved to provide the fast calculation of the
complex error function by means of this method for
large |z| values (Region Q, Fig. 1). For moderate |z|
values (Region R, see Fig. 1) were used the same
continued fractions in combination with the Taylor
expansion in a special direction. For small |z| values
(Region S, see Fig. 1) the Taylor expansion at the zero
point
0
12
2
!)12(
2
1)exp()(
n
n
nn
zi
zzw
was used.
Fig. 1. Three regions of the first quadrant of a complex
plane, where W(z) is calculated on the base three
techniques
It worth to mention that it is enough to calculate this
function in any quadrant due to its properties of
symmetry.
This method provides the accuracy up to 14
significant digits and the average computation time of
the single function value approximately equal to 10
times of the single computation of the exponential
function. It worth to note, the computation time in this
method is the most “expensive” on the real axis and
near it in the Region R.
ISSN 1562-6016. ВАНТ. 2017. №1(107) 77
The main purpose of the present work is to increase
the efficiency of calculation in the most important
region R.
1. EVALUATION OF THE COMPLEX
ERROR FUNCTION OF THE REAL
ARGUMENT
In the theory of plasma waves, solving boundary
value problems requires computation of complex error
function only on the real axis, so it makes sense to
consider this case separately. On the real axis this
function can be defined by the means of the formula [8]
.
1
)(
2
xt
dt
t
e
i
xw
(1)
This function belongs to the class of analytic
functions defined by Cauchy integral, with the integral
density tending to 0 at infinity. In accordance with the
formula of Plemelj-Sokhotskii real part of (1) equals to
the density of the integral at the pole xt multiplied
by i i.e.,
2xe and its imaginary part, named usually
as the Dawson integral, equals to the principal value of
this integral. In principle, for the IIx values in the
Region R, this principal value can be evaluated by
means of the direct numerical integration on the base of
the following integral form which is nonsingular on the
real axis [9]:
2 2 2( 2 )1 1
Im .
xt t t xe e e
w(x) iD(x) dt dt
i t x i t x
(2)
The letter P denotes the principal value of the next
integral. Since the numerator of the integrand in (2) is
expressed in terms of exponentials, it can be estimated
numerically up to the any required precision.
Fortunately, the formulae (2) has an extra interesting
property, which is manifested only in its use in
conjunction with the famous Euler-Maclaurin formula
,
1
)(
)12(
)(
)12(
)!2(
2
2
)()(
!1
)...()( 1
11
m
a
m
fb
m
f
m
mh
m
B
afbf
b
a
hB
yyhdxxf n
(3)
)(,1,...,2,1),(0 bfnynyyyafy , are the
values of function )(xf in successive equally spaced
points with the step h , and mB are the Bernoulli
coefficients. Let’s apply the formula (3) to the right side
of formula (2) by setting
.)(
1
22 )2(
tf
xt
ee
i
xtt
(4)
Then the Dawson integral can be presented as
.)()(
)!2(
)()(
!1
)...()()(Im
1
)12()12(
2
2
1
11
m
mm
m
m
n
x
fxf
m
hB
fxf
hB
yyhdttfxw
(5)
If we present function )(xw as the expansion at the
point x
0
12
1
)(
k
kx
ki
xw
(6)
we can conclude the even properties of )(tf with
respect to the point x: )()( hxfhxf for
every 0h . It’s easy to see, that
)...3,2,1(,0)()12( mf m and 0)( f .
Taking into account that 2/1
1
B and
0)()12( xf m we have
,
2
)...(
1
)(Im 11
)2( 22
nn
x xtt
y
ih
yy
ih
dt
xt
ee
i
xw
(7)
where
2( ) 4(1 )
( ) .
x jh xjh
n j
e e
y f x jh
jh
Since the limiting behavior of Bernoulli coefficients
is m
m mOB )2/(! , this formula is exact for some
1h rather close to 1. Direct calculations along the
formula (7) show that the value 5.0h provides the
accuracy up to 14 significant digits. This allows one
performing computations about ten times faster in
comparison with the direct numerical integration on the
base the single nonsingular form (2).
2. COMPUTATIONAL PROCEDURE, TESTS
After transformation of the final expression (7) for
evaluation of the Dawson integral in order to increase
the efficiency of this algorithm we can have reduction
of the calculation up to only three exponents:
2xe ,
25.0e , xe and some additional power algebra.
Numerical calculations have shown that efficiency of
evaluation of the complex error function on the real axis
is higher than algorithm of Gautschi-Poppe-Wijers
(GPW).
From Table it follows that calculation of )(xw for
most “expensive” in computational time real values of
argument from the region R demonstrates about twice
higher efficiency in comparison with the algorithm
[1, 4].
Comparison of computational times, which are required
for estimation of function )(xw in the points
6,5,4,3,2x from region R for our algorithm (the 2nd
column) and the algorithm of GPW (the 3d column)
with accuracy up to 14 significant digits
x-values of
argument
Time of our
computation
Time of GPW
computation
2.0 5.4990999×10-2 0.160976008
3.0 5.5991001×10-2 0.154975995
4.0 5.3992003×10-2 0.138979003
5.0 5.2992001×10-2 0.123980999
6.0 5.3991001×10-2 9.6986003×10-2
3. JET-APPLICATION
The technique of the plasma dispersion function
evaluation, presented above, was utilized to accelerate
numerical calculations in the ray tracing code as an
example of application in the field of nuclear fusion
78 ISSN 1562-6016. ВАНТ. 2017. №1(107)
[10]. This code was used for the theoretical
investigation of ion cyclotron resonance frequency
plasma heating in a tokamak of JET-type [11]. It was
assumed that the antenna to excite the fast mode (FW)
of the fast magnetosonic wave was located at the inner
side of the plasma column. We investigated the scenario
of plasma heating by means of FW excitation and its
subsequent conversion to the ion Bernstein wave (IBW)
in a plasma depth with subsequent IBW absorption in
the central regions of the plasma column. This scenario
is rather promising for high-power plasma heating in
large size tokamaks because the RF energy is easily
introduced into the FW and then fully absorbed by the
IBW. The energy is easy injected into FW due to
narrow area of FW opacity in the periphery of the
plasma and then after conversion the entire energy is
completely absorbed by the IBW which is linearly
polarized and has a wavelength much smaller than that
of fast wave. Note that direct excitation of IBW is
currently the serious problem due to the wide opaque
barrier on the periphery of the plasma column [12].
Ray tracing code allows a treatment of general
tokamak magnetic field configurations with the poloidal
divertor and realistic wall geometry on the base the
equilibrium field data from the EFIT code [13]. In
calculations was used cylindrical coordinate system
),,( ZRrr
.
The geometrical optics equations describing the
propagation of rays in spatially non-uniform plasma are:
/
/
D
kD
dt
rd
,
/
/
D
rD
dt
kd
, (8)
where k
is the wave vector and ),,,( || rNNDD
is
dispersion relation for waves. Here, /kcN
is
refractive index, BBNN /)(||
and
||NNN
are
longitudinal component of the refractive index along the
magnetic and the perpendicular field, respectively. In
the region of ion cyclotron resonance frequency the
dispersion relation for plasma waves contains the next
components of the dielectric tensor
2
11 Nc ,
i BiBi
Tipi
c
v
)4)((
3
22222
222
, (9)
i Bi
pic
22
2
1 1
,
i BiBi
pi
)( 22
2
2
,
2
2
3 1
pe
))z(Wzi +1(2 0e0e
2
0 ez .
Here
Tee vkz ||0 2 and )(W is the plasma
dispersion function.
Power carried by the ray Q is defined as
eQQ 0
where
0Q is initial power and
dt
t
r
k
T
o
Im . (10)
The following parameters of the plasma and excited FW
were used in the numerical calculations: deuterium
plasma with the addition of hydrogen
( )0(/ dh nn 5…16 %),
313101.3 cmne is the density
of the deuterium in the core, Ti,e = 1…3 keV are the ion
deuterium and hydrogen) and electron temperatures,
B(0) = 1.8 T is magnetic field in the center of column,
8l is toroidal wave number, 810822.1 srad / is
FW frequency (Fig. 2).
Fig. 2. Ray trajectories in JET-like tokamak for cases:
)0(/ dh nn 16% )0(,eiT =1 keV (1); )0(/ dh nn 16%,
)0(,eiT =2 keV (2); )0(/ dh nn 16%, )0(,eiT =3 keV (3);
( )0(/ dh nn 5%, )0(,eiT =3 keV ( 4)
A detailed analysis of that RF plasma heating
scenario in JET-like tokamak is not the subject of
thispaper. Moreover, the full investigation of this
scenario requires as well the use of the 1-D wave model
to solve boundary problem for the system of Maxwell-
Vlasov equations, that also leads to massive
computations of the plasma dispersion function.
However, it is easy to see that the application of the
technique presented in this paper for the plasma
dispersion function evaluation allows approximately
twice to accelerate calculations on the base both models.
CONCLUSIONS
1. On the base of nonsingular form for Cauchy
integrals and the Euler-Maclaurin formula a new method
for computing the complex error function for real values
of argument is given.
2. Comparison of this method with the most
successful one of Gautschi -Poppe-Wijers has shown that
doubled efficiency for real argument values in the most
“expensive” in terms of computational time region R.
3. This method can be used as well for computation
of any Cauchy integrals defined on the real axis,
provided that the density of the corresponding integral
vanishes at infinity.
ISSN 1562-6016. ВАНТ. 2017. №1(107) 79
4. The efficiency of using this method was shown in
application for investigation of ICRF plasma heating in
JET-type tokamaks, where it has been implemented to
ray tracing code. As far as this code usually featured for
the study of the shortwave IBW in tokamak JET, it leads
to a large number of trajectories and a large number of
calculations of the plasma dispersion function during
calculation of each trajectory. For acceleration of such
calculations this method was applied to.
REFERENCES
1. G. Poppe & C. Wijers. More Efficient Computation
of the Complex Error Function // ACM Transactions on
Mathematical Software. 1990, v. 16 (1), p. 38-46.
2. V. Fadeev, N. Terent’ev. Tables of values of the
function ])exp(/21)[exp()(
0
22
z
dttizzw for
complex argument. Moscow: “Gosud. Izdat. Teh.-Teor.
Lit.”, 1954.
3. B.D. Fried, S.D. Conte. The Plasma Dispersion
Function. The Hilbert Transform of Gaussian. New
York: “Academic Press”, 1961.
4. W. Gautschi. Efficient computation of the complex
error function // SIAM Journal on Numerical Analysis.
1969, v. 7, p. 187-198.
5. J.A.C. Weideman. Computation of the Complex Error
Function // SIAM Journal on Numerical Analysis. 1994,
v. 31(5), p. 1497-1518.
6. M. Zaghoul, A. Ali. Algorithm 916 Computing the
Fadeeva and Voight Functions // ACM Transactions on
Mathematical Software. 2011, v. 38(2), p. 1-22.
7. S. Abrarov, B. Quine. Efficient algorithmic
implementation of the Voigt/complex error function
based on exponential series approximation // Appl.
Math. Comp. 2011, v. 218(5), p. 1894-1902.
8. M. Abramowitz, I. Stegun. Handbook of
Mathematical Functions. New York: “National Bureau
of Standards”, 1964.
9. F. Castejón, S. Pavlov. Relativistic plasma dielectric
tensor evaluation based on the exact plasma dispersion
function concept // Physics of Plasmas. 2006, v. 13.
10. E.R. Tracy, A.J. Brizard, A.S. Richardson,
A.N. Kaufman. Ray Tracing and Beyond. Cambridge:
“University Press”, 2014.
11. http://www.ccfe.ac.uk/JET.aspx
12. A.V. Longinov, K.N. Stepano. Radio-frequency
heating in the ion cyclotron frequency range, p. 93-238.
The book: High-Frequency Plasma Heating, Editor
A.G. Litvak. AIP: “New York”, 1992.
13. https://fusion.gat.com/theory/Efit.
Article received 25.09.2016
БЫСТРОЕ ВЫЧИСЛЕНИЕ КОМПЛЕКСНОЙ ФУНКЦИИ ОШИБОК РЕАЛЬНОГО АРГУМЕНТА
C.В. Малко, С.С. Павлов, К.К. Третьяк
Предлагается новый эффективный метод вычисления комплексной функции ошибок реального
аргумента на основе несингулярной интегральной формы Коши и формулы Эйлера-Маклорена. Показано,
что для действительных значений аргумента, которые являются наиболее "дорогими" по времени
вычислений, новый алгоритм примерно в два раза быстрее, чем алгоритм 680 с той же точностью. Алгоритм
был успешно использован в коде лучевых траекторий, который применялся для исследования ИЦР-нагрева
плазмы в токамаке типа JET.
ШВИДКЕ ОБЧИСЛЕННЯ КОМПЛЕКСНОЇ ФУНКЦІЇ ПОМИЛОК РЕАЛЬНОГО АРГУМЕНТУ
C.В. Малко, С.С. Павлов, К.К. Трет’як
Пропонується новий ефективний метод обчислення комплексної функції помилок реального аргументу
на основі несінгулярної інтегральної форми Коші та формули Ейлера-Маклорена. Показано, що для дійсних
значень аргументу, які є найбільш "дорогими" за часом обчислень, новий алгоритм приблизно в два рази
швидший, ніж алгоритм 680 з тією ж точністю. Алгоритм був успішно використаний в коді променевих
траєкторій, який застосовувався для дослідження ІЦР-нагріву плазми в токамаці типу JET.
http://www.ccfe.ac.uk/JET.aspx
https://fusion.gat.com/theory/Efit
|