Improving Solution of Euler Equations by a Gas-Kinetic BGK Method

2009年2月西北工业大学学报Feb.2009第27卷第1期Journal of Northwestern Polytechnical UniversityVol 27 No. 1Improving Solution of Euler equationsby a gas-Kinetic BGK MethodLiu Ya, Gao Chao, F. LiuI National Key Laboratory of Aerodynamic Design and ResearchNorthwestern Polytechnical University, Xi' an 710072, China2. Department of Mechanical and Aerospace Engineering, University of Californiavine, CA 92697-3975, USAAbstract: Aim. The well known JST (Jameson-Schmidt-Turkel) scheme requires the use of adissipation term. We propose using gas-kinetic BGK(Bhatnagar-Gross-Krook)method, which isbased on the more fundamental Boltzmann equation, in order to obviate the use of dissipationterm and obtain, we believe, an improved solution. Section 1 deals essentially with three things:(1)as analytical solution of molecular probability density function at the cell interface has beenobtained by the Boltzmann equation with BGK model, we can compute the flux term byintegrating the density function in the phase space; eqs. (8)and(11)require careful attention;(2)the integrations can be expressed as the moments of Maxwellian distribution with differentlimits according to the analytical solution; eqs. (9)and (10)require careful attention;(3)thediscrete equation by finite volume method can be solved using the time marching methodComputations are performed by the BGK method for the Sods shock tube problem and a two-dimensional shock reflection problem. The results are compared with those of the conventionaljST scheme in Figs. 1 and 2. The BGK method provides better resolution of shock waves andother features of the flow fieldsKey words: Euler equations, Boltzmann equation, gas-kinetic BGK(Bhatnagar-Gross-Krook)scheme, JST(Jameson-Schmidt-Turkel)scheme, shock tube, shock reflectionwhere I and y are the Cartesian coordinates and1 Introduction and the Gas-KineticBGK Scheme2+ppvpv+ pThe Euler equations for two dimensionalpE风HPVH(2)inviscid flow can be written in integral form for aw is the macroscopic state vector, F and G are fluxregion n with boundary an asvectors in the z and y directions, respectively. pa wa +f (dy-Gdr)=0 1 U,V, P, E, and h denote the density, Cartesia中国煤化工Date received: 2007-10-23CNMHGAbout the first author: Liu Ya (1984-), Ph. D Candidate, Northwestern Polytechnical University. Research2西北工业大学学报第27卷velocity components, pressure, total energy, and upwind scheme, one uses Fi+/. =Fi. i.Howevertotal enthalpy, respectively. For a perfect gasfor a system of equations such as the EulerE方(Un+v)equations,the flux vector F contains signals( 3) propagating in both directions. In order to achieveH=E+the needed proper upwinding for the differentIn a finite-volume method, the discretization isvaves, Steger and Warming[a) and Van Leer[3]accomplished by dividing the field into small proposed the so-called flux-vector splittingcontrol volumes and then applying eq (1)to AVschemes. Godunovld, on the other hand, proposedOn a 2D Cartesian grid, this yieldsto evaluate F+1/. by solving a Riemann problemwith Fu. and Fi+l. as its left and right states.dr tatv- F-I/ Gt+1/- G j1/=0 Solution of the Riemann problem automatically(4) reflects the proper amount of upwinding for thewhere Wi,i is the cell-averaged value of w, Ar andvarious waves Inthe system. Since the eAy are the grid sizes in the r and y directions, equations are a nonlinear system, the exactrespectively. Fi+/a, is the Euler flux f at the solution of the Riemann problem is very time-proposed an approximatei+a,j cell interface, and Gi A+i/ is the euleRiemann solver to reduce the computational effortflux at the i,j+a cell interface, eq (4)mayThis type of schemes is usually referred to as theflux-difference splitting scheme.be solved by a time marching method for systemsAll of the abof ordinary differential equationshowever, calculate the flux term based on thea key element is the evaluation of the fluxEuler equations. It is known from gas kineticvectors Fu+1/, and Gi, H+i/. Conventional centraltheory, however, that the Euler equations can bedifference or upwind type schemes use central orderived from the Boltzmann equation for theupwind-biased interpolation schemes to calculatemotion of the gas molecules. Therefore, a moreFi+1. from its surrounding cell-center values. Forfundamental approach to evaluate the euler fluxexample, Jameson, Schmidt, and Turkel in theirvector is by solving the Boltzmann equation withJST scheme i] use the simple average: Fi+/. =0 Fi. and Fi+l. as its left and right equilibrium(F. +Fi+l ) Such a scheme is formally second states. Xut proposed such a method by usingorder. However, it is found to produce numerical BGK model for the Boltzamnn equation. Analyticaoscillations across shock waves and other flow solutions can be obtained to the BGK equation.discontinuities. In order to eliminate suchThe BGK model in the 2-D case isnumerical oscillations, it is necessary to add to eq妥++。影-(5)(4)a numerical dissipation term. The JST schemeuses an adaptive blended 2-nd and 4-th-orderwhere f and g are the molecular probabilityfinite-difference term of the flow variables todensity functions of space T, y, time t, particleprovide the needed numerical dissipation. The keyvelocity u, v and internal variable 4. Because ofthe mass, mto oscillation-free shock solution is to useand energy conservation inupwinding for the flux in the direction of wavethe process of particle collisions, the compatibilitycondition ispropagation, Adding numerical dissipation to acentraI-difference type scheme is equivalent to中国煤化工=1,2,3,4(6)upwinding. The coefficient of this term, however, wherCNMHGneeds to be tuned properly in order to obtainaccurate oscillation-free solution. In a first-order=1,,,(2+v2+)(7)Liu Ya et al; Improving Solution of Euler Equations by a Gas-Kinetic BGK MethodThe general solution of eq (5)at the cell interface results are at t=6. 1 ms. The initial conditions inan be expressed asthe present computation are the followingf+ln,=(1-e-")g+e-"f。(8)p=1.000,u=0,p=10,0≤x<0.5,Here for inviscid flow, we ignore all high order(12)spatial and temporal terms in the expansion of f P=0. 125,u=0, P=10,0.50,discontinthat of the JST(10)method. Decreasing the amount of numericalFrom eg. (5)we can get the relation between dissipation in the JST scheme causes spuriousthe macroscopic variables with microscoposcillation of the computed solutivariables,W=fp. dudvde F=ufD.dudvde00.8G=vfd,dudedFrom eqns(8)and (11)we can construct the fluxat the cell interface. Then the time marching stagein eq(4)can be performed.The above gas kinetic BGK method adds theappropriate amount of dissipation for the wave(b)pressureequations without having to tune certain numericalBGK acheneST acheneparameters needed in the conventional schemes. Itis also expected to provide more accurate androbust solutions for a wider range of flow0.5conditions. In the following section, we compare o0results by the gas kinetic BGK method to those bythe conventional JST scheme for two importanttest cases.c)velocity(a) temperature2 Results and DiscussionsFigure 1 Density, pressure, velocityand temperature distributionsA. Sod shock tube problem. The classic Sodin shock tu be at time t=6. 1 msshock tube problem presents an exact solution toH中国煤化工blique shockthe one-dimensional Euler equations containingwavewave reflsimultaneously a shock wave, a contact discontiCNMHnuity, and an expansion fan. The problem is fromin this test case.The computational domain is a rectangle of lengthcomputed in this test case using 500 cells. All the西北工业大学学报第27卷and height 1, with a uniform grid size of 60X2The boundary conditions produce an incidentThe boundary conditions are that of a reflecting shock with an angle of 29 in a free stream of Machsurface along the bottom surface, supersonicnumberM=2.9.outflow along the right surface, and prescribedFig 2 shows the computational results of thisxed values on the other two sides. At the inletcase as compared to the exact solution.The(P, U,, p)oo. n)=(1, 2.9,0,1/1. 4)(14) density contours on the left are obtained by theop surfaceJST scheme and those on the right are by the BGK(P,U,V,p)x,=(1.69997,2.61934method. Although the differences between the two0.50632,1.52819)(15) solutions are not obvious in the density contourplots, the density and pressure distributions showthat the BGK method provides betterDensity ContoursDensity Contoursresolution of the shock waves, demonstrating theof BGK Schemeof the BGK method over theconventional JST scheme.3 Conclusions8) JST Contours(b)BGK ContoursA computer code based on the gas kineticBGK method is developed for solving the twodimensional Euler equations. The gas kinetic BGKmethod solves the BGK model for the boltzmannequation to obtain the interface flux vectors in afinite-volume method for the eulerinstead of using numerical interpolation and addingnumerical dissipation terms in a conventionaldensityscheme Computations of the one-dimensionalshock tube problem and a two-dimensional shockFigure 2 Density contours and density and pressurereflection problem show that the gas kinetic BGKdistributions along y=0. 5 for the shockmethod provides better resolution of the shockreflection problemwaves and other flow discontinuities than theconventional Jameson, Schmidt, and Turkel (ST)schemeReferences:[1] Jameson A, Schmidt W, Turkel E. Numerical Solution of the Euler Equations by Finite Volume Methods Using Runge-Kutta Time-Stepping Schemes. AIAA-1981-1259[2] Steger J L, Warming R F. Flux Vector Splitting of the Inviscid Gasdynamic Equations with Application to Finite-Difference Methods. J Computational Physics, 1981, 40:263-293[3] Van Leer B. Flux Vector Splitting for the Euler Equations.Procs in Fluid DynamiSpringer Verlag, 1982, 507-51中国煤化工[4] Godunov SK. A Difference Scheme for Numerical ComputationCNMHMath Sb,1959,47:271~306[5] Roe P L. Approximate Riemann Solvers, Parameter Vectors, and Difference Schemes. J of Computational Physics第1期Liu Ya et al Improving Solution of euler Equations by a Gas-Kinetic BGK Method1981,43:357~372[6] Kun Xu. Gas-Kinetic Schemes for Unsteady Compressible Flow Simulations. Von Karman Institute for Fluid DynamicsLecture Series 1998-03, 29 Computational Fluid Dynamics, 1998[7] Michael A Saad. Compressible Flow. Prentice-Hall, 1985[8] Hong Luo, Joseph D. Baum, Rainald Lohner. A Hybrid Cartesian Grid and Gridless Method for Compressible FlowsJ of Computational Physics, 2006, 214: 618-6328AaA机加员碘员加a》》8》822》》2护护》》》《季文美文集》出版2008年9月西北工业大学出版社出版了《季文美文集》第一版,主编叶金福姜澄宇;编辑部主任徐澄成员唐宗焕、张庆恩、季钦伟、张光慎、祁随元、唐国靖自第003页的“季文美教授重要年表”摘录了如下3条(略有微小变动):(1)1912年1月5日生,2001年6月20日病逝;(2)1934年毕业于交通大学同年考取公费留学到意大利都灵大学学习航空工程,1936年获博士学位;(3)194-1945年兼任交通大学航空系主任;1952年后任华东航空学院、西安航空学院西北工业大学教授;1979~1981年西北工业大学副校长,1982~1984年西北工业大学校长,1984~2001年西北工业大学名誉校长。原国防科学技术工业委员会主任、中国航天科技集团总经理张庆伟2008年6月20日于上海为《季文美文集》作序。兹从序文摘录如下:西北工业大学已故老校长季文美先生,是我国著名的航空教育家和现代力学家在季先生诞辰96周年之际,学校决定编辑出版《季文美文集》并嘱我作序。作为学生和晚辈我既感高兴又恐难当此命但每忆昔日师恩风范,思付再三,写下一些粗浅的感受“1978年我考人西北工业大学飞机设计专业。1982年大学毕业后分配到航空工业部西安飞机设计研究所;3年后,再次回到母校攻读飞行器设计控制理论及应用方向的硕士学位。学习期间,季先生虽然没有直接授业于我,但他作为学校主要领导,所制定的教育教学指导思想所提出的教育教学改革措施,却使我们直接受益。这些年来,自己能在工作中取得一些成绩,与在校时打下的坚实基础和形成的‘三实作风是分不开“季先生从20世纪40年代初就执教于上海交大是当时教授会成员之一,许多今天在航空航天和其他领域卓有建树的著名专家、学者都曾是他的学生“在国际交往中,季先生既十分关注最新的学术动态,取他人之长为我所用,又表现出知识分子炽热的爱国情怀。1989年,为了维护国家的声誉争取原定的国际复合材料会议(ICCM)和第18届国际航空科学大会(ICAS'92)如期在我国召开,他不顾77岁的高龄费尽周折,四处奔走,赴深圳、英国和以色列,说服了外国友人;同时他以其对国际航空科学的突出贡献和人格魅力被国际航空科学理事会授予在国际上享有盛誉的莫里斯·鲁瓦奖。”胡沛泉2009年1月中国煤化工CNMHG