Improved linearized model collision operator for the highly collisional regime

The linearized model collision operator for multiple species plasmas given by H. Sugama, T.-H. Watanabe, and M. Nunami [Phys.\ Plasmas {\bf 16}, 112503 (2009)] is improved to be properly applicable up to the highly collisional regime. The improved linearized model operator retains conservation laws of particles, momentum, and energy as well as it reproduces the same friction-flow relations as derived by the linearized Landau operator so that this model can be used to correctly evaluate neoclassical transport fluxes in all collisionality regimes. The adjointness relations and Boltzmann's H-theorem are exactly satisfied by the improved operator except in the case of collisions between unlike particle species with unequal temperatures where these relations and H-theorem still holds approximately because there is a large difference between the masses of the two species with significantly different temperatures. Even in the unequal-temperature case, the improved operator can also be modified so as to exactly satisfy the adjointness relations while it causes the values of the friction coefficients to deviate from those given by the Landau operator. In addition, for application to gyrokinetic simulations of turbulent transport, the improved operator is transformed into the gyrophase-averaged form with keeping the finite gyroradius effect.

A well-established Coulomb collision term is given by the Landau operator 31 which is nonlinear for like-species collisions or bilinear for unlike-species collisions. The linearized Landau operator 32-34 obtained by perturbatively expanding the distribution functions about the local Maxwellian is more tractable than the full Landau operator [35][36][37] and the former is preferred to be used for transport studies when the deviation from the Maxwellian is sufficiently small. Since the field particle part of the linearized Landau operator is not as easy to evaluate as its test particle part, several linearized model collision operators have been proposed, in which simplified versions of the field particle part are used. [17][18][19][20][21][22][23] As an example, Sugama et al. 23 presented a linearized model collision operator for multiple ion species plasmas which conserves particles, momentum, and energy, and satisfies adjointness relations and Boltzmann's H-theorem even for collisions between unlike particle species with unequal temperatures. This model called the Sugama operator has been successfully applied to studies of neoclassical and turbulent transport in relatively low collisional regimes. [38][39][40][41][42][43][44][45] The difference between the field particle part of the Sugama operator and that of the exact linearized Landau operator is anticipated to increase in a highly collisional regime. Even in very-high-temperature fusion plasmas like the ITER plasma, 14,15,46 minority impurity ions such as tungsten are considered to remain in the Pfirsch-Schlüter regime even though bulk ions and electrons are in the banana regime. For such a case, it is necessary to use a collision model which is accurate in all collisionality regimes. In this work, the Sugama operator is improved to present the new linearized model collision operator, which is properly applicable to all cases from low to high collisionality. The improved model is constructed so as to give exactly the same friction-flow relations as those derived from the linearized Landau operator. Therefore, it can be used in drift kinetic simulations to accurately evaluate neoclassical transport fluxes in all collisionality regimes. Then, it is noted that the exact friction-flow relations no longer rigorously keep the symmetry property in the case of collisions between unlike particle species with unequal temperatures, where neither the improved model operator nor the linearized Landau operator is completely self-adjoint. Since the self-adjointness is practically useful for analytical or numerical derivation of the Onsager symmetric neoclassical transport coefficients, [1][2][3][4][5][47][48][49][50] further modification of the improved model for the unequal-temperature case is considered in the present paper to restore the adjointness relations by relaxing the accuracy of the friction-flow re-lations. In addition, the improved collision operator in the form suitable for application to gyrokinetic simulations of turbulent transport is derived by taking the gyrophase average with the finite gyroradius effect taken into account.
It is instructive to note here that Hirshman and Sigmar 16 presented a linearized model collision operator which is similar to ours in that spherical harmonic functions and Laguerre polynomials are used to expand distribution functions as well as key properties of the original linearized Landau collision operator are retained. In their work, 16 an elegant and skillful method of constructing novel basis functions is presented to approximate both test and field particle operators including spherical harmonic functions of all degree numbers (l's) although an explicit expression of their model collision operator is given in their paper only for the case where spherical harmonic functions of degrees l > 2 are dropped. To satisfy conservation laws of momentum and energy, the field particle part of the Sugama collision operator contains the l = 0 and 1 parts which are expressed using the test particle part and take similar forms to those of the Hirshman-Sigmar operator. In the present paper, the improved Sugama collision operator is given by adding the correction terms into only the l = 1 spherical harmonic component of the original Sugama operator in order to correctly reproduce the friction-flow relations which determine collisional transport and influence turbulent transport through controlling micro/macroscopic plasma flow profiles. However, the procedures shown in the present work can be extended to give correction terms to all other spherical harmonic components. It is also shown by Abel et al. 22 that, when the Hirshman-Sigmar model operator is transformed to its gyrophaseaveraged form for application to the gyrokinetic equation, the problematic gyroradius dependence appears in the energy diffusion term in the test particle operator. Therefore, for the gyrokinetic case, our model operator is more favorable than the Hirshman-Sigmar model operator.
The rest of this paper is organized ad follows. In Sec. II, we briefly explain the Landau collision operator and its linearization, from which the associated matrix elements are defined to obtain the friction coefficients entering the friction-flow relations. Then, after reviewing the definition and properties of the original Sugama operator in Sec. III, its improved version is presented in Sec. IV, where we write down the correction term to reproduce the same matrix elements and friction coefficients as given by the linearized Landau operator. In Sec. V, the improved operator is expressed in the form suitable for gyrokinetic equations. Finally, conclusions are given in Sec. VI. In Appendix A, a collisional energy transfer rate between unlike species with unequal temperatures is estimated depending on the ratio between the masses of the two species. In Appendix B, effects of unequal temperatures of colliding particle species on the adjointness relations and matrix elements associated with the linearized Landau operator are discussed. The detailed expressions of the matrix elements are shown in Appendix C. In addition, Appendix D presents a modified version of the improved operator which exactly satisfies the adjointness relations even for collisions between unlike particle species with unequal temperatures although it consequently makes the values of the friction coefficients deviate from those given by the Landau operator.

II. LANDAU COLLISION OPERATOR AND FRICTION-FLOW RELATIONS
The Landau operator for collisions between particle species a and b is written as 31 where and ln Λ is the Coulomb logarithm. The particle mass and charge are denoted by m s and e s , respectively, where the particle species is denoted by the subscript s(= a, b).
The distribution function f s (v) generally depends not only on the velocity v but also on the position and time variables (x, t) although the dependence on (x, t) are not explicitly shown here. Writing the distribution function by the sum of the equilibrium part and the small perturbation part as f s = f s0 + δf s , we obtain where the last term C ab (δf a , δf b ) is neglected hereafter. We now assume the equilibrium distribution functions to take the Maxwellian form where n s is the density, v T s ≡ (2T s /m s ) 1/2 is the thermal velocity, and T s is the temperature. Then, the first term on the right-hand side of Eq. (3) is written as The collision time τ ab is defined by . It is easily seen that C ab (f aM , f bM ) vanishes for T a = T b . We hereafter assume that T a /T b = O(1). When m a /m b = O(1), we have α ab = O(1) and C ab (f aM , f bM ) ∼ −(T a /T b − 1)f aM /τ ab . In this case, as explained in Appendix A, we may consider that collisions cause species a and b to have the equal temperature T a = T b after a time scale longer than τ ab .
The second and third terms on the right-hand side of Eq. (3) are called the test and field particle parts, respectively, and the sum of them gives the linearized collision operator, We now expand the perturbed distribution functions δf s (s = a, b) as where Y m l (θ, ϕ) represent spherical harmonic functions and (v, θ, ϕ) are spherical coordinates in the velocity space. The l = 1 component δf (l=1) s of the distribution function δf s is further expanded in terms of the Laguerre polynomials L where x s ≡ v/v T s . The flow vectors u sj (j = 0, 1, 2, · · · ) are defined by For j = 0 and j = 1, we can write u s0 = u s and u s1 = −(2/5)(q s /p s ), where u s ≡ n −1 ) represent the fluid velocity and the heat flow, respectively. We next use the Laguerre polynomials to expand the l = 1 spherical harmonic component of the collision operator before deriving the friction-flow relations in Eq. (12). The resultant expansion [given below in Eq. (9)] contains the coefficients (denoted by C abj ) as functionals of distribution functions, into which the expression in Eq. (7) is substituted to define the matrix elements M jk ab and N jk ab for representing the the friction coefficients l ab ij later. The l = 1 component of the collision term in Eq. (5) is written as Here, C abj (j = 0, 1, 2, · · · ) are defined by where the matrix elements M jk ab and N jk ab (j, k = 0, 1, 2, · · · ) are given from the test and field particle operators, respectively, as 3 In Eq. (11), v denotes the velocity component parallel to the background magnetic field although it can be replaced with the velocity component in any other direction because of the spherical symmetry of the collision operator.
Using the linear collision operator, the friction forces F ai (i = 1, 2, · · · ) are given by 3 Here, the first two-order friction forces are written as . The friction coefficients l ab ij (i, j = 1, 2, · · · ) are defined by 3 where δ ab denotes the Kronecker delta (δ ab = 1 for a = b and δ ab = 0 for a = b). From the momentum conservation in collisions [see Eq. (B8) in Appendix B], we obtain and a l ab 1j = 0 (j = 1, 2, · · · ).
The adjointness relations for the linearized Landau collision operator is written as Strictly speaking, the linearized Landau operator satisfies the adjointness relations in Eq. (16) rigorously only for the case of T a = T b . In this case, the symmetry properties of M ij ab , N ij ab , and l ab ij are derived from Eq. (16) as As explained in Appendix A, T a and T b are significantly different from each other only when m a /m b ≪ 1 or m a /m b ≫ 1. It is explained in Appendix B that, even for this case of unequal temperatures, the adjointness relations given in Eq. (16), the symmetry properties in Eq. (17), and Boltzmann's H-theorem in the form shown later in Eq. (30) are regarded as approximately valid because of the large difference between m a and m b . We also should note that the Onsager symmetry for collisional transport coefficients is derived from the adjointness relations which are also used to give useful methods for solving drift kinetic equations and evaluating neoclassical transport fluxes. 1-5,47-50

