Computational fluid dynamics basics

In short

Computational fluid dynamics (CFD) solves the governing equations of fluid motion numerically, by dividing a flow domain into millions of small cells and marching the equations to a converged or time-accurate solution. It is the third pillar of aerodynamic design alongside analytical theory and wind tunnel testing, and it dominates modern aerospace practice because it gives full-field data (pressure, velocity, temperature, species) at a fraction of tunnel cost. In Indian university curricula it appears as "CFD" or "Numerical Methods in Fluid Flow and Heat Transfer" (AICTE model curriculum, 6th/7th semester), and GATE Aerospace examines it through fluid mechanics, numerical methods and compressible flow. Employers that test CFD in interviews include HAL, DRDO labs (ADA, ADE, GTRE, DRDL), CSIR-NAL, ISRO/VSSC and LPSC, Tata Advanced Systems, Godrej Aerospace, and the large captive engineering centres (Airbus India, Boeing India, Safran, Collins Aerospace, GE Aerospace Bengaluru, Rolls-Royce India, Honeywell, Ansys and Altair India). Typical graduate roles are CFD analyst, aerodynamics engineer and methods/validation engineer.

Explained

Governing equations and their mathematical character

• CFD solves conservation of mass, momentum and energy for a continuum fluid. For a Newtonian fluid these are the Navier–Stokes equations (strictly, the continuity + momentum + energy set).
• Continuity: ∂ρ/∂t + ∇·(ρV) = 0. Momentum: ρ DV/Dt = −∇p + ∇·τ + ρg, with τ the viscous stress tensor. Energy: in terms of total enthalpy or internal energy, with conduction (Fourier) and viscous dissipation terms.
• Needs a closure: equation of state p = ρRT for a perfect gas, plus μ(T) from Sutherland's law and k = cpμ/Pr.
• Conservative (integral) form is preferred in CFD because it captures shocks correctly; non-conservative (differential) form is mathematically equivalent only for smooth flows.
• Simplified limits examined in exams: Euler equations (μ = 0, inviscid, still nonlinear and shock-capturing), potential flow (∇²φ = 0, incompressible irrotational), boundary layer equations (parabolised), Stokes flow (Re ≪ 1).
• Equation classification (governs the numerical method):
  • Elliptic — e.g. Laplace, incompressible pressure equation, subsonic steady potential flow: information travels in all directions, so boundary conditions are needed on the whole boundary; solved iteratively.
  • Parabolic — unsteady diffusion, boundary layer equations: marching in one direction (time or x).
  • Hyperbolic — unsteady Euler, steady supersonic potential flow: information travels along characteristics at finite speed; permits discontinuities (shocks).
• The transonic Euler/N-S problem is mixed elliptic–hyperbolic, which is why transonic CFD is hard.

Discretisation: turning PDEs into algebra

Three classical methods:

• Finite difference method (FDM): replaces derivatives by Taylor-series differences on a structured grid. Simple, high order easy, but awkward on complex geometry and not inherently conservative.
• Finite volume method (FVM): integrates the conservation law over each control volume and balances fluxes across faces. Locally and globally conservative, works on arbitrary unstructured cells — the method used by ANSYS Fluent, CFX, STAR-CCM+, OpenFOAM, SU2, and all industrial aerospace codes.
• Finite element method (FEM): weighted residual / weak form with shape functions; strong in structures, used in CFD mainly via stabilised formulations (SUPG/PSPG) and in FSI.
• Spectral / spectral-element and discontinuous Galerkin (DG) methods give very high order accuracy; used in research LES/DNS.
• Order of accuracy: truncation error ~ O(Δxⁿ). First-order upwind is very stable but highly numerically diffusive (smears shocks and shear layers); second-order central is accurate but oscillatory in convection-dominated flow.
• Upwinding respects the direction of propagation. Practical schemes: second-order upwind, QUICK, MUSCL, and TVD/limiter schemes (minmod, van Leer, Venkatakrishnan) that are second order in smooth regions but drop to first order at discontinuities to avoid overshoot — this is the Godunov barrier in practice.
• Flux schemes for compressible flow: central + artificial dissipation (Jameson–Schmidt–Turkel), flux-vector splitting (van Leer, AUSM, AUSM+), flux-difference splitting / approximate Riemann solvers (Roe, HLLC).
• Time integration: explicit (Euler, 4-stage Runge–Kutta — cheap per step, CFL-limited) versus implicit (backward Euler, dual-time, LU-SGS — large steps, needs linear solves). Steady runs use local time stepping, multigrid and residual smoothing to accelerate convergence.

