Recommendations of Future Work

In this thesis, we have demonstrated that the QDS-N² method is a very fast numerical method without being subject to convergence problem like conventional CFD methods. However, there are several areas need to be done in pushing the method forward, The areas outlined below should be examined in the future:

• To further reduce the computational time for large-scale multi-dimensional problems, the method should be implemented on multiple extension graphics processing units (GPU).

• To further reduce the numerical diffusion, a high-order stencil in calculating the conservative fluxes may be considered.

• To implement an adaptive mesh function in the region where large gradient of flow properties occurs.

• To extend the QDS-N² method for modelling the Navier-Stokes equation by employing the Chapmann-Engskog expansion theory to account for the viscous effect.

References

[1] Albright B.J. Don S. Lemons, Michael E. Jones, and Dan Winske., 2002, “Quiet direct simulation of Eulerian fluids,” Phys. Rev. E, Vol. 65, pp.1-4.

[2] Albright B.J, W. Daughton, D.S. Lemons, D. Winske, and M.E. Jones, 2002,

“Quiet direct simulation of plasmas,” Phys. Plasma, Vol. 9, Issue 5,pp. 1898.

[3] Anderson J. D.,1990, Modern Compressible Flow. McGraw-Hill, New York.

[4] Bazhenova T. V., Gvozdeva L. G., and Nettleton M. A., 1984, “Unsteady interactions of shock waves.” Prog. Aerospace Sci., Vol. 21, pp. 249–331.

[5] Bird G.A., 1994, Molecular Gas Dynamics and the Direct Simulation of Gas Flows, Clarendon Press, Oxford.

[6] Bird G.A., 1963,“Approach to translational equilibrium in a rigid sphere gas,”

Phys. Fluids, Vol. 6, pp. 1518-1519.

[7] Boyd I.D., Y. Jafry and J.V. Beukel, 1994, “Particle simulations of helium micro thruster flows,” J. Spacecraft Rockets, Vol. 31, Issue 2, pp. 271-277.

[8] Cave H.M., C.-W. Lim, M.C. Jermy, S.P. Krumdieck, M.R. Smith Y.-J. Lin, and J.-S. Wu, 2010, “Multi-Species Fluxes for the Parallel Quiet Direct Simulation (QDS) Method,” in: Proceedings of27^th International Symposium on Rarefied Gas Dynamics.

[9] Čada M., Torrilhon M., 2009,“Compact third-order limiter functions for finite volume methods,” J. Compt. Phys.,Vol.228, pp. 4118-4145.

[10] Coronell D.G. and K.F. Jensen, 1992, “Analysis of Transition Regime Flows in Low Pressure Chemical Vapor Deposition Reactors Using the Direct Simulation Monte Carlo Method,” J. Electro chemical. Soc., Vol. 8, pp. 2264-2273.

[11] Chen W.F., 2003, “Study of hybrid DSMC/EPSM method,” Chin. J. Comput.

[12] Cook G., 1998, “High accuracy capture of curved shock fronts using the method of space-time conservation element and solution element,” in:

Proceedings of 37^th AIAA Aerospace Sciences Meeting and Exhibit.

[13] Culbert B. Laney, 1998, “Computational Gasdynamics,” Cambridge University Press.

[14] Dyke M. van, 1997, An Album of Fluid Motion. Parabolic Press, Stanford.

[15] Gardiner C.W., 1985, “Handbook of Stochastic methods,” 2^nd ed., Springer, New York, pp.106-107.

[16] Gombosi T.I., 1994, Gaskinetic Theory, Cambridge University Press, Cambridge.

[17] Godlewski, Edwige, 1996, “Numerical approximation of hyperbolic systems of conservation laws,” Springer, New York.

[18] Huang J.-C., Lin Herng, Yang J.-Y., 2009, “Implicit preconditioned WENO scheme for steady viscous flow computation,” J. Comp. Phys., Vol. 228, pp.

420-438.

