Applied Mathematics
Vol.05 No.10(2014), Article ID:46526,11 pages

A closed-form solution of a kinetic integral equation for rarefied gas flow in a cylindrical duct

Carmo Henrique Kamphorst1, Patricia Rodrigues2, Liliane Basso Barichello3

1Departamento de Ciências Exatas e da Terra, Universidade Regional Integrada do Alto Uruguai e das Missões (URI), Frederico Westphalen, Brazil

2Departamento de Ciências Agronômicas e Ambientais, Universidade Federal de Santa Maria (UFSM), Frederico Westphalen, Brazil

3Instituto de Matemática, Universidade Federal do Rio Grande do Sul (UFRGS), Porto Alegre, Brazil


Copyright © 2014 by authors and Scientific Research Publishing Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY).

Received 25 March 2014; revised 25 April 2014; accepted 2 May 2014


A spectral method based on Hermite cubic splines expansions combined with a collocation scheme is used to develop a solution for the vector form integral S-model kinetic equation describing rarefied gas flows in cylindrical geometry. Some manipulations are made to facilitate the computational treatment of the singularities inherent to the kernel. Numerical results for the simulation of flows generated by pressure and thermal gradients, Poiseuille and thermal-creep problems, are presented.


Rarefied gas dynamics, integral equation, S-model, collocation schemes

1. Introduction

Research in the field of rarefied gas dynamics (RGD) involving the analysis and description of micro and nanoflows has been fundamental for the development of micro devices and new technologies, specially in the field of micro-electro-mechanical systems (MEMS) [1] -[5] .

The Boltzmann equation is the fundamental mathematical model for dealing with problems in RGD [6] . Investigations of its solution and associated kinetic models have found much attention in the last few years [7] - [13] . Numerical and analytical approaches have been developed, mainly in Cartesian coordinates [14] [15] .

The solution of Boltzmann-like problems is even more challenging in multidimensional and complex geome- tries, particularly, if one seeks for analytical approaches. In cylindrical coordinates, to overcome such difficulties, the integral form of the kinetic equations has been employed. In this context, several approaches have been proposed in the past years [16] - [19] , particularly related to the solution of the BGK model [20] , as, for example, the use of classical Legendre polynomials expansions [21] . A special approach used in an earlier work of Ferziger [22] , is based on a transformation proposed by Mitsis [23] , which allows us to recast the problem in cylindrical coordinates into a much simpler problem in plane Cartesian coordinates. The transformation, however, is restricted to some specific (simpler) physical cases.

Mitsis’ transformation [23] was used by Siewert and Valougeorgis [24] to develop a solution for the integral form of the S-model [25] kinetic equation, which is a vector equation, while the BGK kinetic equation is a scalar problem. It has been addressed in the literature, however, that the S-model is more appropriate than the BGK model to deal, for instance, with nonisothermal flows [24] [26] .

In this work, we use an spectral approach based on Hermite cubic splines expansions along with a collocation scheme to solve the vector integral form of the S-model kinetic equation for flow in cylindrical tube. Numerical results are obtained for the Poiseuille and thermal-creep problems. As in other analytical approaches [24] the solution is explicit, in terms of the spatial variable. However, the technique proposed in this work, that has been also effective to deal with the BGK integral equation [27] , may be extended to deal with the problems which include reflexive boundary conditions.

2. The formulation

The initial formulation of the problems we treat in this work, is the integro-differential formulation associated with the kinetic S-model [24] , for describing the behavior of gas flows in a cylindrical tube of radius, with no variation in the axial direction, such that the distribution of particles and the boundary conditions depends on the spatial variable, written in dimensionless units, and the particle velocity vector, also given in dimensionless units, expressed in the cylindrical coordinates. Still, in the initial formulation, the flow is subject to reflexive boundary conditions at. In Ref. [24] , the authors follow a previous derivation regarding to the BGK model [28] to derive in details the integral form of the equation for the S-model, which describes the flow of a rarefied gas in a straight cylindrical duct of radius, which will be our starting point of this work. Since the derivations were done in details in [24] [28] , we do not repeat it here. We start, already, from considering the integral vector equation [24]