Pressure–velocity coupling and linear solvers

• In incompressible flow there is no pressure equation; pressure is a Lagrange multiplier enforcing ∇·V = 0. Pressure is recovered by a projection/segregated approach: SIMPLE, SIMPLEC, SIMPLER, PISO (good for transient), or a coupled solver.
• Collocated grids need Rhie–Chow interpolation (or staggered grids / the MAC method) to avoid the checkerboard pressure mode.
• Compressible solvers use density-based formulations where pressure follows from the equation of state; low-Mach stiffness is cured by preconditioning.
• Resulting sparse systems are solved by Gauss–Seidel, line-relaxation, TDMA (Thomas algorithm) in 1-D, or Krylov methods (GMRES, BiCGStab) with ILU preconditioning; algebraic multigrid (AMG) is the industrial workhorse.
• Under-relaxation factors (typically 0.3 for pressure, 0.7 for momentum in SIMPLE) trade convergence rate against stability.

Mesh generation and mesh quality

• Structured (H, C, O topology; multi-block) — best accuracy per cell, aligned with flow, ideal for aerofoils and nozzles, but labour-intensive for complex geometry.
• Unstructured — tetrahedra, prisms, pyramids, polyhedra; automatic for complex airframes. Hybrid: prism/hex layers in the boundary layer + tets in the far field. Cartesian cut-cell / immersed boundary for very complex or moving bodies. Overset (Chimera) grids for store separation and rotors.
• Quality metrics: aspect ratio, skewness (< ~0.85 for tets), orthogonality at walls, cell volume growth ratio (< ~1.2 in the boundary layer), non-orthogonality, negative volumes.
• Boundary layer resolution is governed by y⁺: wall-resolved meshes need y⁺ ≈ 1 with 20–40 layers through the boundary layer; wall-function meshes target y⁺ ≈ 30–300 (never 5–30, the buffer layer, where both treatments are poor).
• Grid independence (convergence) study: run coarse → medium → fine with refinement ratio r ≈ 1.3–2, plot the target quantity (CL, CD, Nu, mass flow) and show it asymptotes. Quantify with Richardson extrapolation and the Grid Convergence Index (GCI) — the ASME V&V 20 / AIAA G-077 practice.
• Domain extent matters: far-field boundaries typically 20–50 chords away for an aerofoil, with vortex/far-field corrections for lift-induced blockage.

Turbulence modelling

• Turbulence is three-dimensional, unsteady and multi-scale; the ratio of largest to smallest (Kolmogorov) scales grows with Re, so direct resolution is unaffordable at flight Reynolds numbers.
• DNS — resolves all scales; cost ~ Re³ (roughly Re^9/4–Re³ depending on measure). Research only, low Re.
• LES — resolves large eddies, models the sub-grid scales (Smagorinsky, dynamic, WALE). Needs fine near-wall meshes; wall-modelled LES is the compromise.
• Hybrid RANS-LES — DES/DDES/IDDES, SAS; RANS in attached boundary layers, LES in separated regions. Standard for massively separated flow, cavity noise, jets.
• RANS — Reynolds averaging introduces the Reynolds stress −ρu'ᵢu'ⱼ, which must be modelled (the closure problem). Industrial default. For compressible flow, Favre (density-weighted) averaging is used.
• Boussinesq eddy-viscosity hypothesis: −ρu'ᵢu'ⱼ = μt(∂Uᵢ/∂xⱼ + ∂Uⱼ/∂xᵢ) − (2/3)ρkδᵢⱼ.
• RANS model families:
  • Algebraic/zero-equation: Baldwin–Lomax, Cebeci–Smith (legacy).
  • One-equation: Spalart–Allmaras (S-A) — solves for ν̃t; cheap, robust, the aerospace favourite for attached and mildly separated external flow.
  • Two-equation: standard k–ε (good for free shear flows, poor near walls and in adverse pressure gradients, needs wall functions); RNG and realizable k–ε; k–ω (good near walls, sensitive to freestream ω); Menter SST k–ω — blends k–ω near the wall with k–ε outside plus a shear-stress limiter, the best general-purpose choice for aerodynamics with separation.
  • Reynolds stress model (RSM): solves seven equations for the stresses directly; no isotropic-viscosity assumption, good for swirl and secondary flow, expensive and less robust.