[19] Hirsch C., 1990, “Numerical Computation of Internal and External Flows, Vol.

2,” John Wiley & Sons, New York.

[20] Jameson A., W. Schmidt and E. Turkel, 1981, “Numerical solution of the Euler equations by finite volume methods using Runge-Kutta time-stepping schemes,”

the 14^th Aerospace Sciences Meeting, AIAA-81-12S9.

[21] Jameson A., T. J. Baker and N. P. Weatherill, 1986, “Calculation of inviscid transonic flow over a complete aircraft,” AIAA24th Aerospace Sciences Meeting, AIAA -86-0103.

[22] John D. Anderson Jr., 1995, Computational Fluid Dynamics.

[23] Jiang G., Shu C.-W., 1996, “Efficient implementation of weighted ENO schemes,” J. Comp. Phys., Vol. 126,pp. 202–228.

[24] Keats W.A., Lien F.S., 2004, “Two-dimensional anisotropic Cartesian mesh adaptation for the Euler Equations,” Int. J. Numer. Meth. Eng., Vol. 46, pp. 1099.

[25] Kröner D, 1997, “Numerical Schemes for Conservation Laws.” Wiley Teubner.

[26] Le Beau G. J. and F. E. Lumpkin III, 2001, “Application highlights of the DSMC Analysis Code (DAC) software for simulating rarefied flows,” Comp.

Meth. Appl. Mech. Eng. Vol. 191, Issue 6, pp. 595-609.

[27] Lim Chin Wai, Hadley M. Cave, Mark C. Jermy, Susan P. Krumdieck and Jong-Shinn Wu, 2010, “An Approximate Method for Solving Unsteady Transitional and Rarefied Flow Regimes in Pulsed Pressure Chemical Vapor Deposition Process using the Quiet Direct Simulation Method,” in: Proceedings of 27^thInternational Symposium on Rarefied Gas Dynamics, California.

[28] LeVeque R. J., 2004, “The dynamics of pressureless dust clouds and delta waves,” J. Hyperbolic Differential Equations, Vol 1, pp. 315-327.

[29] Karniadakis G.E. and A. Beskok, 2002, “Micro flows: fundamentals and simulation” Springer, New York.

[30] Kwon M.-K. and Y.-K. Hwang, 2006, “A Study on the Pumping Performance of a MolecularDrag Pump in the Rarefied Gas Flow,” in: Proceedings of 25^th International Symposium on Rarefied Gas Dynamics, St. Petersburg, Russia.

[31] MacCormack R. W. and B. S. Baldwin, 1975, “A numerical method for solving the Navier-Stokes equations with application to shock boundary layer interaction”, the13^th Aerospace Sciences Meeting, AIAA-7s-1.

[32] Macrossan M.N., 1989, “The Equilibrium Flux Method for the calculation of flows withnon-equilibrium chemical reactions,” J. Comp. Phys. Vol. 80, Issue 1, pp. 204-231.

[33] Macrossan M.N., 2001, “A particle only hybrid method for near continuum flows,” in: Proceedings of 22^nd International Symposium on Rarefied Gas Dynamics, Vol. 585, pp. 426.

[34] Majda A., 1984, “Compressible fluid flow and systems of conservation laws in several space variables,” Springer-Verlag, New York.

[35] Pullin D.I., 1980, “Direct simulation methods for compressible inviscid ideal gas flow,” J. Comp. Phys. Vol. 34, pp. 231-244.

[36] Salichs S. Serna, 2005, High Order Accurate Shock Capturing Schemes for Hyperbolic Conservation Laws based on a new Class of Limiters, Ph.D. Thesis, Doctoral Theses in Applied Mathematics, Universitat de Valencia.

[37] Schulz-Rinne C.W., Collins J.P., and H.M. Glaz, 1993, ‘Numerical Solution of the Riemann problem for two-dimensional gas dynamics,” SIAM J. Sci. Compt., Vol. 14, pp. 1394–1414.

