1 /35 . Past, present & future of CFD: a limited review by Brian Spalding

1. 1/35. Past, present & future of CFD: a limited review by Brian Spalding About the past: • focussing on contributions made by Imperial College Mechanical Engineering Department. About the present, observations about: • the FVM versus FEM contest, • computational-grid trends, • turbulence-model trends. About the future: • Population models of turbulence, • FVM for fluid~solid interactions, • The APParition, • Partially-parabolic applications.

2. 2/35. The Past: Before the digital age Airplanes appeared years before digital computers. Yet designers could even then predict their lift and drag. They used a combination of potential-flow theory with boundary-layer theory. This proceeded by iteration: 1. First source-sink distributions were sought which caused streamlines to fit the airplane, which led to 2. distributions of pressure over the surface. 3.They then used boundary-layer theory to calculate the ‘displacement thickness’ of the layer, i.e. the extent to which the airplane seemed bigger than first assumed. 4. Then they repeated steps 1, 2, 3; until convergence.

3. 3/35. The Past: Is the pre-digital method relevant today? Their boundary-layer theory was crude: • two-dimensional, • integral, • with assumed velocity profiles. Therefore wind-tunnel tests were needed in addition. But the principle was sound. And it still is; for computers remain too small to allow adequately fine elliptic grids, despite the use of many levels of sub-division.

4. 4/35. The Past: CFD at IC MED; How we stumbled into it Prior to 1965, IC Mech Eng still used integral-profile methods for 2D boundary-layer flows. Profiles were polynomials, with coefficients deduced from weighted-integral conservation equations. Then the Eureka-insight flash: piece-wise-linear profiles were more flexible; with integration over ‘pieces’; i.e. with unity weighting factors. We had invented (our own kind of) finite-volume CFD. The rest is history. A remark aside: The finite-element community still uses non-unity weighting factors; to their great disadvantage. One day they’ll learn. More about this later.

5. 5/35. The Past at IC MED: Features of the first computer program It first appeared in Patankar’s PhD Thesis of 1967, later published as a book. It simulated 2D parabolic flows, using axial distance and dimensionless stream function as co-ordinates, so minimising false diffusion. The grid width expanded and contracted to cover only the region of interest. So it was ‘self-adaptive’. It handled turbulence via Prandtl’s mixing-length model. Wall functions made their first appearance in it. It used the TDMA, without iteration, for cross-stream solution; and it ‘marched’ in the main-flow direction.

6. 6/35. The Past at IC MED: More about 2D parabolic computer programs A second computer program, GENMIX, applied the same method to more general 2D parabolic flows, e.g. wakes, plumes, wall jets for film-cooling, flames, etc. A major use was for systematic validation studies of the then-emerging two-equation turbulence models. It too was published as a book, with coding; and therefore used by non-IC researchers. However the IC group had already developed and published a stream-function~vorticity program for simulating 2D elliptic flows. We sought therefore a method to escape its restrictions viz. to two dimensions; and to uniform density.

7. 7/35. The Past at IC MED: SIVA, SIMPLE and 3D Harlow and Welch (1965) had published 3D methods for unsteady compressible flows; but our 2D methods handled practically more-important steady and incompressible ones. We wanted to continue. Our first success was with SIVA (= Simultaneous Variable Adjustment) ( Caretto et al 1971). But it worked point-by-point and converged slowly. So SIMPLE came into existence, purloining elements from predecessors, but surpassing them all (being later surpassed in its turn by SIMPLER, SIMPLEST, SIMPLEC, etc.) Later it was extended to two-phase, free-surface, magneto-hydrodynamics, and much more.

8. 8/35. The Past at IC MED: 3D parabolic and partially parabolic The first publication of SIMPLE (1972) was for 3D parabolic flows. This is fact seldom remembered. Eager to show SIMPLE’s elliptic capability, we too quickly exemplified it. The world followed, and scarcely noticed its parabolic capability. Who uses it nowadays? However, finding our computers too small (they still are; and may be forever), we later created the partially-parabolic method.This stores 3D only pressures, but velocities 2D, saving memory at the expense of time. There were numerous publication; but little world-wide following. Our fault, no doubt; but the world’s loss (I believe). More about this below.

9. 9/35. The Past at IC MED: More about partially-parabolic IC’s publications concerned flows in curved and coiled tubes, rotating ducts, 2D turbines, and around ship’s hulls. Too few! And with wrong emphasis. Their authors (I was one) compared partially-parabolic with fully-parabolic (and therefore incorrect) solutions, rather than with more accurate fully-elliptic ones. Another distraction was: did the turbulence models used fit the experiments? Beside the point. What should have been stressed was: Replacing fully-elliptic by partially-parabolic, 1. scarcely affected accuracy or computer time, but 2. greatly reduced computer memory requirement. When authors forget the point, most readers will miss it.