III. SUGAMA OPERATOR
The linearized model collision operator for collisions between species a and b given by Sugama et al. 23 is written here as The test-particle part C T S ab (δf a ) of the Sugama operator is defined by where C T 0 ab is defined by Eq. (B2) in Appendix B and the operator Q ab is given for an arbitrary distribution function g by with the dimensionless parameter θ ab , The projection operators P 1a and P 2a is defined by where where C T ab represents the test particle part of the linearized Landau collision operator given in Sec. II. We here note that C T S ab is defined such that the self-adjointness condition, holds exactly even if T a = T b . The field particle part C F S ab (δf b ) of the Sugama operator is given by where and In Eqs. (25)- (26), and are used. We see γ ab = γ ba and η ab = η ba from Eqs. (27) and (28), respectively. It can be easily verified that the test-particle operator C T S ab and the field particle part C F S ab defined in Eqs. (19) and (24) obey conservation laws for particles, momentum, and energy. In addition, C F S ab satisfies the adjointness relation, It is shown in Ref. 23 that the Sugama operator satisfies Boltzmann's H-theorem, We also find that, for the case of m a /m b ≪ 1, C T S ab and C F S ab coincide with C T ab and C F ab of the linearized Landau collision operator to the lowest order in (m a /m b ) 1/2 . For the case of m a /m b ≫ 1, C F S ab equals C F ab to the lowest order in (m b /m a ) 1/2 while, for T a = T b , C T S ab differs from C T ab by the non-self-adjoint part which remains to the lowest order [see the paragraph including Eq. (B16) in Appendix B].
The matrix elements M (S)jk ab and N (S)jk ab (j, k = 0, 1, 2, · · · ) are defined by replacing C T ab and C F ab with C T S ab and C F S ab , respectively, in Eq. (11). Similarly, the friction coefficients l in Eq. (13). Then, from the momentum conservation law satisfied by C T S ab and C F S ab , we obtain and a l (S)ab 1j We also have From the adjointness relations shown in Eqs. (23) and (29), the symmetry properties of the matrix elements and the friction coefficients are derived as In addition, it is found from replacing (C T ab , C F ab ) with (C T S ab , C F S ab ) in Eq. (11) and using Eqs. (24), (31), (33), and (34) that N (S)ij ab (i, j = 0, 1, 2, · · · ) are represented by In Appendix C, we find the detailed expressions of the matrix elements (M ) as well as (M ij ab , N ij ab ). In the moment method, 3,50 the neoclassical transport coefficients, with which the radial particle and heat fluxes and the parallel current are linearly related to the radial density and temperature gradients and the parallel electric field, can be expressed in terms of the viscosity coefficients and the friction coefficients. The friction coefficients l (S)ab ij derived from the Sugama operator do not all coincide with l ab ij given by the Landau operator even for the case of T a = T b (see Sec. III.A below). The dependence of the neoclassical transport coefficients on the friction coefficients becomes relatively strong in the highly collisional regime where accurate values of the friction coefficients need to be derived from the model collision operator for correctly describing the neoclassical transport. In Sec. IV, the improved Sugama operator is presented to reproduce such accurate friction coefficients.
A. Equal temperature case When T a = T b , the test particle part of the Sugama operator is equivalent to that of the Linearized Landau collision operator, which can be easily verified from Eq. (19) with θ ab = 1, Q ab (g) = g, and C T ab = C T 0 ab for that case. In this equal temperature case, we have from which we see and We find from Eqs. (33) and (39) that the deviations of the friction coefficients l (S)ab ij evaluated by the Sugama operator from l ab ij by the Landau operator appear only for i ≥ 2 and j ≥ 2. It is also noted that, for collisions between particles of like species (a = b), the Sugama operator is equivalent to the linearized model collision operator given in Refs. 18,19,22 .