[38] Serna S., 2006, “A class of extended limiters applied to piecewise hyperbolic methods,“ SIAM J. Sci. Compt., Vol. 28, pp. 123–140.

[39] Shu C.-W., Osher S., 1988, “Efficient implementation of essentially non-oscillatory shock-capturing schemes,” J. Comp. Phys., Vol. 77, pp. 439–471.

[40] Skews B.W., 1967, “The perturbed region behind a diffracting shock wave.” J.

Fluid Mech., Vol. 29, Issue 4, pp.705–719.

[41] Smith M.R., 2003, “Hybrid methods in near-continuum flows,” M.Phil. Thesis, University of Queensland, Australia.

[42] Smith M.R., M. N. Macrossan, M. M. Abdel-Jawad, 2008, “Effects of direction decoupling in flux calculation in finite volume solvers,” J. Comp. Phys. Vol. 227, pp. 4142-4161.

[43] Smith M.R., 2008, “The true direction equilibrium flux method and its application”, PhD Thesis, University of Queensland, Australia.

[44] Smith M.R., H.M. Cave , J.-S. Wu , M.C. Jermy , Y.-S. Chen, 2009, “An improved Quiet Direct Simulation method for Eulerian fluids using a second-order scheme,” J. Comp. Phys. Vol. 228, pp. 2213–2224.

[45] Smoller, Joel, 1983, “Shock waves and reaction-diffusion equations,” Springer-Verlag, New York.

[46] Sun M. and K. Takayama., 1997, “The formation of a secondary shock wave behind a shock wave diffracting at a convex corner.” Shock Waves, Vol. 7, pp.

287–295.

[47] Takayama K. and O. Inoue., 1991, “Shock wave diffraction over a 90 degree sharp corner.” Shock Waves, Vol. 1, pp. 301–312.

[48] Teshima K. and M. Usami, 2001, “DSMC simulation of axisymmetric supersonic free jets,” Comp. Fluid Dyne. J., Vol. 10, pp. 525-530.

[49] Toro E.F., 1999, Riemann Solvers and Numerical Methods for Fluid Dynamics, second ed., Springer, Berlin, Germany.

[50] Toro E. F., 1991, “A linearized Riemann solver for the time-dependent Euler equations of gas dynamics,” Proc. Roy. Soc. London Ser. A, 434, pp. 683-693.

[51] Toro E. F., 2009, “Riemann solvers and numerical methods for fluid dynamics:

a practical introduction.” Springer-Verlag, New York.

[52] Tutkun, B., Edis, F., O., 2010, “A GPU Application for High-Order Compact Finite Difference Scheme,” 22nd Int. Conference on Parallel CFD 2010, Kaohsiung, Taiwan.

[53] Van Leer B., 1977, "Towards the ultimate conservative difference scheme III.

Upstream-centered finite-difference schemes for ideal compressible flow", J.

Comp. Phys., Vol. 23, Issue 3, pp. 263–275.

[54] Visbal Miguel R., Gaitonde Datta V., 1999, “High-Order Accurate Methods for Complex Unsteady Subsonic Flows,” AIAA Journal, Vol. 37, pp. 10.

[55] Wagner W., 1992, “A convergence proof for Bird's direct simulation Monte Carlo method for the Boltzmann equation,” J. Stat. Phys., Vol. 66, pp. 1011-1044.

[56] Woodward P., Colella P., 1984, ”The numerical simulation of two-dimensional fluid flow with strong shock,” J. Comp. Phys., Vol. 54, pp.115.

[57] Wu M.G., 2003, “Study on hybrid EPSM/EPSM method for chemical reacting flow,” Chin. J. Comp. Phys., Vol. 20, Issue 3, pp. 564.

Appendix A

Flow Structures of the Shock Wave Diffraction

Keats [Keats et al., 2004] describes the theory of the shock wave diffraction that including experimental and computational result from Skews [1967], and summarized the secondary physical phenomena in the perturbed region behind the shock wave.