10. 10/35. The Past at IC MED: starting the CFD software industry As engineers our aim was to be useful. So we offered our services to industry; with success. Imperial College was ill-suiteded to industrial contracts; so CHAM Ltd was founded as a pioneering ‘spin-off’. At first, each task was treated as a ‘start from scratch’; and only CHAM personnel could use the software. Soon came another ‘Eureka’ moment: Why not create a general-purpose package, with a closed-off core and open-to-users outside? Then the customers’ own staff could use it. Hence: PHOENICS (1981); followed by many emulators: FLUENT, Star-CD, Flow-3D, CFX, etc., etc.

11. 11/35. Present: some trends: FVM versus FEM Starting with FIDAP, finite-element-based codes appeared in the CFD-software market; and multplied. The burgeoning FEM literature implied that differential equations not multiplied by non-unity weighting functions (NUWFs) could’t be solved. Else why do it? Yet Finite-Volume codes use unity weighting factors (UWFs), i.e. they use no weighting at all. Those who never understood the FEM literature can now take heart: FEM-based CFD codes are no more. The UWFists have at last prevailed …. … for CFD.

12. 12/35. Present trends: FVM for solid stress? FEM still dominates the solid-stress field; but FVM can solve the problems just as well, calculating displacements in place of velocities. So a single FVM computer code can solve fluid-solid-interaction problems, e.g. thermal and mechanical stresses in a gas-turbine blade. Slides extracted from an earlier lecture will now illustrate this. The computer code used is PHOENICS which has a built-in SSFT, i.e.simultaneous-solid-fluid-thermal capability.

13. 13/35. A typical SSFT problem: blade in hot gas, cooledinternally Hot gas flows outside an internally-cooled blade-like solid. The un-structured grid which is used is shown below. The picture above shows the whole calculation domain, with gas inlet on the left and outlet on the right. Also visible is the central tube, which introduces the cooling air. The problem is illustrative, with idealised geometry. The smallest cells are placed near the curved solid-fluid interfaces.

14. 14/35. FVM solution for velocities Velocity vectors in the gas stream. Red is fast and green slow.

15. 15/35. Displacement and thermal strain in solid Displacement vectors computed at same time and in the same (SIMPLE) way as velocities. Thermal-strain distributions here shown as contours are also computed simultaneously. Repetition for sake of emphasis: The simulation is performed by a singleFVM-based code as part of a single calculation, with full compatibility between conditions in solid and fluid? Can any FEM code do the same?

16. 16/35. Pressure distribution in gas Pressure contours in the flowing gas. Red is high; blue is low. Pressure is computed by SIMPLE for gas only.

17. 17/35. Thermally and mechanically- induced stresses X-, y- and z-direction thermally-induced-stress contours within the ‘blade’. Red is compressive, blue – tensile. Notetheir strongly three-dimensional variation.

18. 18/35. Summary of experiences with FVM applied to solid-stress problems Many comparisons with both analytical and finite-element solutions have been made. They confirm that FVM for stress-in-solids problems is practicable, accurate and economical; it is at least as good as FEM. This is a fertile field for research, still almost explored. SIMPLE works well for both fluid flow and solid stress; but surely better SSFT-specific algorithms can be found. Other questions remaining to be answered concern: relative advantages of staggered and collocated structured grids; and of (various kinds of) unstructured grids. Extensions are also required to time-dependent phenomena: • to large (enough to influence the flow) displacements; and • to non-linear and plastic deformations.

19. 19/35. Present trends: grids for arbitrary body shape First used were body-fitted-coordinate grids which were topologically Cartesian, i.e. still ‘structured’; but sometimes hard to create. Therefore unstructured grids with tetrahedral cells (copied from FEM) were popular for many years. Polyhedral cells followed. These too present creation difficulties; and the current trend is back to Cartesian, albeit sub-divided as in just-shown solid-stress example. Then surface curvature may be allowed for by use of the ‘Immersed-Boundary Method’ (IBM).

20. 20/35. Present trends: early IBM examples The PHOENICS IBM viz. PARSOL, simulates, on the right, flow through a louvred wall. It looks realistic. But quantitative accuracy is improbable. The same is true of the football stadium on the left. Plausibility is easy to get. Reliable accuracy - much harder. The.

21. 21/35. Present trends: SPARSOL (Structured PARSOL) Some versions of the IBM perform poorly for solids which are thin compared with grid cells; careful calculation of the intersection locations is necessary. However realism can indeed be procured, with sufficient care (see below).

22. 22/35. Present grid trends: Various storage locations 1. Staggered: pressures and scalars at cell centres; velocities on cell boundaries. This is the natural choice. 2. Collocated: all variables at cell centres. Sometimes (unwisely?) preferred. 3.Other: Xcell (various). Seldom used, but having merit. In grid on right, scalars , e.g. temperatures are stored at triangle centroids. So they are 4 times (8 in 3D) as numerous as pressures and velocities.

23. 23/35. Present grid trends: How Xcell reduces false diffusion Blue fluid flows in from left, red from below. Grid is Cartesian staggered. Interface is blurred, with diagonally- directed flow, as seen on right. This is well known false diffusion due to upwind differencing. Less well known is that upwind differencing with Xcell grid causes no blurring at all for flow angles at 0, 90 or 45 degrees, as shown on right. False diffusion does exist at other angles, but less than without Xcell.