IV. IMPROVED SUGAMA OPERATOR
The improved Sugama operator C

LS(imp) ab
is defined by adding the correction part ∆C LS ab to the original one C LS ab , In order for C LS(imp) ab to reproduce the same friction-flow relations and friction coefficients as those in Eqs. (12) and (13) derived from the linearized Landau collision operator C L ab , the correction term ∆C LS ab (δf a , δf b ) is defined by where the matrix elements M jk ab and N jk ab (M  11) and (13), we can confirm that C LS(imp) ab still gives the same values to M ij ab , N ij ab , and l ab ij as C L ab does, and accordingly, the improved operator correctly reproduces the friction-flow relations in Eq. (12) derived from the Landau operator as well as it retains conservation laws of particles, momentum, and energy. Therefore, we can expect that the classical and neoclassical transport fluxes are accurately evaluated using the improved operator up to the highly collisional regime. As shown in the literature, 4,51 in order to correctly describe the neoclassical transport for the case where all particle species belong to the Pfirsch-Schlüter collisionality regime, we need accurate values for at least the part of the friction coefficients l ab ij with i, j = 1, 2, 3. Accordingly, in this highly collisional case, truncation of the summation j and k in Eqs. (41) and (42) should not be done unless the terms with j ≤ 2 and k ≤ 2 are retained. When the truncation is done such that ∆M 0k ab and ∆N 0k ba with k ≤ k max (k max : an arbitrary integer number) are included, the matrix elements associated with the improved operator satisfy the relations given in the same form as in Eq. (14) or Eq. (31) so that the momentum conservation law still holds as well as the conservation laws of particles and energy. We note here that the correction term ∆C LS ab (δf a , δf b ) is given for only the l = 1 spherical harmonic component of the Sugama operator to correctly reproduce the friction-flow relations which are regarded as an important factor in accurate evaluation of collisional transport and flow profiles influencing turbulent transport. As pointed out in Ref. 3, since the spherical harmonic functions of degree l is an eigenfunction of the pitch-angle-scattering operator included in the test particle collision part with an eigenvalue proportional to −l(l + 1), the test particle part tends to be more dominant over the field particle part as l is larger. Also, in the highly collisional regime, anisotropic components of the distribution function represented by the spherical harmonic functions of higher l's are considered to be stronger damped. Thus, without correction terms in the l ≥ 2 spherical harmonic components of the field particle part, the improved operator is expected to work accurately for describing the distribution function at high collisionality. Besides, in principle, we can extend our procedures to add correction terms to all other spherical harmonic components similarly with the approximation method of Hirshman and Sigmar. 16 Then, the Landau field particle operator is recovered by using spherical harmonic functions and Laguerre functions of all degrees.
As described in Appendix B, the adjointness relations, the resultant symmetry properties for M ij ab , N ij ab , and l ab ij in Eq. (17), and Boltzmann's H-theorem in the form of Eq. (30) are not exactly but only approximately satisfied by the linearized Landau operator and the improved Sugama operator for the case of unequal species temperatures T a = T b . Also, it should be recalled here that the two species need to have very different masses for their temperatures to be significantly different from each other. When the summations j and k in Eqs. (41) and (42) are truncated at the same maximum number j max = k max , the adjointness relations of the improved Sugama operator are still satisfied for T a = T b because the matrix elements M jk ab and N jk ab evaluated by the improved operator still keep the symmetry properties. On the other hand, the H-theorem is not guaranteed by this truncation even for T a = T b because not all but only some fraction of the matrix elements N jk ab (j, k = 0, 1, 2, · · · ) of the Sugama field particle operator are replaced with those of the Landau field particle operator. [Note that the friction-flow relations including all matrix elements M jk ab and N jk ab (j, k = 0, 1, 2, · · · ) given by the linearized Landau operator is equivalent to the l = 1 spherical harmonic part of that Landau operator which satisfies the H-theorem for the l = 1 parts of the distribution functions.] It is shown in Ref. 23 that the H-theorem for the original Sugama operator can be derived from the fact that its field particle part can be completely expressed in terms of the test particle part although the same technique of the derivation of the H-theorem cannot be used for the truncated version of the improved Sugama operator. However, we can still expect that the H-theorem is approximately satisfied by the truncated model if the l = 1 parts of the distribution functions are well represented by the linear combinations of only low-order Laguerre polynomials. As shown in Ref. 4 , sufficiently accurate evaluations of collisional (classical and neoclassical) transport fluxes can be made using the friction-flow relations including the Laguerre polynomial moments up to the order of j = 2. This appears to be because higherorder Laguerre-polynomial components of the distribution functions are stronger suppressed by the energy diffusion operator [see Eq. (B4) in Appendix B]. Therefore, except for the case where the field particle distribution takes a special form due to some external sources, we don't generally expect that higher (j ≥ 3) moments need to be retained in Eq. (41).
It is easily found from the definition of the improved operator in Eqs.  6) and (7)]. Therefore, if δf s (s = a, b) are given by the perturbed Maxwellian with the perturbed densities δn s and temperatures δT s as for which δT a /T a = δT b /T b is assumed, then C