Skews performed experiments for a variety of Mach numbers and convex corner angles, and has outline the structure of the perturbed region; the structure is shown in Figure 3-5.1. Skews determined experimentally and tabulated the following correlations:

• The slipstream angle variation with the shock Mach number Ms.

• The terminator angle variation with Ms.

• The relationship between Ms and the velocity of the secondary shock.

• The contact surface velocity variation with Ms.

• The variation of the vortex angle and velocity with Ms. The flow structures can be described as follows:

• Incident shock: Diffracts in a similar way to a sound wave: its radius of curvature is approximately Wt.

• Reflected sound wave: Propagates upstream and marks the start of the curvature of the incident shock.

• Slipstream: Due to separation, it separates high-velocity gas on the upper side from almost stationary gas on the lower side. It represents the outermost characteristic of the Prandtl-Meyer expansion fan.

• Terminator: The first characteristic of the Prandtl-Meyer expansion; the single separating the terminator from the horizontal increases with rising Ms.

• Second shock: The region between the slipstream and the terminator is a uniform flow region parallel to the slipstream, and the second shock is a normal shock caused when the flow in this region exceed Mach 1.0 [Sun et al., 1997].

• Vortex: Located just below the slipstream, its location is well defined for Ms≤ 1.5. The angle between the vortex and the slipstream decreases as Ms increases.

• Contact surface: Originates at the intersection point of the reflected sound wave and the incident shock, but is highly diffuse in this region. It becomes better-defined as it nears the region containing the rest of the flow structures.

Tables

Table 2-1 The value of weight and abscissas for the Gaussian quadrature.

Number of QDS particles Weight (w_J) Abscissas (q_J)

2 1

2 2

± 1

2 π

3 0 2

3 π 1 6

±2 1

6 π

4 1 3 6

2 2

± − ^{4 3}

(

^{−π 6}

)

1 3 6

2 2

± + ^{4 3}

(

^{+π 6}

)

 

   

Table 3-1 Comparison of computational expenses for QDS schemes using 2N and N² dimensional reconstruction.

Number of cells

QDS solvers

2N N²

300x100 8.41 min 23.15 min

450x150 28.9 min 78.3 min

600x200 68.56 min 183.6 min

1000x500 478.4 min 1282.6 min

   

Table 3-2 QDS scheme time cost in Euler-4-shocks interaction case.

Number of cells QDS solvers

2N N²

1000x1000 13.29 hours 55.6 hours

100x100 45 (s) 189 (s)

200x200 375(s) 1520 (s)

300x300 1307 (s) 5199 (s)

         

   

Table 3-3 pre- and post-shock fluid initial conditions.

Mach number

1.5 2.4

ρ2/ρ

1 1.862 3.212

T₂/T₁ 1.32 2.04

Up 0.8215 1.956

   

Table 3-4 Parallel Performance for a 2D shock-bubble problem with 2.4 millioncomputational cells.

Number of processors Computation time (sec.) Number of cells

1 361.613 20000

4 348.995 80,000

9 428.812 180,000

16 410.294 320,000

25 425.229 500,000

36 444.135 720,000

49 472.599 980,000

 

Table 3-5 Parallel computation times for shock-bubble problem with 500,000 cells at 2000 time steps in simulation time 0.2.

Number of processors Computation time (s) Number of Cells for one processor

1 8220.19 500,000

2 4299.89 250,000

4 2761.14 125,000

8 1210.08 62,500

16 762.245 31,250

25 425.229 20,000

32 347.833 15,625

64 158.845 7,813

80 132.77 6,250

100 106.164 5,000

256 42.9935 1,953

   

   

 

Table 3-6 Parallel computation times for shock-bubble problem with 2 million cells in 2000 time steps.

Number of processors Computation time (sec.)

Number of cells for one processor

1 29740.59 2,000,000

4 8744.45 500,000

8 4978.56 250,000

16 2890.74 125,000

25 1827.92 80,000

32 1348.79 62,500

64 679.745 31,250