• Transition modelling: γ–Re_θt (Langtry–Menter) correlation-based transition SST, or e^N / amplification-factor-transport methods — important for UAV and sailplane aerofoils at Re ~ 10⁵–10⁶ and for laminar-flow wing design.
• Known RANS weaknesses: shock-induced separation, strongly curved/rotating flow, buffet onset, separated maximum-lift prediction — the reason wind-tunnel and flight validation remain mandatory.

Boundary conditions, solution process and errors

Standard CFD workflow (exam answer):

1. Define the problem — objective, flow regime (Mach, Reynolds), 2-D/3-D, steady/unsteady, incompressible (M < 0.3) or compressible.
2. Geometry and cleanup — CAD import, defeature, close gaps, create the fluid domain by subtraction, use symmetry.
3. Mesh — surface mesh, boundary-layer prisms from a y⁺ target, volume mesh, check quality.
4. Physics setup — solver (pressure-based/density-based), turbulence model, fluid properties, reference values.
5. Boundary and initial conditions — see below.
6. Discretisation and solver controls — schemes, relaxation, Courant number, convergence criteria.
7. Solve and monitor — residuals plus integrated quantities (CL, CD, mass imbalance).
8. Post-process — contours, streamlines, Cp, Cf, y⁺, forces, surface integrals.
9. Verification and validation — grid study, compare with tunnel or flight data, state uncertainty.
10. Boundary condition types: velocity-inlet / mass-flow-inlet, pressure-far-field (characteristic-based, for compressible external flow), pressure-outlet (with back-pressure or supersonic extrapolation), wall (no-slip; adiabatic or fixed temperature or heat flux; moving/rotating), symmetry (zero normal velocity, zero normal gradients), periodic/translational-rotational (turbomachinery passages), and non-reflecting/characteristic BCs for acoustics.
11. Numerical BC discipline: total pressure/total temperature at a subsonic inlet; all variables specified at a supersonic inlet, none at a supersonic outlet; never place a boundary in a region of recirculation.
12. Convergence means (i) residuals dropped 3–5 orders, (ii) monitored integral quantities flat, and (iii) global mass/energy imbalance < ~0.1%. Residual drop alone is not convergence.
13. Error taxonomy (V&V vocabulary): round-off, iterative (convergence), discretisation (truncation) errors — these are *verification*; modelling error (turbulence, chemistry, geometry simplification) — this is *validation*. "Verification = solving the equations right; validation = solving the right equations."
14. Stability: explicit schemes must satisfy the CFL condition; the von Neumann method gives the stability limit of a linear scheme; false diffusion and dispersion are the two faces of truncation error.
15. Practice and tools in India: commercial — ANSYS Fluent/CFX/ICEM, Siemens STAR-CCM+, Cadence Fidelity, Hexagon Cradle; open source — OpenFOAM, SU2, Code_Saturne; Indian indigenous codes — VSSC's PARAS-3D, CSIR-NAL and VSSC in-house solvers, and Project PARAM/AIRAWAT and the National Supercomputing Mission machines (PARAM Siddhi-AI, Pratyush at IITM Pune and Mihir at NCMRWF Noida) on which large aerospace and weather runs are performed. There is no Indian Standard that certifies CFD results; aerospace work follows customer and certification guidance (CEMILAC/DGAQA for military airworthiness in India, DGCA/EASA/FAA CS-25/FAR-25 for civil), where CFD supports but does not replace test evidence.

Key points