LS(imp) ab
(δf a , δf b ) vanishes as C LS ab (δf a , δf b ) does. 23 However, when δf s (s = a, b) are written as the shifted Maxwellian δf s = f sM (m a /T s )(u s ·v) with the same flow velocity u a = u b and the different equilibrium tempera- (δf a , δf b ) does exactly. This is related to the fact that the symmetry properties l ij ab = l ji ba (i, j = 1, 2, · · · ) are slightly broken when T a = T b (see Appendix D). When the above-mentioned adjointness relations and resultant symmetry properties are satisfied, they provide useful techniques for calculating the neoclassical transport coefficients. 1-5,47-50 Therefore, it will be beneficial for such applications if we can have a linearized collision model which satisfies the adjoint relations exactly even for T a = T b while giving small inaccuracies to the values of the matrix elements and the friction coefficients. Such a model is presented in Appendix D where the correction part of the improved Sugama operator is symmetrized.

V. COLLISION OPERATOR FOR GYROKINETIC EQUATIONS
There are two types of gyrokinetic equations. One is the gyrokinetic equation derived by using the WKB rep-resentation for the perturbed distribution function [52][53][54][55][56][57] which has a high wavenumber in the direction perpendicular to the equilibrium magnetic field B. The other is derived by using the Lie transform technique to properly define the gyrocenter coordinates for description of the total distribution function. 58,59 The collision operator for the former type of the gyrokinetic equation is considered in this section as well as in the literature. [20][21][22][23] On the other hand, several studies have been done to represent the collision operator for the latter type in the gyrocenter coordinates. [24][25][26][27][28][29][30] When applying the improved Sugama operator to the gyrokinetic equation for the perpendicular wavenumber vector k ⊥ , the collision operator is transformed into the following form, where , Ω a ≡ e a B/m a c) and dξ/(2π) represent the gyroradius vector and the gyrophase average, respectively, and h ak ⊥ is given from the nonadiabatic part of the perturbed particle distribution function δf ak ⊥ = −(eφ k ⊥ /T a )f aM + e −ik ⊥ ·ρ a h ak ⊥ . The detailed expression of the first term on the right-hand of Eq. (48) is shown in Ref. 23 while the second term is written as ) denote the zeroth-and first-order Bessel functions of the normalized perpendicular wavenumber k ⊥ v ⊥ /Ω s , respectively, and are used.
In the case of application to the drift kinetic equation for studying neoclassical transport, we neglect the finite gyroradius effects and take the limit k ⊥ · ρ s → 0 (s = a, b). Then, we put J 0s → 1, J 1s → 0, and u ⊥sk [h sk ⊥ ] → 0 in Eqs. (49) and (50).