128 340.572 15,625

256 172.599 7,813

   

   

 

Table 3-7 Parallel computation times for shock-bubble problem with 12.5 million cells in 2000 time steps.

Number of processors Computation time (sec.)

Number of cells for one processor

1 187130.1 12,500,000

8 27532.8 1,562,500

16 19929.3 781,250

32 8613.3 390,625

64 4463.82 195,313

100 2806.29 125,000

128 2225.91 97,656

256 1115.02 48,828

   

Figures

Figure 2-1. Schematic showing the way fluxes of conserved quantities between source and destination cells are calculated using the “overlap” function in QDS [Smith et al.,

2009].

Destination cell

Source cell (Area, AS)

vyk∆t A

vxj∆t

Overlap area, A= vxjvyk∆t² Mass transported from source to destination cell is proportional to A/AS

       

  Figure 2-2.Flowchart describing QDS-2N particle computation with gradient inclusion.

                 

  Figure 2-3. The special reconstruction convention for current amount of conserved

quantity Q in one cell.

Figure 2-4. QDS flux procedure within a general (arbitrary) spatial reconstruction of conserved quantity Q.

Figure 2-5. Flowchart describing QDS particle computation with gradient inclusion.

Figure 2-6. Two-dimensional motion of a single QDS particle showing “sub-particle”

contributions.

Figure 2-7. Three-dimensional motion of a single QDS particle showing “sub-particle”

contributions.The green parallelogram is presented the concept in two-dimensional.

DESTINATION A

DESTINATION B

DESTINATION C

DESTINATION D

DESTINATION E

DESTINATION F

DESTINATION G

SOURCE

 

     

Figure 3-1. Waves generated in shock tube following rupture of diaphragm [Anderson, 1990].  

     

Contact surface moving at the velocity of the gas behind the shock when diaphragm ruptured.

Normal shock wave propagating

Expansion wave propagating

outflow outflow

Figure 3-2. The shock tube problem as computed by pre-QDS method with a uniform grid of 200 zones. The results were discussed thedifference to the QDS 1^st ~3^rdmethod

and Riemann solver using MINMOD limiter at time 0.1.

Figure 3-3. The interaction of two blast wave computed by the QDS method with 400 grids at t = 0.0038. The solid black line is WENO (fifth order) scheme with 10,000

grids.

Figure 3-4. Density profile of the shock-acoustic-wave case at t = 1.8. The solid black line is WENO-3 (fifth order) with 2000 grids compared with QDS method which

without limiter form 1^st order to 3^rdorder.

0.5

               

  Figure 3-5. The structure of shock bubble interaction.

 

(a) (b)

(c)

Figure 3-6. Zoom of shock-bubble Schlieren image with 1000×500 cells at time of 0.2.

QDS 2^nd order (a) 2N method with van Leer’s limiter, (b) N² method, and (c) 2^nd order TVD result presented in [Čada, 2009] using the same resolution.

 

Figure 3-7. Zoom of Schlieren image of shock bubble problem at time of 0.2; (a) QDS-N² method with 300×100 cells; (b) QDS-2N method with 300x100 cells; (c) QDS-2N

method with 450×150 cells; (d) QDS-2N scheme with 600×200 cells.

   

 

         

 

Figure 3-8. The initial conditions for the first problem of Euler-4 shocks interaction.

     

   

         

(a)             (b)

 

(c)             (d)

Figure 3-9. Zoom of density contour line of Euler-four-shocks problem. Comparing the second-order QDS-N² method (a) using 100×100 grids with MC limiter and 2N method using 100×100 grids (b) and 200×200 grids (c), 300×300 grids (d) with MC limiter at

time of 0.4.

QDS 2n 2nd 100x100 cells ,CFL=0.5 Density, MC limiter

QDS 2N 2nd 200x200 cells, CFL=0.5 Density, Mc limiter

QDS 2N 2nd 300x300 cells, CFl=0.5 Density, MC limiter

 