for. In Equation (1) the inhomogeneous term is given by


where the constants and are defined so that and for the so-called problem of Poiseuille flow while and for the thermal-creep problem. Still, the kernel of integral equation is






Here and are the modified Bessel functions of zero order, first and second kind [29] , respectively. Furthermore, in Equation (3),







Now, if we define [24] ,


we can evaluate the quantities of physical interest for the gas flow, such as the velocity profile [24]


and the heat-flow profile


Consequently, we write the particle-flow rate as


and the heat-flow rate


The Equations (11) to (14) are used for both problems, Poiseuille and thermal-creep, depending on the specified values of the constants and which define the expression used for the source term given in Equation (2).

Finally, we note that, in Equation (1), the term brings the information from the boundary condition, and it is not really an inhomogeneous term when reflexive boundary conditions are applied. For this reason, the only case considered in Ref. [24] was the diffuse reflection, where, at the end. At this point, in this work, we also disregard this term in our formulation and will comment on that later on in this text.

3. The solution

Firstly, we rewrite Equation (1) as,


for, and for the case with no specular reflection,.

We then propose a solution to Equation (15) in terms of a truncated expansion of Hermite cubic splines functions [30] , as


where are vectors


whose components (constants) and are to be determined. In regard to the splines functions,

first of all, we consider there to be knots defined on the interval, such that

(see Appendix for details).

If we now substitute Equation (16) into Equation (15) we obtain


which still can be expressed in the form






At this point, we introduce a collocation scheme, such that, if we evaluate Equation (19) at the collocation points


we end up with


for, which is a linear system to the constants and used in expansion

given by Equation (16).

4. Computational Aspects

The procedure described in Section 3 may be considered quite straightforward. However the solution of the

linear system given by Equations (23), depends on the explicit evaluation of the integrals and

, which are written in terms of the kernel, Equations (3) to (9), which clearly includes singularities

that introduce numerical difficulties. In this way, a series of computational steps have to be taken into account.

As a first step, as used in previous works [24] [31] , for computational reasons, the modified Bessel functions are rewritten as


such that Equations (4) and (5) are now expressed as


for and


for. To be consistent, from now on, we consider the previously given integral equation kernel, Equation (3), expressed in the form


Continuing, if we look back to Equations (6) to (9), we can see that







where, again, the (constant components) matrices, , are defined in Equations (7) to (9). More importantly, in writing the kernel as above, it is to point out that in these definitions, the only expression which effectively involves the term which induces a singular behavior, is the first one in the right hand-side. We note, however, that this same expression was part already of the BGK model kernel we had to deal in a previous work [21] .

Finally, for future use, later on in this text, we can also express the kernel of the integral equation, in a matrix form,


following, of course, the previous definitions given in Equations (7) to (9) and Equations (29) to (31). Here,





Going back to the definition of the linear system we need to solve, Equations (23), we now focus our

attention in the definition of the functions and given in Equations (20) and (21).

Considering that is the support of the spline function, in other words, ,

we note that


if and


if with


In an analogous form,


if and


if with


In fact, to be more specific, we can express the vector, in terms of its components


such that




where and are the constants to be defined in the proposed expansion, Equation (16), obtained as the solution of the linear system Equations (23). The terms, for, are defined in terms of expressions given in Equations (33) to (35),


if and


if. Similarly, the vector, defined by Equations (39) and (40), is rewritten in a matrix form








if and


if. Considering the explicit derivations above, at this point, we can rewrite the linear system defined in Equations (23), in the following form




for, where


for, and are given by the expressions shown in Equations (41)-(46) and

(50)-(51). Still, the terms and are defined by the components of the source term, so that




depending on the problem, as mentioned earlier: and for the Poiseuille flow; and for the thermal-creep problem.