VI. CONCLUSIONS
In this paper, the improved linearized model collision operator which can be applied up to the highly collisional regime is presented. The improved operator is constructed by adding the correction part to the previous model by Sugama et al. so as to reproduce the same friction-flow relations as those given by the linearized Landau collision operator. In the improved model, conservation laws of particles, momentum, and energy are retained while the adjointness relations and Boltzmann's H-theorem are approximately valid for collisions between unlike particle species with unequal temperatures and very different masses. It is also shown that the improved operator can be modified to satisfy the adjointness relations exactly even in the unequal-temperature case. This modification causes the friction coefficients to deviate from those given by the Landau operator although the influence of the deviations is made small by the very different masses.
Performing the gyrophase average with keeping the finite gyroradius effect, the improved operator is represented in the suitable form for gyrokinetic equations. In the zero-gyroradius limit, the gyrophase-averaged improved operator can be used in drift kinetic equations to accurately evaluate neoclassical transport in all collisionality regimes. It is considered that only the terms with j ≤ 2 in the Laguerre polynomial expansion of the correction part of the operator need to be kept even for the most collisional case where all particle species are in the Pfirsch-Schlüter collisionality regime. The present model is expected to be useful for simulation studies of neoclassical and turbulent transport processes in plasmas including multi-species of particles in various collisional regimes.