• CFD — numerical solution of the mass, momentum and energy conservation equations over a discretised domain.
• Navier–Stokes — momentum equations for a Newtonian viscous fluid; Euler equations are their inviscid limit.
• Conservative form — integral flux-balance form; essential for correct shock capture.
• Equation type — elliptic (subsonic, needs BCs all round), parabolic (marching), hyperbolic (supersonic, characteristics); transonic is mixed.
• FVM — integrates fluxes over control volumes; conservative and geometry-flexible, hence the industry standard.
• Order of accuracy — first-order upwind is stable but numerically diffusive; second-order schemes are the working minimum for aerodynamics.
• Riemann solvers — Roe, HLLC, AUSM: upwind flux schemes that resolve shocks sharply.
• CFL number — Co = uΔt/Δx; explicit schemes need Co ≲ 1, implicit schemes can exceed it.
• SIMPLE/PISO — segregated pressure-correction algorithms for incompressible flow; Rhie–Chow interpolation prevents checkerboard pressure.
• y⁺ — non-dimensional wall distance; ≈1 for wall-resolved, 30–300 for wall functions, avoid 5–30.
• Closure problem — Reynolds averaging creates unknown Reynolds stresses, modelled via eddy viscosity (Boussinesq).
• SST k–ω — best general-purpose RANS model for aerodynamics with separation; Spalart–Allmaras is the cheap aerospace workhorse.
• DNS > LES > DES > RANS — decreasing cost and decreasing fidelity; DNS cost scales roughly as Re³.
• Grid independence — refine until the answer stops changing; quantify with Richardson extrapolation and GCI.
• Verification vs validation — solving the equations right vs solving the right equations.
• Convergence — residuals *and* integral quantities *and* mass balance, together.
• Limitations — turbulence and transition modelling, separated and buffet flows, mesh/user dependence, cost of unsteady runs.
• Tools — Fluent, STAR-CCM+, OpenFOAM, SU2; HPC under India's National Supercomputing Mission.

Formulas and examples

Example 1 — First cell height from a y⁺ target (classic interview/viva question). Aerofoil chord c = 0.5 m in air at U = 50 m/s, sea level: ρ = 1.225 kg/m³, μ = 1.789 × 10⁻⁵ Pa·s. Target y⁺ = 1 at x = c.

1. Re_c = ρUc/μ = 1.225 × 50 × 0.5 / 1.789 × 10⁻⁵ = 30.625 / 1.789 × 10⁻⁵ ≈ 1.712 × 10⁶.
2. Cf = 0.058 Re^−0.2: Re^0.2 = exp(0.2 × ln 1.712 × 10⁶) = exp(0.2 × 14.353) = exp(2.8706) ≈ 17.65 → Cf ≈ 0.058/17.65 ≈ 3.286 × 10⁻³.
3. τw = ½ρU²Cf = 0.5 × 1.225 × 2500 × 3.286 × 10⁻³ = 1531.25 × 3.286 × 10⁻³ ≈ 5.032 Pa.
4. u* = √(τw/ρ) = √(5.032/1.225) = √4.108 ≈ 2.027 m/s.
5. y = y⁺μ/(ρu*) = 1 × 1.789 × 10⁻⁵ / (1.225 × 2.027) = 1.789 × 10⁻⁵/2.483 ≈ 7.20 × 10⁻⁶ m ≈ 7.2 µm.
6. The first cell centre must sit at this height, so the first cell height is ≈ 1.44 × 10⁻⁵ m if the solver uses cell-centre distance. With growth ratio 1.15, reaching a ~5 mm boundary layer needs about n layers where 1.44e-5 × (1.15ⁿ − 1)/0.15 ≈ 0.005 → 1.15ⁿ ≈ 53 → n ≈ ln53/ln1.15 ≈ 28 layers.

Example 2 — Time step from a CFL limit. 2-D unsteady simulation of vortex shedding, smallest cell Δx = 2 mm, local velocity u = 60 m/s, explicit scheme requiring Co ≤ 1.

• Δt_max = Co·Δx/u = 1 × 0.002/60 = 3.33 × 10⁻⁵ s.
• For 20 shedding cycles at St = 0.2 behind a D = 0.05 m cylinder at U = 60 m/s: f = St·U/D = 0.2 × 60/0.05 = 240 Hz, period = 4.167 ms, total time = 20 × 4.167 ms = 83.3 ms.
• Steps needed = 0.0833 / 3.33 × 10⁻⁵ ≈ 2,500 steps — which is why an implicit scheme at Co ≈ 5 (with ~25 steps per cycle, i.e. ~1.7 × 10⁻⁴ s) is used instead.

Example 3 — Grid convergence index for drag. Three grids with refinement ratio r = 2 give CD: coarse f₃ = 0.0412, medium f₂ = 0.0381, fine f₁ = 0.0373.