The linear system given by Equations (52) and (53), is of order (or, which solution, the constants and, allow us to evaluate the quantities of physical interest such as velocity profile


and the heat flow profile


Numerical quadrature schemes

Searching for a general procedure, as much as possible, the basic approach for evaluating the integrals involved in this derivation was based on the use of the Gauss-Legendre quadrature scheme [29] with nodes and weights defined in the interval. In this form, convenient mappings, according to the integration interval, were used.

In this way, we define


and we consider a quadrature scheme with nodes and weights for evaluating the integrals defined in Equations (27) to (29) as





At this point, however, a special remark has to be done in regard to the evaluation of Equation (60). While

the general procedure described above works well as long as is not small―in fact, using

one obtains results which agree in seven digits with results obtained by the software Maple―the case has to be treated more carefully, if ones search for more accurate results as a general case.

To deal with that, we subtract (and add) from Equation (29) the expression


In the resulting integral expression, we introduce the change of variables


and, subsequently, once more, we add (and subtract) the integral term


To be clear, when (we considered typically) we evaluate as






The reasons for doing that is, firstly, the use of these additional functions results in a smoother function when approaching. In addition to that, can be expressed as


where is the elliptic integral of the first order, that can be evaluated numerically by applying the algorithm given by [29] . Still, , can be evaluated by using Maple, such that. The integral, given by Equation (67), can be then evaluated without much difficulty by applying the scheme of Gauss-Legendre quadrature, by mapping into.

Therefore, in regard to the computational scheme, we note that is evaluated either by Equation (60) if or Equation (66), if. This procedure is, in fact, analogous to the one we used to deal with the BGK model [21] [27] and it was implemented once we seek for increasing the accuracy of the final results, as we will comment later on in this text.

To all other integrations processes needed in this formulation, more specifically, the expressions given by Equations (46) and (51) as well as Equations (13) and (14), linear mappings from the integration interval to the interval were performed, to the use of the Gauss-Legendre quadrature scheme. Clearly, in the case of Equations (46) and (51) the dependence of the integration limits on the interval definition of the splines functions have to be taken into account to improve the computational results. In this way, we evaluate these equations as


if and


if. Here we used






With this, we conclude the evaluation of the elements of the system given by Equations (52) and (53).

5. Numerical results

As input data in our program we need to specify the radius of the duct, , the order of the quadrature schemes and the number of knots (associated to the splines functions) which will define in the linear system given by Equations (52) and (53). We implemented a FORTRAN program and subroutines from LINPACK [32] were used to solve the linear system. The solution was then used to obtain the physical quantities of interest, the velocity and heat-flow profiles, Equations (57) and (58), and subsequently Equations (13) and (14).

As a first test, to check and to have confidence in our approach, we reproduced all the results listed in Ref. [24] , for the velocity and heat-flow profiles (case), as well as velocity slips and particle flow rates for values of ranging from to. We do not repeat the results here since they were already reported in [24] .

As we mentioned earlier, the approach used in that work can be very accurate, since it works with a Cartesian coordinates problem associated with a very accurate discrete ordinates solution, after using Mitsis’ transfor- mation [23] . In order to get agreement in six to seven digits, for all cases, with Siewert and Valougeorgis’ results

(listed with seven digits) [24] , we used, , and, for

to. If we use, for instance, lower order approximations as, ,

and results can also be reproduced, in general, with two to three digits of agreement with [24] , particularly for.

We then investigated additional cases. To provide the results we list in Table 1 and Table 2, we used the first set of parameters listed above. In regard to Table 2, we remark that, the Onsager reciprocity relation [33] [34] , is verified. Still, the results listed in Table 2 have some agreement with Ref. [26] , within the accuracy

Table 1. Velocity and heat flow profiles in a duct with R = 1.

Table 2. Flow rates.

declared in that paper (0.5%). In addition to provide accurate results, it is important to remark that, the approa- ch proposed in this work may be used to simulate the case which includes specular reflection at the surface, as we show in a subsequent work, and, a better description of the boundary effects is a relevant issue when it comes to the modeling of microflows.

6. Conclusion

In this work, we have used a classical spectral approach, based on the use of Hermite cubic splines combined with a collocation scheme, to obtain a solution, in a closed-form, to the S-model kinetic integral vector equation relevant to describe rarefied gas flow in cylindrical coordinates. The main purpose of this work was to verify the performance of this approach as an extension of an earlier work related to a simpler integral model equation, the BGK model, trying to derive, as much as possible, a general approach which can be now, in a subsequent work, extented to include the treatment of reflexive boundary conditions. The numerical results obtained showed good agreement with accurate results available in the literature.


The authors would like to thank CNPq of Brazil for partial financial support of this work. One of the authors (LBB) would like to thank the Mathematics Department of North Carolina State University for the kind hospitality during a visit when this work was written.


  1. Karniadakis, G., Beskok, A. and Aluru, N. (2005) Microflows and Nanoflows. Springer, New York.
  2. Kakaç, S., Vasiliev, L.L., Bayazitoglu, Y. and Yener, Y. (2005) Microscale Heat Transfer: Fundamentals and Applications. Springer, New York.
  3. Li, D. (2008) Micro and Nanoscale Gas Dynamics. Springer, New York.
  4. Wang, M., Lan, X. and Li, Z. (2008) Analyses of Gas Flows in Micro- and Nanochannels. International Journal of Heat and Mass Transfer, 51, 3630-3641.
  5. Gad-el-Hak, M. (2005) The MEMS Handbook. Mechanical Engineering Handbook Series, CRC Press, New York.
  6. Cercignani, C. (1988) The Boltzmann Equation and Its Applications. Springer-Verlag, New York.
  7. Sone, Y., Ohwada, T. and Aoki, K. (1989) Evaporation and Condensation on a Plane Condensed Phase: Numerical Analysis of the Linearized Boltzmann Equation for Hard-Sphere Molecules. Physics of Fluids A, 1, 1398-1405.
  8. Sharipov, F. and Seleznev, S. (1998) Data on Internal Rarefied Gas Flows. Journal of Physical and Chemical Reference Data, 27, 657-706.
  9. Williams, M.M.R. (2001) A Review of the Rarefied Gas Dynamics Theory Associated with Some Classical Problems in Flow and Heat Transfer. Zeitschrift für Angewandte Mathematik und Physik, 52, 500-516.
  10. Barichello, L.B., Camargo, M., Rodrigues, P. and Siewert, C.E. (2001) Unified Solutions to Classical Flow Problems Based on the BGK Model. Zeitschrift für Angewandte Mathematik und Physik, 52, 517-534.
  11. Barichello, L.B. and Siewert, C.E. (2003) Some Comments on Modeling the Linearized Boltzmann Equation. Journal of Quantitative Spectroscopy and Radiative Transfer, 77, 43-59.
  12. Siewert, C.E. (2003) The Linearized Boltzmann Equation: Concise and Accurate Solutions to Basic Flow Problems. Zeitschrift für Angewandte Mathematik und Physik, 54, 273-303.
  13. Scherer, C.S., Prolo Filho, J.F. and Barichello, L.B. (2010) An Analytical Approach to the Unified Solution of Kinetic Equations in the Rarefied Gas Dynamics. I. Flow Problems. Zeitschrift für Angewandte Mathematik und Physik, 60, 70-115.
  14. Sone, Y. (2002) Kinetic Theory and Fluid Dynamics. Birkhäuser, Boston.
  15. Ivchenko, I.N., Loyalka, S.K. and Tompson, R.V. (2007) Analytical Methods for Problems of Molecular Transport. Springer, New York.
  16. Lang, H. and Loyalka, S.K. (1984) Some Analytical Results for Thermal Transpiration and the Mechanocaloric Effect in a Cylindrical Tube. Physics of Fluids, 27, 1616-1619.
  17. Valougeorgis, D. and Thomas Jr., J.R. (1986) Exact Numerical Results for Poiseuille and Thermal Creep Flow in a Cy- lindrical Tube. Physics of Fluids, 29, 423-429.
  18. Siewert, C.E. (2000) Poiseuille and Thermal-Creep Flow in a Cylindrical Tube. Journal of Computational Physics, 160, 470-480.
  19. Loyalka, S.K. and Tompson, R.V. (2009) The Velocity Slip Problem: Accurate Solutions of the BGK Model Integral Equation. European Journal of Mechanics-B/Fluids, 28, 211-213.
  20. Bhatnagar, P.L., Gross, E.P. and Krook, M. (1954) A Model for Collision Processes in Gases. I. Small Amplitude Processes in Charged and Neutral One-Component Systems. Physical Review, 94, 511-525.
  21. Rodrigues, P., Kamphorst, C.H. and Barichello, L.B. (2009) A Spectral Method for Rarefied Gas Dynamics Problems in Cylindrical Geometry. International Journal of Pure and Applied Mathematics, 51, 181-187.
  22. Ferziger, J.H. (1967) Flow of a Rarefied Gas through a Cylindrical Tube. Physics of Fluids, 10, 1448-1453.
  23. Mitsis, G.J. (1963) Transport Solutions to the Monoenergetic Critical Problems. Argonne National Laboratory, Argonne, Illinois, Report ANL-6787, 370-382.
  24. Siewert, C.E. and Valougeorgis, D. (2002) An Analytic Discrete-Ordinates Solution of the S model in Rarefied Gas Dinamics: Poiseuille and Thermal-Creep Flow in a Cylindrical Tube. Journal of Quantitative Spectroscopy and Radia- tive Transfer, 72, 531-550.
  25. Shakhov, E.M. (1974) Method of Investigation of Rarefied Gas Flows. Nauka, Moscow.
  26. Sharipov, F. (1996) Rarefied Gas Flows through a Long Tube at Any Temperature Ratio. Journal of Vacuum Science & Technology A, 14, 2627-2635.
  27. Kamphorst, C.H., Rodrigues, P. and Barichello, L.B. (2009) A Spectral Approach to Compute Rarefied Gas Flows in Cylindrical Geometry. International Nuclear Atlantic Conference, Rio de Janeiro, 27 September-2 October 2009. (CD-Rom)
  28. Barichello, L.B., Camargo, M., Rodrigues, P. and Siewert, C.E. (2002) An Integral Equation Basic to the BGK Model for Flow in a Cylindrical Tube. Zeitschrift für angewandte Mathematik und Physik, 53, 769-781.
  29. Abramowitz, M. and Stegun, I.A. (1965) Handbook of Mathematical Function. Dover Publications, New York.
  30. Schultz, M.N. (1973) Spline Analysis. Prentice Hall, Englewood Cliffs.
  31. Barichello, L.B., Rodrigues, P. and Siewert, C.E. (2002) An Analytical Discrete-Ordinates Solution for Dual-Mode Heat Transfer in a Cylinder. Journal of Quantitative Spectroscopy and Radiative Transfer, 73, 583-602.
  32. Dongarra, J.J., Bunch, J.R., Moler, C.B. and Stewart, G.W. (1979) LINPACK User’s Guide. SIAM, Philadelphia.
  33. Sharipov, F. (1994) Onsager-Casimir Reciprocity Relations for Open Gaseous Systems at Arbitrary Rarefaction. I. General Theory for Single Gas. Physica A, 203, 437-456.
  34. Sharipov, F. (1994) Onsager-Casimir Reciprocity Relations for Open Gaseous Systems at Arbitrary Rarefaction. II. Application of the Theory for Single Gas. Physica A, 203, 457-485.


The Spline Functions

We follow here the Hermite cubic splines functions as defined by Schultz [30]. In this way, we first consider the knots, defined in the interval as,


where, in this work, we use. To express a function, for, in terms of the spline functions, we write


where are constants to be determined and. In fact, associated with each one of the knots, there are two spline functions, which are defined differently for even or odd values of,


for. Here, we introduce



and we consider that the spline functions are zero unless otherwise, to write



for, and


Similarly, we write



for, and