ACKNOWLEDGMENTS
The authors sincerely thank Dr. Y. Idomura and Dr. K. Obrejan for beneficial discussion on kinetic simulation of impurity transport. This work is supported in part by JSPS Grants-in-Aid for Scientific Research Grant No. 19H01879 and in part by the NIFS Collaborative Research Program NIFS18KNTT045.

Appendix A: COLLISIONAL ENERGY TRANSFER BETWEEN UNLIKE SPECIES WITH UNEQUAL TEMPERATURES
Using Eq. (4), the collisional energy transfer from species b to a, which equals the opposite sign of that from species a to b, is given by where each species is assumed to be in the local equilibrium state represented by the Maxwellian distribution function. Thus, if T a = T b , collisions cause the temperatures of the two species to approach to each other, and the characteristic rate ν et ab of the collisional energy transfer from species b to a is given by We now assume that |e Therefore, when m a /m b ≪ 1 or m a /m b ≫ 1, collisional energy exchange between species a and b is so slow that f aM and f bM are permitted to have unequal temperatures T a = T b .
We now consider the case where m a /m b ≫ 1, |e a /e b | ≫ 1, n a /n b ≪ 1, and T a /T b = O(1). This can happen when a and b represent heavy minority impurity ions with high charge number and bulk hydrogen isotopes (or electrons), respectively. Then, we have α ab ≪ 1 and For example, using Eq. (A3) in realistic cases as in the JET and ASDEX Upgrade experiments shown in Ref. 12 , the characteristic rates ν et W i and ν et W e of the collisional energy transfer from tungsten impurities (W) to bulk hydrogen isotopes (i) and to electrons (e) are estimated to satisfy and respectively. Under these conditions, it is reasonable to assume T W = T i although it is not to assume T W = T e .