(a) (b)

(c)

Figure 3-10. Zoom of the density contour lines of Euler four shocks problem. (a) the 2^nd order TVD-MUSCL method taken from Čada [Čada et al., 2009] using 1000x1000 points, CFL=0.8. (b) The third-order QDS-N² method used 1000x1000 grids with MC limiter at time of 0.8. (c) The third-order QDS-2N method used 1000x1000 grids with

MC limiter.

             

 

Figure 3-11. The initial conditions for the second problem of Euler-four-shock interaction.

ρ =1.0, ux = 0.75 p=1.0, u_y = −0.5

ρ= 2.0, ux = 0.75

p=1.0, u_y = 0.5

ρ =1.0, u_x = −0.75 p=1.0, uy = 0.5

ρ = 3.0, u_x = −0.75 p=1.0, uy = −0.5

Figure 3-12. Density profile of the four contacts problem for second-order TVD-MUSCL method taken from [Čada et al., 2009].

 

  (a)  

  (b)  

Figure 3-13. Density contour obtained from QDS N² solver (a) and 2N solver (b) by using 1000×1000 cells, 2^nd order method with MINMOD limiter. The CFL number is

0.5. Level form 0 to 2.4 at 0.05 interval of line.

(a) (b)

(c) (d)

Figure 3-14. Density contour obtained from QDS-2N solver with 5 particles (a) and 9 particles in each direction (b);QDS-N²method with 5 particles (c) and 9 particles in each

direction (d) by using 1,000×1,000 cells, 2^nd order method with MINMOD limiterat time of 0.8. The CFL number is 0.5. Level form 0 to 2.1 at 0.05 interval of line.

(a)

(b)

Figure 3-15. Density contour obtained from QDS-N²solverusing (a) 2,000×2,000and (b) 3,000×3,000 cells, 2^nd order method with MINMOD limiterat time of 0.8. The CFL

number is 0.5. Level form 0 to 2.1 at 0.05 interval of line.

Figure 3-16. Geometry and boundary conditions for the Mach 3 flow over a forward facing step in a wind tunnel. All boundaries with exceptions to the inflow and outflow

are secularly reflective. The outflow boundary is calculated through interpolation of states of interior cells.

  2.4

0.6 0.2 1.0

Inflow boundary

Outflow boundary Reflect

Reflect

 

Figure 3-17. Contour of density at 4.0s for Mach 3 flow over a foeward facing step in a wind tunnel. Compare the 2^nd order QDS-2N method (top) and QDS-N² method (middle) for 600×200 grids. (buttom) The result of Keats and Lien [Keats et al., 2004]

X

Y

0 0.5 1 1.5 2 2.5 3

0 0.5 1

[Mch 3]t=4.0

QDS 2N 2nd 600x200 cells,CFL=0.1 Density, MINMOD limiter.

30 contours:0.2568 ~ 6.067

Figure 3-18. Structure of the perturbed region behind a diffracting shock wave, defined by from Skews [1967].

Figure 3-19. The output for compulsory figure for shock wave diffraction (by Takayama [1991]).

Figure 3-20. The initial geometry of the shock wave diffraction over degree sharp corner.

Wall

Inlet ^state2 M

_s

_state1

Outlet

Wall

Figure 3-21. Schematic of moving shock waves [Anderson, 1990].

(a)

(b) (c)

Figure 3-22. The density contours of the shock wave diffracting over 90 degree sharp corner with 400 × 400 grid, Ms=1.5. (a) the second-order TVD extension of Godunov method [Takayama et al., 1991]. (b) the second-order QDS-2N method and (c) the

second-order QDS-N² method with MC linter, CFL=0.5.

(a)

(b) (c)

Figure 3-23. Schlieren image of the shock wave diffracting over a 90 degree sharp corner, Ms=1.5.(a) the experimental result made form Ritzerfeldet al. [Takayama et al., 1991]. (b) second-order QDS-2N method and (c) QDS-N² method with 400× 400 cells,