1. Differences: f₂ − f₁ = 0.0008; f₃ − f₂ = 0.0031.
2. Observed order p = ln(0.0031/0.0008)/ln 2 = ln(3.875)/0.6931 = 1.3545/0.6931 ≈ 1.954 (≈ 2, as expected for a second-order scheme).
3. Richardson estimate: f_exact ≈ f₁ + (f₁ − f₂)/(r^p − 1) = 0.0373 + (−0.0008)/(2^1.954 − 1) = 0.0373 − 0.0008/2.871 = 0.0373 − 0.000279 ≈ 0.03702.
4. GCI(fine) = 1.25 × |(f₂ − f₁)/f₁| /(r^p − 1) × 100 = 1.25 × (0.0008/0.0373)/2.871 × 100 = 1.25 × 0.02145/2.871 × 100 ≈ 0.93%.
5. Interpretation: CD = 0.0370 with roughly ±1% numerical uncertainty — reportable; a GCI of 15% would not be.

Example 4 — Turbulence inlet quantities. Wind-tunnel-style inlet: U = 30 m/s, turbulence intensity I = 1% (0.01), turbulent length scale ℓ = 0.07 × D with D = 0.3 m → ℓ = 0.021 m.

1. k = 1.5 (U I)² = 1.5 × (30 × 0.01)² = 1.5 × 0.09 = 0.135 m²/s².
2. ε = Cμ^0.75 k^1.5/ℓ = 0.09^0.75 × 0.135^1.5 / 0.021. 0.09^0.75 ≈ 0.1643; 0.135^1.5 = 0.135 × √0.135 = 0.135 × 0.3674 ≈ 0.04960. → ε ≈ 0.1643 × 0.04960/0.021 ≈ 0.388 m²/s³.
3. ω = ε/(Cμ k) = 0.388/(0.09 × 0.135) = 0.388/0.01215 ≈ 31.9 s⁻¹.
4. μt/μ = ρCμk²/(εμ) = 1.225 × 0.09 × 0.0182/(0.388 × 1.789 × 10⁻⁵) ≈ 0.002008/6.94 × 10⁻⁶ ≈ 289 — a plausible eddy-viscosity ratio for a duct inlet.

Example 5 — Cost of DNS versus RANS (why RANS is used). Wing boundary layer at Re = 10⁷.

• DNS cell count scales roughly as Re^9/4 for a boundary layer: (10⁷)^2.25 = 10^15.75 ≈ 5.6 × 10¹⁵ cells — impossible (a large HPC run today handles 10⁹–10¹⁰ cells).
• A RANS mesh for the same wing is ~2 × 10⁷ cells, i.e. about 8 orders of magnitude cheaper. Hence RANS for design sweeps, DES/LES only for specific separated-flow questions.

Where it is used

• Airframe aerodynamics: cruise drag, high-lift (flap/slat) performance, store separation and transonic buffet at HAL, ADA (LCA Tejas and AMCA), CSIR-NAL (SARAS Mk2, HANSA-NG) and ADE (UAV/Rustom, Tapas).
• Propulsion: compressor and turbine blade rows, combustor flow and spray, nozzle and afterburner analysis at GTRE (Kaveri), GE Aerospace Bengaluru, Rolls-Royce India, Safran, Honeywell Aerospace Bengaluru.
• Launch vehicles and re-entry: base flows, plume–jet interaction, aero-thermal heating and TPS sizing for PSLV/GSLV/LVM3, RLV-TD and Gaganyaan crew module work at VSSC, LPSC and IPRC; hypersonic and scramjet CFD at DRDL.
• Rotorcraft and propellers: overset/sliding-mesh rotor CFD at HAL's Rotary Wing R&D Centre for ALH Dhruv, LCH Prachand and LUH.
• Services and product companies: CFD analyst roles at Ansys India, Altair, Siemens, Cyient, L&T Technology Services, Tata Elxsi, Tata Consultancy Services, Capgemini Engineering and Quest Global — supplying analysis to global OEMs.
• Adjacent sectors hiring the same skill: automotive external aero and underhood thermal (Maruti Suzuki, Tata Motors, Mahindra), HVAC and electronics cooling, wind energy (Suzlon), process and nuclear flows (BARC, NPCIL), and weather/ocean modelling at IITM Pune, NCMRWF and INCOIS.