Appendix B: EFFECTS OF UNEQUAL TEMPERATURES ON ADJOINTNESS RELATIONS
Based on the Landau collision operator defined in Eq. (1), the test and field particle parts are written as where Here, C T 0 ab (δf a ) consists of the pitch-angle-scattering part ν ab D (v)L(δf a ) and the energy diffusion part C ab v (δf a ). The pitch-angle-scattering operator L is defined by where I denotes the unit tensor and (v, θ, ϕ) represent spherical coordinates in the velocity space. The energy diffusion operator C ab v is defined by The collision frequencies for pitch-angle scattering and energy diffusion are given by s = a, b). We can easily confirm that ν D L, C ab v , and accordingly C T 0 ab are all self-adjoint so that holds for arbitrary functions δf a and δg a of v. It can also be shown that C F 0 ab satisfies the adjointness relation written as The remaining test and field particle operators C T N ab and C F N ab do not keep adjoint relations such as Eqs. (B5) and (B6) satisfied by C T 0 ab and C F 0 ab , respectively, although C T N ab and C F N ab vanish for T a = T b . We also note that the two pairs of the operators (C T 0 ab , C F 0 ab ) and (C T N ab , C F N ab ) independently satisfy the particle, momentum and energy conservation laws, which are written as and respectively. From the Galilean invariance and spherical symmetry of the Landau collision operator, we have an identity, , f bM (v)] = 0, for an arbitrary vector u which is independent of v. Then, taking the u → 0 limit of the above identity and using the particle and momentum conservation laws, we can derive another type of relations, Then, the momentum conservation law shown in Eq. (B8) is used to find respectively. Also from Eq. (B10), we obtain It should be noted that the symmetry properties of N (0)ij ab and N 00 ab take different forms with respect to the way the temperatures enter.
In the case of m a /m b = O(1), the temperatures T a and T b are expected to be close to each other because of the relatively fast energy exchange due to collisions. Therefore, only when m a ≪ m b or m a ≫ m b , T a can be significantly different from T b . In the limiting case m a ≪ m b , it is shown that the pitch-anglescattering term ν D (v)L(δf a ) is dominant in the test particle operator C T ab (δf a ) where the energy scattering term C ab v (δf a ) and the non-adjoint part C T N ab (δf a ) are negligible in the lowest order of the expansion with respect to (m a /m b ) 1/2 . However, when T a = T b , C F N ab (δf b ) is not negligible but it is necessary to keep contributions from both C F 0 ab (δf b ) and C F N ab (δf b ) for accurately evaluating collisional momentum transfer. Then, it can be shown that, to the lowest order in (m a /m b ) 1/2 , the test and field particle parts of the Sugama operator C S ab (δf a , δf b ) = C T S ab (δf a ) + C F S ab (δf b ) correctly approximate C T ab (δf a ) and C F ab (δf b ) of the linearized Landau operator, respectively.
We next consider the case in which m a ≫ m b and T a = T b hold. In this case, C T N ab (δf a ) is not negligibly small compared with C T 0 ab (δf a ) while C F N ab (δf b ) does not contribute to C F ab (δf b ) in the lowest order of the expansion with respect to (m b /m a ) 1/2 . Then, C F ab (δf b ) is well approximated by either C F 0 ab (δf b ) or C F S ab (δf b ) although the difference of C T ab (δf a ) from C T 0 ab (δf a ) or C T S ab (δf a ) is significant. However, this difference doesn't cause serious errors in solving the kinetic equation for δf a as far as C T ab (δf a )/C T aa (δf a ) ∼ (e b /e a ) 2 (n b /n a )(m b /m a ) 1/2 becomes very small. [This ratio (e b /e a ) 2 (n b /n a )(m b /m a ) 1/2 can be large in such a case of tungsten impurity as mentioned in Appendix A although, for that case, T a = T b is expected so that C T ab (δf a ) = C T 0 ab (δf a ) = C T S ab (δf a ) holds.] Except for this limiting case of m a ≫ m b and T a = T b , we can suppose that the matrix elements M ij ab evaluated by C T ab (δf a ) = C T 0 ab (δf a ) + C T N ab (δf a ) satisfy the symmetry relations of the same form as those for M When m a /m b ≪ 1 and T a = T b , C F N ab (δf a ) makes a significant contribution to C F ab (δf a ) = C F 0 ab (δf a ) + C F N ab (δf a ). In this case, we can show that, to the lowest order in (m a /m b ) 1/2 , are satisfied by the matrix elements N ij ab associated with C F ab (δf a ). Note that the second relation in Eq. (B15), which holds exactly, is included as a special case in the symmetry relations shown by Eq. (B17) and that they take a different form from those for N In summary, the adjointness relations of the linearized Landau operator C L ab are not satisfied rigorously in collisions between unlike species with unequal temperatures although significantly different temperatures occur in the case where the two species have so different masses that the adjointness relations and symmetry properties of the matrix elements and the friction coefficients can still be used as approximately valid formulas. On the other hand, the Sugama operator C LS ab in Sec. III and the operator C LS * (imp) ab in Appendix D are constructed so as to exactly keep the adjointness relations which can be useful in formulating efficient methods of evaluating Onsager symmetric collisional transport coefficients. [1][2][3][4][5][47][48][49][50] [see Eqs. (B11) and (B12) in Appendix B] which are associated with the test part C T ab = C T 0 ab + C T N ab and the field part C F ab = C F 0 ab +C F N ab of the linearized Landau operator C L ab = C T ab +C F ab are expressed in terms of α ab ≡ v T a /v T b , T a /T b , and m a /m b . In addition, it is shown how to evaluate M  (D5) We can now use the test and field particle part of C LS * (imp) ab (δf a , δf b ) to obtain the matrix elements M * ij ab and N * ij ab in the same way as shown in Eq. (11). Then, the friction coefficients l ab * ij can be derived from M * ij ab and N * ij ab [see Eq. (13)]. Since ∆C F * ab (δf b ) defined in Eq. (D2) gives the correction only in the field particle part, we immediately see that