MC limiter, CFL=0.5.

(a)

(b) (c)

Figure 3-24. Vorticity magnitude contours compared (a) exact solution and two result using 2^nd order (b) QDS 2N method and (c) QDS N² in 800×800 uniform cells. All

results are taken the CFL number to 0.1.

Figure 3-25. The vorticity profiles along the central line passing through the vortex. The comparison contained the exact solution (blue squeal-symbol line), the QDS-N²method

using160×160 cells (red line), 800×800 cells (black dash-dot line), and 2N method using 800×800 cells (purple long-dash line),1600×1600 cells (green doted line). Two

method s are computed in MC limiter and CFL=0.1 at time 8.0.

Figure 3-26. The three-dimensional geometry of the Mach 2 flow over a pillar.

(a)

(b)

Figure 3-27. The Density contour of the Mach 2 flow over a pillar obtained using the second-order QDS-N² method (a) in two-dimension with 200 × 200 cells; (b) in three-dimension with 200 × 200 × 100 cells. The CFL factor is 0.5 using MINMOD limiter.

Figure 3-28. The three-dimensional geometry of the Mach 2 flow over a square block.

(a)

(b)

Figure 3-29. The Density contour of the Mach 2 flow over a square block obtained using the second-order QDS-N² method with 200 × 200 × 100 cells (a) in x-y surface;

(b) in x-z surface. The CFL factor is 0.5 using MINMOD limiter.

Figure 3-30. Parallel Performance of a 2D shock-bubble interaction with 2.4 million computational Cells.

       

  Figure 3-31.Strong scaling performance in the QDS-N² method with 500,000 cells on

various massively parallel systems.

100%

75%

Strong scaling for case of shock bubble interaction Total number of cells:

500,000

Time step =2000 Simulation time: 0.2

       

  Figure 3-32.Strong scaling performance in the QDS-N² method with 2 million cells on

various massively parallel systems.

100%

68.5%

Strong scaling for case of shock bubble interaction Total number of cells:

2000,000

Time step =2000 Simulation time: 0.1

       

  Figure 3-33.Strong scaling performance in the QDS-N² method with 12.5 million cells

on various massively parallel systems.

65.5%

100%

Strong scaling for case of shock bubble interaction Total number of cells:

12500,000 Time step =2000 Simulation time: 0.04

(a)

(b)

Figure 4-1 The value of mass flux for the difference of 2N and N² method. (a) The case 1 with the gradient 1.0e^-5; (b) case 2 with the gradient 1.0e^-6.

(a)

(b)

Figure 4-2.The value of momentum flux in x-direction for the difference of 2N and N² method.(a) The case 1 with the gradient 1.0e^-5; (b) case 2 with the gradient 1.0e^-6.

(a)

(b)

Figure 4-3.The value of momentum flux in y-direction for the difference of 2N and N² method.(a) The case 1 with the gradient 1.0e^-5; (b) case 2 with the gradient 1.0e^-6.

(a)

(b)

Figure 4-4.The value of energy flux for the difference of 2N and N² method.(a) The case 1 with the gradient 1.0e^-5; (b) case 2 with the gradient 1.0e^-6.

(a)

(b)

Figure 4-5. The value of energy flux for the difference of 2N and N² method at the low density range.(a) The case 1 with the gradient 1.0e^-5; (b) case 2 with the gradient 1.0e^-6.

(a) (b)

(b) (d)

Figure 4-6. Contour profile of Shock-bubble interaction. QDS-N² 2^nd order method using 1700×500 cells with MC limiter at time of 0.2. (a) Density, (b) temperature, (c)

velocity in x-direction and (d) velocity in y-direction.

(a)

(b)

Figure 4-7Shock-bubble Schlieren image with 1700×500 cells at time of 0.2. QDS 2^nd order (a) 2N scheme with van Leer’s limiter, (b) N² scheme.

(a)

(b)