Interview questions

• Q: What are the Navier–Stokes equations and why can't we usually solve them analytically? A: They are the conservation equations of mass, momentum and energy for a viscous Newtonian fluid; the convective term V·∇V makes them nonlinear and coupled, so closed-form solutions exist only for a handful of simplified cases.
• Q: Why does nearly every commercial CFD code use the finite volume method? A: Because it integrates fluxes over control volumes, so it conserves mass, momentum and energy exactly at the discrete level and works on arbitrary unstructured cells for complex geometry.
• Q: What is the CFL number and why does it matter? A: Co = uΔt/Δx compares the distance information travels in a time step with the cell size; explicit schemes become unstable above about Co = 1, while implicit schemes allow larger values at the cost of solving a linear system.
• Q: What is y⁺ and what value should you target? A: y⁺ is the non-dimensional wall distance based on friction velocity; aim for y⁺ ≈ 1 for wall-resolved low-Re models, or 30–300 for wall functions, and avoid 5–30 where neither treatment is valid.
• Q: Difference between RANS, LES and DNS? A: RANS models all turbulence after time-averaging; LES resolves large eddies and models sub-grid scales; DNS resolves everything — cost rises steeply in that order and only RANS is affordable for routine design at flight Reynolds numbers.
• Q: Which turbulence model would you pick for a wing near stall and why? A: SST k–ω (or Spalart–Allmaras with care), because SST behaves correctly near walls under adverse pressure gradients and limits eddy viscosity in separating shear layers; standard k–ε would delay separation.
• Q: What is the difference between verification and validation? A: Verification checks that the code and mesh solve the chosen equations correctly (grid study, order of accuracy); validation checks that those equations represent reality, by comparison with experiment.
• Q: How do you know a steady run has converged? A: Residuals have dropped three to five orders, the monitored integral quantities (CL, CD, mass flow) are flat over many iterations, and the global mass and energy imbalance is below about 0.1 per cent.
• Q: Why is first-order upwind discretisation discouraged for final results? A: Its truncation error acts like extra diffusion, smearing shocks, wakes and boundary layers, so it gives stable but artificially diffused and inaccurate answers.
• Q: What is numerical/false diffusion and how do you reduce it? A: Artificial smearing from low-order convection schemes and mesh misalignment with the flow; reduce it with second-order or limiter-based schemes, finer meshes and grids aligned to the flow direction.
• Q: Why does an incompressible solver need SIMPLE or PISO? A: There is no independent equation for pressure when density is constant, so pressure must be iterated as a correction that drives the velocity field to satisfy continuity.
• Q: CFD gives beautiful pictures that don't match the tunnel. What do you check first? A: In order: boundary conditions and reference values, mesh quality and y⁺, residual and integral convergence, discretisation order, then the turbulence/transition model — and confirm the experiment's own corrections (blockage, wall interference, support tares).
• Q: Can CFD replace wind tunnel and flight testing for certification? A: No — certification authorities accept CFD as supporting evidence within a validated envelope; test data remains the primary proof, because turbulence, transition and separation modelling still carry significant uncertainty.

Common mistakes

• Reporting a single-mesh result without any grid-independence study or uncertainty estimate.
• Treating a residual drop as convergence while forces or mass imbalance are still drifting.
• Using wall functions on a y⁺ ≈ 1 mesh (or a low-Re model on a y⁺ ≈ 100 mesh) — the two must match.
• Leaving first-order upwind on from the start-up phase and quoting those numbers as the answer.
• Placing an outlet boundary inside a recirculation zone, or a far-field boundary too close, and then blaming the turbulence model.
• Running a steady solver on an inherently unsteady flow (vortex shedding, deep stall) and averaging the oscillating residuals.
• Confusing the Reynolds-averaged mean field with instantaneous flow, or quoting Mach-number-independent "incompressible" results above M ≈ 0.3.
• Forgetting unit and reference-value consistency (reference area, length and dynamic pressure) so CL and CD come out scaled wrongly.

Related topics

Wind tunnel testing · CAD and 3D modelling · Compressible flow and shock waves · Boundary layer theory and viscous flow · Numerical methods and linear algebra · Aircraft design and aerodynamic optimisation