24. 24/35. Present grid trends: More-advanced Xcell Cartesian sub-divided grids can also be ‘triangularised’: see right. Not all cells need to be ‘triangularised’; only those where scalar gradients are large. In another version of Xcell, velocities are also stored at triangle centroids. This is ‘semi-collocated’ Xcell. Because pressures and velocities are not stored at the same points, it is free from the ‘checker-boarding’ ailment of conventionally collocated grids. What might be termed a ‘smart-grid’ technology is emerging which is ‘solution adaptive’.

25. 25/35. Present: turbulence-model trends in terms of how many population members Variables of popular turbulence models, e.g. k-e, LES, are local averages; implied population has 1 member. Eddy-break-up model for combustion (1971) is most-used 2-member example. Population theory is not new. Only 2-(or more)-member-population models can represent chemical reaction, swirling-flows, and un-mixing. The next three slides concern an un-mixing experiment, first performed in 1978, which no 1-member (i.e. conventional) model has ever been able to simulate. Will any one accept the challenge?

26. 26/35. Turbulence trends: The Stafford experiment Fill the lower half of a glass-sided vessel with coloured salty water, and the top half with clear fresh water. Connect electrodes to a battery. The salty water heats more rapidly than the fresh. The consequent Rayleigh-Taylor instability causes mixing. Soon the vessel appears to be filled with coloured fluid. Quickly switch off the current; then the two fluids start to un-mix! In the end, the original sharp interface is restored.

27. 27/35. Turbulence trends: Mixing followed by unmixing (Sapozhnikov and Mitiakov, 2010)

28. 28/35. Turbulence trends: A 2-member population model can do it Each member has its own vertical- direction velocity. One is +ve the other –ve. Values are calculated from Navier Stokes. At the start (on the left), the volume fraction is unity in the bottom half and zero in the top half. Later (in the middle) fragments of salty fluid rise, and even begin to concentrate at the upper surface. Later still (on the right), the heating has stopped; so the salty fragments, lose heat to the fresh water and fall down to the bottom again. Just as the video showed. Two-member models can simulate both mixing and un-mixing.

29. 29/35. The future – perhaps. Models of turbulence will use multi-member populations One-member models in effect represent temperature by one unity ordinate at calculated abscissa. They know nothing about PDFs. Multi-member models: 1. focus on several calculated ordinates at arbitrary abscissae. 2. simulate inter-member-physics; 3. calculate PDFs; 4. solve more equations; 5. generate much more information. Computation is cheap: ignorance is expensive. So surely multi-member models must become the norm.

30. 30/35. The future – perhaps Use one (FVM) code for both solid-stress and fluid/heat-flow problems? Common sense says ‘Yes!’; for the world-wide cost of FEM-for-solid-but-FVM-for-fluid is enormous; and all because FEM carried needless pre-computer baggage (the NUWFs) into the computer age: and others have been gulled into using it. The picture answers. And why?

31. 31/35. The future – perhaps Will general-purpose CFD codes survive? My answer? Yes, but underground. Apps will dominate. Apps, aka SimScenes, apply CFD to special classes of equipment, i.e. Simulation Scenarios,via app-specific menus and buttons. App users may know no more about CFD than apple eaters about arboriculture. Apps and apples can be equally healthy if the tree-roots are well nourished by the underlying CFD code.

32. 32/35. The future – perhaps: Revival ot the partially-parabolic method? The revived method would solve the simple potential-flow elliptic equation outside boundary layer, wake and jet. Inside each of these it would solve 3D parabolic Navier-Stokes equations on as fine a grid as needed (easy because only 2D storage is required). Elliptic and parabolic solutions would alternate, exchanging domain-boundary information each time. Why should this not work? Is it perhaps already used?

33. 33/35. The future – perhaps: Partially parabolic for automobiles? Early (1988) PHOENICS needlessly solved elliptic Navier-Stokes far from the vehicle surface where the flow is inviscid. Near much of surface the flow is 3D parabolic. But elliptic Navier-Stokes must be used for the wake. And behind wheels and wing mirrors.

34. 34/35. The future – perhaps: Implementation of partially- parabolic for automobiles What is needed, in order to implement such a hybrid solution procedure, is: • a CFD code with 2D and 3D, elliptic and parabolic capabilities (Note that PHOENICS is an acronym for Parabolic Hyperbolic Or Elliptic Numerical Integration Code Series; so it will do); • a flexible module for grid-to-grid transfer and interpolation of solved-for variables (PHOENICS has an embryonic one); and • a user-friendly module for problem set-up and run-cycle control (coming soon). Nothing of significant difficulty is involved.

35. 35/35. The future – perhaps: Terrestrial applications of partially-parabolic Urban-air-flow and wind-farm simulations have a predominant flow direction. Analysis can be parabolic over most of volume, with embedded elliptic sub-domains of recirculation. Grid fineness can be varied according to accuracy needs of each region. This and many other possible extensions of the partially-parabolic method promises a highly profitable future. Let’s help bring that about. The End