Figure 4-8. The density contour obtianed using (a) the 2N method; (b) the N² method with MC limiter, CFL=0.5, 1700×500 cells.

(a)

(b)

Figure 4-9. The temperature contour obtianed using (a) the 2N method; (b) the N² method with MC limiter, CFL=0.5, 1700×500 cells.

(a) (b)

(c) (d)

Figure 4-10. contour profile of N² method. (a)Density, (b) temperature, (c) velocity in x-direction and (d) velocity in y-direction.

Figure 4-11. Contour of density at 4.0s for Mach 3 flow over a foeward facing step in a

wind tunnel. Compare the 2^nd order QDS-2N method (top) and QDS-N² method (middle) for 600×200 grids. (buttom) The result of Keats and Lien [Keats et al., 2004]

 

Figure 4-12. Contour of temperature obtained using second-order QDS-2N (top) and QDS-N² method using 4 simulation particles for Mach 3 flow over a forward facing

step in a wind tunnel.

(a)

(b) (c)

Figure 4-13. Schlieren image of the shock wave diffracting over a 90 degree sharp corner, Ms=2.4.(a) The experimental result made form Ritzerfeld et al. [Dyke, 1997]. (b) the second-order QDS-2N method, and (c) QDS-N2 method with 1000×1000 cells, MC

limiter, CFL=0.2.

Figure 4-14. Density contours of Mach 2.4 shock diffraction using the second-order 2N method (top) and the N² method with 1000×1000 uniform grids. The CFL number uses

0.2 with MC limiter.

Figure 4-15. Temperature contours of Mach 2.4 shock diffraction using the second-order 2N method (top) and the N² method with 1000×1000 uniform grids. The CFL

number uses 0.2 with MC limiter.

Figure 4-16. Velocity contours of Mach 2.4 shock diffraction using the second-order 2N method (top) and the N² method with 1000×1000 uniform grids in x-direction. The CFL

number uses 0.2 with MC limiter.

List of Publications

Journals: (* corresponding author)

1. Y.-J. Lin, M. R. Smith, F.-A. Kuo, H. M. Cave, M. C. Jermy and J.-S. Wu*, “A True-direction Reconstruction to the Multi-dimensional Quiet Direct Simulation Method,” Computer Physics Communications, 2013 (in press).

http://dx.doi.org/10.1016/j.cpc.2013.05.007

2. Y.-J. Lin, M. R. Smith and J.-S. Wu*, “Theoretical Analysis ofFlux Deviation of Quiet Direct Simulation Method: 2N vs. N²,” Computer Physics

Communications, 2013 (under preparation).

3. Y.-J. Lin, M. R. Smith and J.-S. Wu*, “Parallel Implementation of Quiet Direct Simulation Monte Carlo Method,” Computers and Fluids, 2013 (under

preparation).

International Conference Papers: (* corresponding author)

1. Y.-J. Lin, M. Smith, Fang-An Kuo, H. Cave and J.-S. Wu*, “Two-dimensional quiet direct simulation (QDS) using GPUs,” Conference on Computational Physics, Kaohsiung, Taiwan, December 15-19, 2009.

2. Y.-J. Lin, M. R. Smith, H. M. Cave, J.-C. Huang and J.-S. Wu*,, “General Higher Order Extension to the Quiet Direct Simulation Method,” The 27th International Symposium on Rarefied Gas Dynamics, California, July 10-15, 2010.

3. H.M. Cave, C.-W. Lim, M.C. Jermy, S.P. Krumdieck, M.R. Smith, Y.-J. Lin and J.-S. Wu*, “Multi-Species Fluxes for the Parallel Quiet Direct Simulation

3. H.M. Cave, C.-W. Lim, M.C. Jermy, S.P. Krumdieck, M.R. Smith, Y.-J. Lin and J.-S. Wu*, “Multi-Species Fluxes for the Parallel Quiet Direct Simulation

在文檔中 應用真實流向概念重建通量法則於靜直接模擬法速解尤拉方程式的研究 (頁 83-158)