\makeatletter
\@ifundefined{HCode}                  
{\PassOptionsToPackage{russian,english}{babel}
\documentclass[screen, CRMECA, Unicode, Thematic,published]{cedram}
\newenvironment{noXML}{}{}%
\def\tsup#1{$^{{#1}}$}
\def\tsub#1{$_{{#1}}$}
\def\sfrac#1#2{{#1}/{#2}}
\def\thead{\noalign{\relax}\hline}
\def\endthead{\noalign{\relax}\hline}
\def\tbody{\noalign{\relax}\hline}
\def\xtbody{\noalign{\relax}}
\def\tabnote#1{\vskip4pt\parbox{.51\linewidth}{#1}}
\def\xxtabnote#1{\vskip4pt\parbox{.65\linewidth}{#1}}
\def\mathbi#1{\boldsymbol{#1}}
\def\xmathop#1{\mathop{\text{#1}}}
\def\rmDelta{\Delta}
\usepackage{etoolbox}
\usepackage{upgreek}
\usepackage{amsfonts}
\usepackage{amssymb} 
\csdef{Seqnsplit}{\\}
\def\botline{\\\hline}
\newcounter{runlevel}
\usepackage{multirow} 
\def\nrow#1{\@tempcnta #1\relax%
\advance\@tempcnta by 1\relax%
\xdef\lenrow{\the\@tempcnta}}
\def\morerows#1#2{\nrow{#1}\multirow{\lenrow}{*}{#2}}
\let\MakeYrStrItalic\relax
\def\jobid{crmeca20220900}
%\graphicspath{{/tmp/\jobid_figs/web/}} 
\graphicspath{{./figures/}} 
\def\rmmu{\upmu}
\def\mn{\phantom{$-$}}
\def\refinput#1{}
\def\back#1{}
\newcommand{\bnabla}{\boldsymbol{\nabla}}
\newcommand{\bgamma}{\boldsymbol{\gamma}}
\newcommand{\bdelta}{\boldsymbol{\delta}}
\newenvironment{inftab}{}{}
\def\Custref#1#2{\hyperlink{#2}{#1}}
\def\og{\guillemotleft}
\def\fg{\guillemotright}
\makeatletter
\g@addto@macro{\UrlBreaks}{\UrlOrds}
\gappto{\UrlBreaks}{\UrlOrds}
\makeatother
\DOI{10.5802/crmeca.240}
\datereceived{2022-09-16}
\daterevised{2024-01-17}
\dateaccepted{2024-01-18}
\ItHasTeXPublished
\def\oneunient#1{\cyrins{#1}}
\def\twounient#1{\cyrins{#1}}
\renewcommand\cyr
{
\renewcommand\rmdefault{wncyr}
\renewcommand\sfdefault{wncyss}
\renewcommand\encodingdefault{OT2}
\normalfont
\selectfont
}
\DeclareRobustCommand{\cyrins}[1]{%
  \begingroup\fontfamily{erewhon-TLF}%
  \foreignlanguage{russian}{#1}%
  \endgroup
}
}
{\documentclass[crmeca]{article}
\usepackage[T1]{fontenc} 
\usepackage[T2A,T1]{fontenc}
\usepackage{upgreek}
\def\CDRdoi{10.5802/crmeca.240}
\def\newline{\unskip\break}
\def\sfrac#1#2{{#1}/{#2}}
\def\xxtabnote#1{\tabnote{#1}}
\def\xmathop#1{\mathop{{#1}}}
\def\xtbody{\tbody}
\def\hypertarget#1#2{\label{#1}}
\let\beqimg\relax
\let\eeqimg\relax
\def\cyrins{}
\def\oneunient#1{\unient{0428}}
\def\twounient#1{\unient{042E}}
}
\makeatother

\usepackage{moreverb}

\dateposted{2024-06-13}
\datepublished{2024-11-05}
\begin{document}

\begin{noXML}

\makeatletter
\def\TITREspecial{\relax}
\def\cdr@specialtitle@english{More than a half century of Computational Fluid Dynamics}
\def\cdr@specialtitle@french{Plus d'un demi-si\`ecle de m\'ecanique des fluides num\'erique}
\makeatother

\CDRsetmeta{articletype}{history-of-sciences}

\title{Evolution of CFD as an engineering science. A~personal
perspective with emphasis on the finite volume method}

\author{\firstname{Akshai Kumar} \lastname{Runchal}\CDRorcid{0000-0002-8226-5234}} 
\address{Analytic \& Computational Research, Inc., 1931 Stradella Road,
Los \'{A}ngeles, CA~90077, USA}
\email{runchal@ACRiCFD.com}
\email{runchal@gmail.com}
\urladdr{www.ACRiCFD.com}

\begin{abstract}
Computational Fluid Dynamics---CFD for short---is a comparatively
recent development that has made a significant impact in engineering
sciences. The foundations of CFD were laid by developments in physics,
numerical analysis and matrix theory over the last 200 years or so. It
was, however, the electronic computer in the middle of the 20{th}
century that led to its birth and its widespread use. It has now become
an invaluable tool for almost every sphere of human activity. A number
of researchers contributed to the tools and technology that came to be
called CFD but the prime movers were two brilliant scientists---one at
Los Alamos National Laboratory and the other at Imperial College,
London. The author was a part of the Imperial College group and this is
a personal perspective on the evolution of CFD.

Over the last few decades, CFD has made a significant impact in a wide
spectrum of engineering and environmental sectors. These include both
traditional disciplines such as aerospace, automotive, mechanical and
thermal sciences, and newer disciplines and emerging fields of human
and social relevance such as biomedical, sports, entertainment, food
processing, fire safety, HVAC and energy efficiency. However, a number
of challenges remain. For many important applications, we lack both the
physics and the mathematical tools to adequately understand the
behavior of fluids. The current generation of CFD tools are resource
intensive and difficult to use that require a long period of training.
Most importantly, except for a very limited set of applications, one
cannot rely on the predictions with a high level of confidence.

Computational Fluid Dynamics is now poised on the threshold of a
revolution. Recent developments in machine learning, AI, computer
hardware, big data, IOT and Virtual Reality tools will lead to design
tools with CFD engines that are robust, reliable and easily accessible
to a practicing engineer rather than be the domain of a CFD expert. CFD
will be a tool that is ubiquitous by its absence---buried inside
devices and applications into diverse areas of human activity. This
would make it possible for a practicing engineer to effectively
visualize a fluid system, interact with its components, conduct CFD
simulation, and control its behaviour through dynamic intervention.
\end{abstract}

\keywords{\kwd{CFD}\kwd{Finite volume methods}\kwd{Machine
learning}\kwd{PINN}\kwd{AI}}

\maketitle

\end{noXML}

\section{Background}\label{sec1}

We inhabit fluids, and fluids inhabit us. We're immersed in air and two
thirds of the earth surface that we live on is water. Two thirds of the
human body is fluid. The solid parts of the earth and that of the human
body are actually highly deformable and often behave like fluids.

Computational Fluid Dynamics---CFD for short---predicts the behaviour
of fluids, plastics and deforming solids. CFD, in a broader sense, is
the science of ``\textit{Things That Move}''. To understand the
physical, biological and environmental world that we live in, we need
CFD. Every human endeavour involves fluids and CFD applications span
human existence: be it, industry, environment, human body, or the
biosphere. CFD tools are critical to aerospace, energy transport,
national security, engineering, design, food and chemical processing,
consumer electronics, environment, disaster management, life sciences,
earth sciences, astronomy, virtual reality, animation and increasingly
entertainment. Over the years it has become more a question of what CFD
can't do, rather than that of what it can do.

This is a rather subjective and personal perspective of the history of
CFD. In that context, my focus is that of a practicing engineer.
Therefore, I have emphasised the physics and engineering aspects,
rather than the significant developments that have occurred, in
parallel, in the fields of mathematics and numerical analysis. I grew
up in the remote hills of Himalayas, and never dreamt of participating
in the birth of a new branch of engineering science. But I was lucky
enough to witness the emergence of CFD at Imperial College. As a
footnote, McLeod Ganj, my hometown, became the home to the Dalai Lama
and a symbol of rejuvenation of Buddhism. It then played an important
role in the human rights movement, and became a prominent place on the
political maps of India, China and the US.

The mathematical foundations of fluid dynamics were laid by the 18{th}
Century Swiss mathematician, Leonhard Euler---one of the most prolific
mathematicians in history. He developed the Euler Equations that
describe the flow of an ideal frictionless fluid in the early 1750's
but formally presented to the Swiss Academy later in 1757
\cite{Euler1757}. Though this was a monumental theoretical
achievement, most flows of practical interest stayed outside the
purview of these equations since they ignored friction that plays a
very important role in the behaviour of real fluids. In 1821, the
French engineer Claude-Louis~\cite{Navier1822} introduced the concept
of viscosity to generalize the Euler equations to deal with the flow of
real fluids. Though Navier realized that there were forces between the
molecules of fluids, he had no formal concept of shear stress or its
relation to the state of the fluid. The next major breakthrough was in
the middle of the 19{th} century by the British physicist and
mathematician, Sir George Gabriel~\cite{Stokes1845}, who derived a
formal mathematical expression that related the stress due to friction
(viscosity) to the strain of the fluid. These developments led to the
Navier--Stokes Equations that put fluid dynamics on a firm and formal
mathematical basis. One unexpected result of these theoretical advances
was that whereas the Euler equations were, in many cases, somewhat
amenable to analytic solutions, the Navier Stokes equations were
intractable except for trivial and highly idealized one or
two-dimensional flows in simple geometrical configurations. 

The Navier--Stokes equations are now recognized to be one of the most
intractable mathematical problems in physics. Despite their broad
applicability, no proof is available that a solution always exists in
three dimensions and that, if it does exist, it is smooth and unique.
The Clay Mathematics Institute \cite{CI2017} has called this
challenge one of the seven most important open problems in mathematics
and has offered a prize of US\$1 million for a solution or a
counterexample.

In the absence of analytic solutions, for the rest of the 19{th}
century, the only option available was to focus on certain flows where
some terms in the Navier--Stokes equations could be ignored. This then
allowed either exact or approximate solutions to the simplified
equations. One early example was the case of creeping flows where the
frictional force is much greater than the \mbox{convective} inertial force.
Stokes~\cite{Stokes1851} formulated equations where the convective terms are
completely ignored and the simplified equations are amenable to exact
or approximate solutions. Another example is that of the Shallow Water
Equations where it is assumed that the horizontal extent of a body is
much larger than its vertical depth. This is the case, for example, for
large water bodies or the atmosphere at the meso- or global-scales.
Order of magnitude analysis then shows that the pressure at any point
is hydrostatic in the vertical direction. The vertical acceleration can
be dropped and the resulting simplified equations allow solutions for
many problems of{\unskip\break} interest.

The next major advance came at the start of the 20{th} century when
Prandtl~\cite{Prandtl1904} introduced the concept of a boundary layer. Many
practical flows consist of a dominant direction of flow. Flow in a pipe
or that on a flat plate, a jet in ambient surrounding, or the wake
behind a cylinder, are common examples. For such flows, he visualized
the flow field as consisting of two regions. A boundary layer region
where flow gradients are dominant and friction plays an important role,
and a far-field region where gradients are negligible and friction can
be ignored. He then derived a simplified form of the Navier--Stokes
equations on the basis of the order-of-magnitude analysis that applied
in the friction-dominated boundary layer. The far-field is then
described by the Euler equations. These simplified equations enabled
both analytical and approximate solutions for a wide range of
applications. Pohlhausen~\cite{Pohlhausen1921} made a major contribution by
visualizing the boundary layer to be described by self-similar analytic
profiles where the profile parameters can be determined from the known
initial and boundary conditions. This led to a rejuvenation of interest
in fluid dynamics and great progress was made over the next few decades
in solving many problems of practical interest with acceptable
approximations.

Many real-life flows of great importance occur in 3 or 4-dimensional
space-time domain. They do not have a dominant direction of flow.
Examples include flows with separation or recirculation, impinging
flows with strong streamline curvature, and flows in complex
geometries. Such flows, actually, define the majority of flows
encountered in industry and in the environment, and are governed by the
complete Navier--Stokes equations. These stayed intractable and evaded
all attempts to obtain an analytic or approximate solution, and till
the middle of the 20{th} century, the solution of the complete
Navier--Stokes equations stayed a distant dream.

It is in this context that Computational Fluid Dynamics arrived on the
scene. What made it possible was the availability of commercial
electronic computers starting with the IBM 701 in 1952. It made
possible generalized, albeit approximate numerical, solutions of the
complete Navier--Stokes equations for a much broader class of practical
problems without having to resort to dropping one term or another in
the governing equations.

Newton famously stated that If I have seen further, it is by standing
on the shoulders of Giants. CFD also stands on the shoulders of
mathematicians who, over a period of 200 years, developed the tools
that would later form the foundation of CFD. These tools belonged to
Numerical Analysis and Matrix Theory. At the start of the 18{th}
century, Brook Taylor~\cite{Taylor1715} introduced the Finite Difference
Method (FDM) that laid the foundations of the differential calculus and
provided a very convenient tool to convert the governing differential
equations into algebraic equations that can then be solved by numerical
methods. Next milestone can be credited to Carl Friedrichs Gauss who,
over a period of 30 years starting in 1810~\cite{Smith1928}, made
seminal contributions that played an important role in CFD. He
introduced both a direct elimination and an iterative method for
solution of a system of linear algebraic equations---a matrix. He also
introduced the method of weighted residuals. His other important
contribution was to develop the Gauss (or Divergence) Theorem that,
among other uses, provides an equivalence between surface and volume
integrals for a closed volume. This can be used to decrease the order
of a differential equation by one. William Rowen Hamilton
\cite{Hamilton1833} made important contribution by defining the
concept of a Hamiltonian---in the form of minimum energy state for a
system---that proved useful in analysis of a dynamical system.
Cauchy~\cite{Cauchy1929} introduced the method of steepest descent for
non-linear equations that is the precursor to the popular Krylov and
Conjugate Gradients methods~\cite{Amritkaretal2015} for solving
algebraic equations. Jacobi~\cite{Jacobi1845} developed a relaxation type of
method for solving a system of algebraic equations. A more important
contribution of Jacobi was to introduce the concept of pre-conditioning
to speed up convergence. Seidel~\cite{Seidel1874} proposed a modification of
the Gauss method, today known as the Gauss--Seidel method, that allowed
improved convergence in many cases. Chebyshev introduced his
polynomials in 1854 that play an important role in iterative matrix
solution~\cite{Varga1962}. 

In the first half of the 20{th} century, rapid progress was made in the
field of numerical analyses. Variational method and the concept of
trial functions for solving partial differential equations, the
Rayleigh--Ritz method, was introduced independently by
Ritz \cite{Ritz1909} and Lord Rayleigh
\cite{Leissa2005}. Richardson~\cite{Richardson1910} made some important
contributions including the method of polynomial iteration, and a
method for improving the accuracy of a solution by extrapolation of
solutions obtained with different grids. Boris Galerkin
\cite{Galerkin1915} introduced the method of converting a
differential equation to its algebraic equivalent through basis
functions which forms the basis of the Finite Element Method (FEM).
Liebman~\cite{Liebmann1918} appears to be the first one to use an
iterative method for solving a PDE. His method is essentially a
Gauss--Seidel method which had previously been applied to matrices
resulting from ordinary algebraic equations. Krylov~\cite{Krylov1931}
introduced the concept of Krylov subspace that forms the basis of the
currently most successful iterative methods for solution of the
iterative equations. Southwell~\cite{Southwell1935} introduced the
successive over relaxation method to solve the system of linear
equations with significant improvement in speed for a certain class of
matrices. Young~\cite{Young1950} devised an efficient algorithm to
determine the over-relaxation factor for a class of elliptic equations.
George Forsyth made important contributions to Numerical Analysis. An
excellent summary of methods for solving linear algebraic equations is
available, for example, in Forsyth~\cite{Forsythe1953}. 

The 2{nd} half of the 20{th} century saw further development of both
the direct and the iterative methods for solving a matrix. Much
progress was made in developing efficient use of direct matric solvers
and Krylov subspace methods for iterative solution reached a point of
maturity. Saad~\cite{Saad2020} has compiled an excellent history of the
development of iterative methods. 

The field of CFD has mushroomed since the introduction of the
electronic computer. However, manual computations of fluid dynamics
problems had existed since a few decades earlier.
Richardson~\cite{Richardson1922} used a FDM to predict weather.
Thom~\cite{Thom1928,Thom1933} obtained a FDM solution for a
two-dimensional flow and a flow past a circular cylinder at a Reynolds
number of 50. He used the two-dimensional stream-function-vorticity
version of Navier--Stokes equations and employed street children as
human computers to compute the resulting matrix of equations.
Southwell~\cite{Southwell1946} derived the FDM solution for flow past a
cylinder up to a Reynolds number of 100.
Kawaguti~\cite{Kawaguti1953,Kawaguti1961} explored the problem of flow
around a circular cylinder and that in a driven square cavity. Allen
and Southwell~\cite{AllenSouthwell1955} added the innovation of
conformal coordinate transformation to map the infinite domain behind a
cylinder into a finite domain. Thom and Apelt~\cite{ThomApelt1961} made
an interesting observation that the range of Reynolds number for flow
behind a cylinder could be extended by the use of under-relaxation. 
Other notable computations were those of Simuni~\cite{Simuni1964} who
reported numerical solutions for a range of viscous flow problems. All
these researchers reported numerical instability as the Reynolds number
was increased. They had employed the central difference scheme, because
of its second order accuracy, to discretize their equations. It was
left to Courant \etal~\cite{Courantetal1952} to make a significant
contribution by recognizing that the central difference scheme is not
the correct choice to get a stable solution for high-speed flows when
the flow shows a strong hyperbolic behaviour. They discovered the
important role of the upwind, or one-sided, differences, for hyperbolic
equations.  However, this was not widely known to the rest of the
research community working on Navier--Stokes equations till well into
the 1960's.
                  
Dorodnytsyn~\cite{Dorodnitsyn1959}, Godunev~\cite{Godunov1959},
Gilinski \etal~\cite{Gilinskiietal1964}, Rusanov~\cite{Rusanov1966},
Yanenko~\cite{Yanenko1967} and a number of other Russian scientists,
many at the Institute of Applied Mathematics of the USSR (now Russian)
Academy of Sciences, contributed immensely to the understanding of the
numerical aspects of CFD. The work of Yanenko had particular influence
in later years as the Fractional Step Method that he proposed became a
very useful option to deal with complex partial differential equations.
Unfortunately, much of this work was unknown in the west till well into
the 1970's.

Before the advent of the electronic computer the use of numerical
methods for fluid flow was severely restricted due to the sheer number
of arithmetic computations. However, these explorations made important
contributions that highlighted the basic problems encountered in the
application of the FDM to Navier--Stokes equations. One interesting
fact is that before the advent of the electronic computer a ``person''
who performed the calculations was called a ``computer''. 

It should be mentioned here that turbulence plays a key, and often
controlling, role in the application of CFD. The formal foundations of
turbulence were laid by the seminal work of Reynolds~\cite{Reynolds1895}
half a century after the development of the Navier--Stokes equations.
Turbulence is, in itself, a very exhaustive and complex field of
research. Significant developments have taken place since the advent of
the electronic computers that have resulted in a wide range of
approaches to deal with the impact of turbulence on fluid flow.
However, the theme of this paper is CFD, hence these important
developments in turbulence are not reviewed here. Many exhaustive
reviews of turbulence and turbulence models have now been published and
these may be referred to by the reader for further information. 

\section{The basic physics and equations of CFD}\label{sec2}

The cornerstone of CFD is the law of conservation of a property of the
fluid. Consider a closed control volume, $V$, with a bounding area
$\mathbi{A}$ and a unit normal $n$ as shown in Figure~\ref{fig1}. 

\begin{figure}[h!]
\includegraphics{fig01}
\caption{\label{fig1}Illustrative control volume.}
\end{figure}

With $\mathbi{J}$ as the net outward flux vector through the area
$\mathbi{A}$, and $G$ as the net rate of generation, the balance of any
given property, ${\Phi}$, is given by:
{\begin{equation}\label{eq1}
 \int _{V} \dfrac{\partial \rho \Phi }{\partial t} \,\mathrm{d}V = {-}
 \int _{A} \mathbi{J} \boldsymbol{\cdot} \mathbi{n}\,\mathrm{d}A + 
 \int _{V} G\,\mathrm{d}V
\end{equation}}\unskip
${\Phi}$ here is a property per unit mass of the fluid and ${\rho}$ is
the mass density of the fluid. If we now use Gauss Theorem and, as
$\delta V \rightarrow 0$, this equation transforms to:
{\begin{equation}\label{eq2}
\dfrac{\partial \rho \Phi }{\partial t} = {-}(\nabla
\cdot \mathbi{J}) + G
\end{equation}}\unskip
{\vspace*{-1pc}}
{\[
\mbox{Rate of Accumulation} = \mbox{Net Influx through surface} + 
\mbox{Net Rate of generation}.
\]}\unskip
Fluxes associated with transport processes can be of many types. At
this stage we will consider only the most common ones: those due to
Convection and Diffusion. Convective flux is due to the motion of the
fluid at the continuum scale where the property moves with the velocity
of the fluid. The diffusive flux is due to the random Brownian motion
of individual fluid molecules that occurs at the molecular or
sub-continuum level. At the continuum scale, a net average of this flux
is expressed in terms of a diffusivity constant and the gradient of the
property of concern. The flux vector, $\mathbi{J}$, is then given by:
{\begin{equation}\label{eq3}
\mathbi{J} = \rho u\Phi - \Gamma \Delta \cdot \Phi.
\end{equation}}\unskip
In this equation, $\mathbi{V}$ is the velocity vector and $\Gamma$ is
the diffusion coefficient. The balance equation then can be written as:
{\begin{equation}\label{eq4}
\dfrac{\partial \rho \Phi }{\partial t} = {-}(\nabla \cdot V\rho \Phi)
+ (\nabla \cdot \Gamma \nabla \Phi) + G.
\end{equation}}\unskip
So, in essence, CFD is really very simple. It is the balance of a
property that moves with the fluid. To that extent, it is no different
than keeping an account of the bank balance; what is important is to
remember the physics behind this process and not the mathematical
details.

An important contribution made by Spalding at Imperial College was to
generalize this equation to the index format given below to represent
any entity that moves:
{\begin{equation}\label{eq5}
\dfrac{\partial }{\partial t} (\alpha \Phi ) + 
\dfrac{\partial }{\partial x_{i}} (\beta V_{i}\Phi)  = 
\dfrac{\partial }{\partial x_{i}} \left(\Gamma _{\Phi}
\dfrac{\partial \Phi}{\partial xj}\right) + 
S_{\Phi} - \alpha \Phi s_{\Phi}.
\end{equation}}\unskip
Here:
{\begin{eqnarray*}
&& \alpha \mbox{ is the accumulation or Storage Coefficient} \\
&& \beta \mbox{ is a specific property of fluid, such as mass density} \\
&& \Gamma _{\Phi} \mbox{ is the effective diffusivity Tensor} \\
&& S_{\Phi} \mbox{ is the source of property } \Phi \\
&& s_{\Phi} \mbox{ is the decay rate of property } \Phi.
\end{eqnarray*}}\unskip
This equation expresses the balance of any property that can move (in
continuum approach). This can be fluid itself, heat, mass. electrons,
cars, humans, fishes, birds, stars---anything in large enough numbers
to be considered a continuum.

\section{The origin of modern CFD}\label{sec3}

Modern CFD started with the advent of the Electronic Computer that
became commercially available in the early 1950's. Though there were a
few researchers who started solving the fluid flows problems, the birth
of what came to be known as Computational Fluid Dynamics (CFD), as we
know it today, can be attributed to two brilliant scientists---Frank
Harlow at Los Alamos National Laboratory and Brian Spalding at Imperial
College. Both worked independently of each other and did not know much
about each other's work till the 1970's. 

An IBM 701, the first commercially available scientific computer, also
known as the Defence Calculator, arrived at Los Alamos in 1953. Soon
thereafter, Harlow, leader of the T3 group, started a dedicated
exploration of fluid dynamics equations with the Finite Difference
Method. The primary interest of his group, was in otherwise intractable
fluid dynamics problems related to the US atomic weapons program known
as the Manhattan Project. Much of this work was classified and most of
the early work appeared as internal reports. Harlow, a physicist by
training, had no interest, in commercial applications of CFD. One of
the first hints of this work to the outside research community was a
paper on the particle-in cell method~\cite{Harlow1957}.

An IBM 7090, later upgraded to IBM 7094, arrived at Imperial College,
London, in 1962. Unaware of the work at Los Alamos, Brian Spalding at
Imperial college, had been trying to get a handle on complex practical
problems faced by practicing engineers. His focus was on industrial
applications and academic publications. By late 1950's Spalding was
already recognized as a major researcher and an expert in the field of
combustion. He had made innovative contributions in theoretical and
experimental combustion that would lead to the Spalding number or the B
factor that is widely used in mass transfer applications today. He was
also active in combining the key hydrodynamic, heat and mass transfer
concepts of von Karman~\cite{VonKarman1921}, Kruzhilin~\cite{Kruzhilin1936} and 
Eckert~\cite{Eckert1950} into a unified framework. He went on to write a
seminal book on Mass transfer~\cite{Spalding1963} that has influenced
many subsequent researchers. In 1964, he outlined a comprehensive
Unified Theory~\cite{Spalding1964} for boundary layers, jets and wakes
and with his students
\cite{RicouSpalding1961,Jayatillaka1966,EscudierNicoll1966,Escudier1967}
undertook an extensive in-depth survey of the experimental data to
determine the ``Universal'' empirical constants needed to account for
turbulence. He then explored generalized profile methods for numerical
solutions of boundary layer flows. His dedication over two decades to
explore these concepts and analyse experimental data, proved crucial in
developing the Finite Volume Method (FVM) for predicting flows of
practical interest that is the most popular method for solving complex
fluid flow problems. 

It must be emphasised that during these early years, there were other
researches who made significant contributions to the development of the
technology. These included Barakat and Clark~\cite{BarakatClark1965},
Burgraff~\cite{Burggraf1966}, Chorin~\cite{Chorin1967}, and a few
others. However, their focus was mostly on research and on solving
specific problems and not on developing generalised CFD methods for
distribution to other researchers.

It was Spalding who made CFD a common and popular tool for engineers
through publications, training courses and industrial consultancy. He,
with his students, published the first books on CFD and incorporated
the first consulting company in CFD. Most of the processes that impact
practical problems---turbulence, combustion, radiation, wall
interactions, multi-phase flows etc. were (and still are) somewhat
poorly understood. There is no generally agreed unique mathematical
formulation of such processes. Approximations, called models---actually
hypotheses---are therefore used to capture some essential aspects of
the impact of such processes so that acceptable answers can be obtained
for practical design purposes. Spalding boldly ventured into this
field---with the philosophy ``if it's good enough for design, it's good
enough for an engineer''---and derived methodologies that enabled an
engineer to gain some insight into complex problems.

In developing the Finite Volume Method (FVM), his group drew expertly
and widely on the work by other researchers and added new innovations,
explored new avenues and gave the CFD practice a coherent shape that
was easy and practical to use for an engineer.
Patankar and Spalding~\cite{PatankarSpalding1966} had already developed a generalized
profile method for solving boundary layer problems that fall in the
category of parabolic flows that have a dominant direction of flow.
However, a broad range of flows of practical importance fall into the
category of elliptical flows with no dominant flow direction. Such
flows include those in complex geometries, flows with separation or
recirculation, those with strong pressure gradients, and flows with
high stream line curvature. Spalding asked two of his students to
concentrate on such flows. Runchal started with the separated flows
behind a backward facing step and that in a driven-lid square cavity,
and Wolfshtein was assigned the problem of the impinging jet on a flat
plate. By the end of 1965, Runchal and Wolfshtein came to the
conclusion that the parabolic profile method could not be extended to
elliptic flows. The only general method then available to tackle the
complete form of the Navier--Stokes equations was the method of finite
differences. Further development of the profile methods was abandoned
at Imperial College and FDM became the method of choice for both
parabolic and elliptic flows. Runchal and Wolfshtein~\cite{RunchalWolfshtein1969} published
the first general purpose Navier--Stokes solver for two-dimensions
elliptic flows and Patankar~\cite{PatankarSpalding1967} published a general-purpose
solver for parabolic flows.

In 1967, Spalding used his insight to invent a physical rather than a
mathematical approach to the FDM. He envisaged the FD grid as a series
of Tanks and Tubes as shown in Figure~\ref{fig2}{\unskip\break} below.

\begin{figure}
\includegraphics{fig02}
\caption{\label{fig2}Tank and tube visualization of finite volume
method.}
\end{figure}

With this insight each node of a finite difference grid becomes an
independent tank which exchanges fluxes with other tanks by tubes
aligned with the grid lines. Though this was a simple conceptual
change, it had a profound effect in shifting the focus to physical
quantities such as fluxes and amounts of a property rather than
mathematical artifacts such as variables, derivatives and basis
functions. Runchal~\cite{Runchal1969} realized that the FVM can be generalized
by using an integral approach based on Gauss Theorem rather than a
differential approach based on the Taylor series which is the basis of
the FDM. With these insights, the Finite Volume Method (FVM) was born
and the governing equations were reformulated in integral form in terms
of accumulated quantities, fluxes, distances, areas and volumes.

Spalding had realized the commercial potential of CFD and a post
experience course was organized at Imperial College in 1969, which
targeted both academic and industrial communities. He also incorporated
the first commercial company in CFD called CHAM in 1969. Based on the
work of Runchal, Wolfshtein and Spalding 
\cite{Runchaletal1967,Runchaletal1969}, the first ever book on CFD
\cite{Gosmanetal1969} was published by Academic Press. Unlike the
focus on intellectual property today, this book contained the full
computer code (ANSWER) for solution of the 2-D Navier--Stokes equations
which was freely distributed. Spalding had a strict rule that all joint
publications carry the names in alphabetic order; Gosman, the editor of
the book is therefore listed first. 

The extension of FVM to three-dimensional flows came in 1972. Spalding
and Patankar~\cite{PatankarSpalding1972} developed a three-dimensional
velocity-pressure based method---SIMPLE---that revolutionized the field
of CFD and established the Imperial College approach as the de-facto
standard for solution of the incompressible version of the
Navier--Stokes equations. Many researchers have since developed
variations of improvements of the SIMPLE approach that have advantages
for specific applications. The SIMPLE approach is part of the family of
pressure-projection methods which are perhaps the most widely used
methods, especially for low Mach number flows
\cite{Kwaketal2005,Wangetal2018}. For High Mach number flows the so-called
density-based methods are generally preferred
\cite{HirschelKrause2009,vanLeerPowell2010,RizziLuckring2021}.
                   
Another technology that started soon thereafter was the extension of
the Finite Element Method (FEM) to fluids. FEM was originally developed
for analysis of structural problems. It employed the minimization
principles and Hamiltonian to find the minimum-energy stable state of a
system. It was originally applicable only to linear governing equations
such as those for structures. Zienkiewicz~\cite{Zienkiewiczetal1967,Zienkiewicz1974}
and his students had considerable success in extending it to fluid flow
problems by switching to variational principals. However, FEM, at that
time, could not deal with high Reynolds number flows or flows with
discontinuities. The path-breaking work in this aspect by Godunov
\cite{Godunov1959} was in the Russian literature and was not yet
widely known in the West. FEM also lacked a clear theoretical basis for
fluids because a unique minimum energy state does not exist for
nonlinear systems. A bigger drawback was that FEM uses global error
minimization. There are many systems, such as those with combustion,
where a small change or a small error, say in temperature, can lead to
orders of magnitude larger error in some other quantity, such as the
reaction rate. Local mass balance, therefore, can be very important for
these highly nonlinear systems. Today, the differences are semantics
since both FEM and FVM have appropriated and borrowed from each other.
However, on the practical side, most current commercial codes still use
FVM.\looseness=-1

As mentioned in the introductory section, turbulence plays a very
important role in practical flows. Before the 1960's, the dominant
method to deal with turbulence was the Prandtl mixing length model
\cite{Prandtl1925} that was based upon the concept of a universal
constant (von Karmon constant) and a mixing length that characterized
the scale of turbulence eddies. Theoretical foundations of more
advanced turbulence models were laid by Chou~\cite{Chou1945},
Prandtl~\cite{Prandtl1945},
Rotta~\cite{Rotta1951} and Davidov~\cite{Davidov1961}. However, before the
electronic computer, it was not possible to make any practical use of
these equations. It was but natural that much of the early work on
turbulence models occurred at Los Alamos and Imperial College lock-step
with CFD. Based upon the Reynolds-Averaged Navier--Stokes (RANS)
equations, Harlow and Nakayama~\cite{HarlowNakayama1968} developed a two-equation model
that allowed characterization of turbulence by a parameter for
turbulence kinetic energy ($k$) and one for a length scale of eddies
represented by the rate of dissipation of turbulence (${\varepsilon}$).
This was the so-called $k$-${\varepsilon}$ model.
Launder and Spalding~\cite{LaunderSpalding1972} generalized this approach and developed
``universal'' constants that are required for the application of this
model. The $k$-${\varepsilon}$ model became the de-facto standard for a
while. Saffman~\cite{Saffman1976} and Wilcox~\cite{Wilcox2006} developed a
$k$-${\omega}$ model where the length scale is now represented by the
vorticity (${\omega}$)---essentially a time scale---of the turbulence
of eddies. This offered some advantages in representing turbulence near
walls. However, turbulence is characterized by a very diverse spectrum
of kinetic energy and length scale of eddies. Representation by an
``average'' energy and an ``average'' length scale does not do it full
justice in most practical applications. An inherent assumption in
two-equation models is that turbulence is isotropic or that its
anisotropy can be defined by a single constant. This is obviously not
the case for most flows near walls or flows with strong pressure
gradients or stream-line curvatures. Over the years many turbulence
models have been developed to account for complexities of turbulence in
practical flows. These include Reynolds-stress models, hybrid
$k$-${\varepsilon}$-${\omega}$ models, LES, ILES, VLES, DES and DNS
models, and many other variations. It is outside the scope of this
paper to go into the details of these. Excellent reviews of turbulence
models are available elsewhere 
\cite{Kadivaretal2023,Bonnet2022,Cambonetal2022,SchiestelChaouat2022,
Visonneauetal2022,Lietal2021,
Yeetal2020,Yusufetal2020,
Duraisamyetal2019,Durbin2018,ArgyropoulosMarkatos2015,
Schiestel2010,Jaramilloetal2007}.

It must be stated here that in the 1960's and 1970's, with the wide
availability of electronic computers, computation of flows (and
structures) was a very active field. Many researchers made a number of
important contributions that made it possible for the Los Alamos and
Imperial College groups to develop their generalized approaches. A few
of these have been mentioned earlier. In addition to further work on
the FDM and FEM, considerable progress has been made in developing the
Spectral, and the Lattice Boltzmann approaches to CFD. Further, most of
the earlier work in CFD was based on the 1{st} and 2{nd} order schemes
which stressed stability over accuracy. Since then, many higher order
schemes have been developed to improve accuracy with limited
computational resources. Also, many significant developments have taken
place in numerical analysis and solution of the matrix of equations
that result from the application of these numerical approaches. A
complete anthology of these is considered beyond the scope of this
subjective and personal account.

One interesting fact is that the acronym, CFD, was not coined by either
Harlow or Spalding. It was first used by Pat Roache in his famous book
called Computational Fluid Dynamics~\cite{Roache1972}.

\section{CFD today}\label{sec4}

Since its introduction only a few decades ago, the field of CFD has
mushroomed into practically every aspect of human existence. IMARC
\cite{IMARC2019} estimated that the global commercial CFD market
will surpass US\$16.0 billion by 2024. Commercial CFD is only a small
part of the total CFD activities. The actual CFD market is at least an
order of magnitude larger since most of the CFD research and problem
solving occurs in government organizations and private R\&D departments
of large corporations and research organizations that is not directly
monetized.

The period since 1970's, has been a period of significant innovation
and consolidation. Early use of CFD was mostly for aerospace and
defense applications. However soon its potential was realized in other
fields. Most flows of practical interest are turbulent. Combustion and
chemical reactions play a very important role in energy conversion,
chemical processing and many other applications. Thus, significant
developments took place in developing models for turbulent flows,
interactions of turbulence with heat and mass transfer, and flows with
combustion and chemically reactions. These turbulence models made it
possible for wide use and acceptability of CFD. However, our focus here
is on CFD and not on turbulence which is a discipline in itself.
Excellent reviews of turbulence models are available which have been
referred to in the earlier section.

CFD has already made a significant impact on the hitherto conventional
applications in the aerospace, automotive, industrial and chemical
engineering. The literature is full of such examples and we will not
summarize these here. But among the non-traditional applications, there
are some areas where CFD is making a significant impact. Some of these
are summarized below. 

Prediction of weather has long been an area of major concern and CFD
has been the cornerstone of these predictions. Many developments in CFD
have had their origins in weather and climate related research.
Richardson~\cite{Richardson1922} undertook one of the first numerical
prediction exercises for any type of flow with a finite-difference
based approach to weather prediction. The widely used LES model for
turbulence has its origins in the work of Smagorinsky
\cite{Smagorinsky1963} to model atmospheric circulation. Lynch
\cite{Lynch2006}, Edwards~\cite{Edwards2011} and Pu and Kalnay
\cite{PuKalnay2018}, among others, have summarized the history and the
basics of weather prediction models.

Wind energy is playing an increasingly larger role in moving towards a
zero-emission future and CFD has been applied in this area for various
aspects including assessment of wind energy potential
\cite{RichardsHoxey1993,Ayotte2008,Sumneretal2010,YamadaKoike2011} 
and flow behind wind turbine wakes \cite{Vermeeretal2003,Sanderseetal2011}. 
Blocken~\cite{Blocken2014} has summarized the developments in Wind
Energy over the last 50 years.

The application of CFD to Environmental flows has a long history.
Starting in the early 1970s, such studies have contributed immensely to
the understanding of our physical environment consisting of the surface
water, ground water and the atmosphere. CFD has been applied
extensively to predict the behaviour of flow in these systems. This has
involved problems related to pollution as well as management of these
valuable resources. Analysis of Environmental flows is now a
burgeoning field with its own literature, professional societies and
government organizations looking after various environmental aspects.
CFD models now play a major role in decision making of a number of
government agencies and private organizations around the world.

Early examples of atmospheric and air quality models are available in
Lange~\cite{Lange1973}, Joynt and Blackman~\cite{JoyntBlackman1976},
Runchal \etal~\cite{Runchaletal1978}, and Pepper and Cooper
\cite{PepperCooper1983}. Since then, a vast literature concerning
development, applications, guidelines and reviews for various aspects
of atmospheric CFD models has been published. This includes, 
Di Sabatino \etal~\cite{DiSabatinoetal2007}, Tominaga \etal
\cite{Tominagaetal2008},  Pepper and Carrington
\cite{PepperCarrington2009}, Franke \etal~\cite{Frankeetal2011}, 
Leelossy \etal~\cite{Leelossyetal2018}, Toja-Silva \etal
\cite{TojaSilvaetal2018} and  Kadaverugu \etal
\cite{Kadaveruguetal2019}.

The surface water flow applications have been across the board to
rivers, estuaries, coastal waters and offshore areas. Effluent
discharges from industrial and urban activities into rivers, lakes and
coastal waters have been focus of CFD studies from the early days of
CFD (e.g.,~\cite{Rotty1974,LeeSengupta1976,Runchal1978}). Zhao \etal
\cite{Zhaoetal2011} and  Mohammadian \etal~\cite{Mohammadianetal2020}
provide a recent review of such effluent discharge applications. Miller
\etal~\cite{Milleretal2013a} have given an overview of water resource
modelling for both surface and ground water flows. 

Storm surges and Tsunamis have been an active area of application of
CFD since the 1970's~\cite{Leendertse1970,Runchal1975}. However, the
devastation caused by the 2004 Indian Ocean Tsunami and the risk posed
by 2011 Tsunami that affected the nuclear reactor in Fukushima has led
to an increased focus on numerical modelling of the runup and
inundation due to Tsunamis and storms~\cite{Morietal2011}, and in the
uncertainties associated with such forecasts~\cite{Godaetal2014}. 
Marras~\cite{MarrasMandli2020} has reviewed various tsunami modelling
approaches and recent developments due to the advent of GPUs and
exa-scale computing.

Groundwater has always played an important role in our water resources.
In the recent past, the quality of groundwater has been seriously
impacted due to leaching of chemicals and pesticides from agricultural,
industrial and urban activities. Study of these resources is therefore
a major incentive for proper management and utilization of these
resources. From the mid-1960's CFD models were developed to analyse
groundwater resources. About the same time, the Oil and Gas industry
realized the potential of these CFD models for optimization of
production. Later, starting mid 1970's, it was realized that the
continued use of nuclear power required an option for disposal of
nuclear waste. Various options were explored but, finally, the
consensus was that the most practical option was near-surface (for
low-level waste) or deep underground (for high-level waste) disposal of
Nuclear Waste. The heat generated by this waste can drastically modify
the ground water flow patterns and any leaching of this waste over time
can create a hazard as well as degrade the ground water quality. All
this has made CFD model development and application a very active
field. Early examples of such models are those from
Freeze and Witherspoon~\cite{FreezeWitherspoon1966}, 
Neuman and Witherspoon~\cite{NeumanWitherspoon1970},
Bredehoeft and Pinder~\cite{BredehoeftPinder1973}, 
Pinder~\cite{Pinder1973},
Trescott \etal~\cite{Trescottetal1976},
Dillon \etal~\cite{Dillonetal1978}, 
Narasimhan \etal~\cite{Narasimhanetal1978},
Runchal and Maini~\cite{RunchalMaini1980},
Reeves \etal~\cite{Reevesetal1986} and 
Runchal \cite{Runchal1985}. More recently numerous
reviews and applications of such models have appeared in the literature
\cite{ZhouLi2011,Milleretal2013b,PatilChore2014,Aderemietal2022,
Yanetal2022}.

Ground source heat pumps are an extremely efficient source for HVAC and
water heating that make use of the Earth's relatively constant
temperature throughout the year. CFD is being used increasingly to
explore the potential for ground source heating and cooling and
optimizing the design of such heat
pumps~\cite{Congedoetal2012,Luoetal2015,Wangetal2018b}.

With global warming leading to higher average temperatures, urban
planners and designers are using CFD to understand and mitigate the
urban heat island effect~\cite{ChungChoo2011} that requires increasing
capacity of air conditioning to make urban space habitable. Here CFD is
playing an increasingly important role in HVAC and energy efficiency
\cite{LawsonLawson2001}, and in fire safety
\cite{YeohYuen2009,McGrattanetal2013,Wangetal2014}, of buildings. CFD
is also being used by urban planners in novel ways to design evacuation
strategies during catastrophic events~\cite{Epsteinetal2011}.

Over the years, the food processing industry has made extensive use of
CFD \cite{Sastrohartonoetal1992,Anandharamakrishnan2013,Jaluria2018,
ManoharanRadhakrishnan2022}. Processes such as spray drying and spray
freezing of food use CFD to analyze and improve the product quality.
CFD is also being used to analyze heat transfer and air flow during
thermal processing of canned solid and liquid foods, and for the
analysis of processes such as baking, pasteurization, sterilization,
refrigeration, crystallization, and even for understanding of the
digestion process itself.

Competitive sport started making use of CFD in the early nineties.
Motor racing was among the first to utilize CFD for competitive
advantage~\cite{Hanna2012}. Apart from external aerodynamics of the
vehicle, CFD has been applied to engines, intakes, exhausts, radiators,
wheels, etc. Future applications will probably include maneuvers such
as overtaking and complex motion such as aquaplaning etc. CFD is also
used in sailing yacht design
\cite{Hedgesetal1996,Shyyetal2007,Peters2009}.  Apart from sail
design, CFD can be used to analyze the hull shape, the interaction
between the various underwater and above water components, sailing
tactics, etc. CFD has also been used to understand the aerodynamics of
ball sports and this trend will continue in the future
\cite{Hanna2012,Jalilianetal2014}. Wei \etal~\cite{Weietal2014} have provided a detailed
discussion of the fluid mechanics of competitive swimming.
Goff~\cite{Goff2013} has provided a review of research into the
aerodynamics of sport projectiles. Immonen~\cite{Immonen2022} has used CFD and
machine learning for the optimal design of golf balls.

In the last decade or so, CFD has started to revolutionize the
pharmaceutical, biomedical and healthcare fields. As may be expected,
the advent of the COVID pandemic immediately led to the application of
CFD to the fluid mechanics of disease transmission 
\cite{Bourouiba2021,Magaretal2021}.
CFD has been used extensively in the aftermath of the Covid-19 pandemic
to investigate the risks of disease transmission in indoor settings and
examine the proposed social distancing rules~\cite{Bourouiba2020}.
The study of airborne transmission of respiratory
pathogens under different ventilation conditions in indoor spaces is of
particular interest. CFD has been used to examine the dynamics of
virus-bearing aerosol droplets including growth in dry and humid air
\cite{Ngetal2021}, in confined spaces such as elevators
\cite{Biswasetal2022}, buses~\cite{Pavansaietal2021}, trains
\cite{Wangetal2022a}, aircraft~\cite{MohamadiFazeli2022} and in devising
mitigation prevention strategies via airflow modulation
\cite{Jeongetal2022}. CFD has also been used to ensure safe
transportation of SARS-COV-2 vaccines by refrigerated containers
\cite{Zhangetal2022a}. 

In a similar manner, computational hemodynamics is making critical
contributions to understanding the dynamics of the cardiovascular and
cerebrospinal fluids. CFD is now an important part of the
patient-specific cardio-vascular modelling, brain imaging, drug
delivery and transport in the human body 
\cite{KathawateAcharya2008,
Haaretal2014,
RahimiGorjietal2022,
Khanietal2022}. CFD has also made significant contribution to diagnosis
and treatment of diseases, development of medical devices, and surgical
planning
\cite{Ku1997,VandeVosseStergiopulos2011,vanBakeletal2018,TaylorFigueroa2009}.
CFD is already an important tool for research, risk assessment,
diagnosis, and treatment planning for aneurysms
\cite{Sforzaetal2012,Turjmanetal2014,Roietal2019,Qiuetal2022}. CFD is
also being used to estimate the wall shear stress in blood vessels
\cite{Wangetal2022b} and to analyze and predict the hemodynamics post
treatment of aneurysms with flow diverting stents
\cite{Zhangetal2022b}. 

Other applications of CFD have been in the clinical investigation and
surgical planning for the respiratory system~\cite{Basrietal2016,
Ayodeleetal2021,Aoyagietal2022,Cordaetal2022}, bone tissue
\cite{Piresetal2022}, the eye
\cite{Kumaretal2007,SiggersEthier2012,Wangetal2022c} and the inner
ear~\cite{Obrist2019}. 

In summary, CFD finds pervasive use in all spheres of human activity
and also helps in understanding and managing of the environment that
has become a critical issue due to the impact of human social and
economic activities.

\section{Primary challenges}\label{sec5}

Though CFD is now a very mature science, many challenges remain on the
horizon for it to become a practical and reliable tool in the hand of
an engineer. CFD today is a combination of three components: Science,
Art and Guesswork. Science is in the Navier--Stokes and the Boltzmann
equations, along with the equations for thermodynamics, heat and mass
transfer, that are based on accepted principals of mathematics, physics
and chemistry. Art comes in designing an efficient software tool
because, due to limited computer resources, compromises have to be made
in designing an efficient computing algorithm, selecting an economical
matrix solver and wrapping the CFD tool in an effective user interface.
What is really problematic about CFD, though, is the guesswork that is
required today to solve all but a subset of practical problems. The
fact is that, for many practical problems, we do not understand many
aspects of real fluids and the phenomena that govern the flow. It will
be more accurate to say that we understand less than we pretend to. For
many flows, CFD rests on uncertain scientific basis---we only have
hypotheses or empirical data to complete our simulations. This is the
Knowledge challenge.

This lack of knowledge---and the uncertainty associated with it---falls
under two categories---the Epistemic uncertainty that relates to our
lack of knowledge of a phenomenon and the Aleatoric uncertainty that
arises from inherently random and stochastic phenomena that we cannot
fully define or account for.

In the first category we have many aspects of real fluids where
underlying physics and chemistry is poorly understood with no agreed
mathematical basis. Most important of these is perhaps the lack of
understanding of turbulence. It has been 200 years since Navier
recognized linear and non-linear flows and almost 150 years since
Reynolds formally studied turbulence. It is turbulence models such as
those described in~\cite{LaunderSpalding1972} that made possible the
wide use and acceptability of CFD. Yet, there is no clear universal
acceptance of how to deal with turbulence. Neither has a widely
acceptable model of turbulence emerged after research spanning more
than 60 years. There are multiple models of turbulence---all with
limited applicability. What is worse is that, except for a very limited
class of problems, it is practically impossible to be \textit{a-priori}
certain that a given turbulence model will provide reliable
predictions. One has to conduct expensive empirical tests, or use
multiple simulations, or use multiple models, to obtain a measure of
confidence in the predictions.

Another area of lack of knowledge is that of chemical or biological
reactions. These play a very important role in many aspects of
industry, environment, and life itself. Most of our energy needs are
met today with processes that depend on chemical reactions. The
environment that we live in and breathe is a result of complex chemical
and biological reactions. Our body depends on chemical and biological
reactions to sustain itself. Yet, for most reactions of importance we
do not fully understand the reaction system or the rate of reaction. A
simple reaction between hydrogen and oxygen can be described as a 1- or
a 50-step reaction based on the operating conditions and the
objectives of the simulation. Further, the reaction rates for many
reactions of practical interest, such as for complex hydrocarbons that
are the mainstay of our energy generation, are not known except under a
limited set of thermodynamic conditions. The situation becomes even
more complicated when the flow is turbulent. The turbulence-chemistry
interaction is a poorly understood phenomena---again with multiple
models with limited applicability. There is no clear understanding of
many biological reactions that occur in living organisms such as those
in a liver, blood flow or the blood--brain barrier. Another example is
that of the reactions that occur deep in the earth such as those in
molten rocks at very high temperatures and pressure. There are many
other examples of Epistemic uncertainty that impact application of CFD
to problems of practical interest.

The second kind of uncertainty, the Aleatoric uncertainty, also plays a
very important role in CFD simulations of real-life systems. Often, we
don't even know the domain of influence for such simulations and many
of the variables that affect simulations are inherently stochastic or
unknown. For example, the synchronic initial and boundary conditions
for atmospheric, ocean or groundwater simulations, are almost never
known with acceptable certainty. For ground water simulations, the
spatial or temporal distribution of the properties, such as the
permeability or porosity of the soil, are known only at a few selected
locations though they are required throughout the domain of study. All
such quantities have to be either guessed as ``mean'' or ``most
probable'' values or assumed to be stochastic functions with assumed
probabilities. Often, therefore, multiple simulations are needed to
gain any confidence in the predictions. When multiple variables are
stochastic, the problem becomes orders of magnitude more difficult.
There is no exact mathematical methodology to deal with multivariate
uncertainty. Resort has to be made to assumptions about the interaction
between such variables with decreasing level of confidence as the
number of stochastic variables increases.

Another group of challenges involve multiple coupled phenomena. Almost
all practical systems involve multiple physical phenomena such as
fluid-structure interaction, heat or mass transfer, chemical reactions
and electromagnetism. Thus, multi-physics is, or will become, the norm
for many CFD simulations. Many problems involve temporal evolution and
structural deformations. Examples are the folding and flapping of wings
of aircrafts or the blood flow through a pulsating heart. These add an
order of magnitude of complexity and uncertainty due to multi-physics
interactions and phenomena that are poorly understood. 

Many problems of great practical interest involve phase change and/or
interaction between multiple fluids or fluids and solids. Phase change
plays a very important role in energy generation as well as
refrigeration. Droplet formation and evaporation are often an integral
component of energy production. The processes involved in such flows
are so complex and dependent on operating conditions, that currently
the only reliable method of simulation is through empirical relations
that apply to distinct regimes of flow. For example, there are many
gaps in our knowledge of how, for example, a jet of water breaks down
into droplets in the ambient air or when condensation or evaporation
occurs such as in the case of water and steam. When multiple fluids are
present, our lack of understanding is even more limiting. Processes
such as atmospheric interaction with the ocean, or the flow of water,
air and immiscible oil in a typical groundwater situation are poorly
understood. In such cases, the only option is to resort to empirical
data or models that are approximations with limited applicability. The
last half century has been rich in attempts to apply CFD to understand
the physics and processes involved in such complex coupled processes.

Another complex and poorly understood phenomenon is that of crystal
formation and growth which plays an important role in energy generation,
semiconductor manufacturing, computing, space research,
pharmaceuticals, and other emerging technologies. The complex processes
involved in this field are multidisciplinary in nature and involve
multiple order of magnitude scaling. For example, understating of
striation and epitaxial growth modes are challenges that are poorly
understood. In 1902, Verneuil developed his well-known growth method,
also called, flame-fusion, that formed the basis of modern crystal
growth technology. In 1915, Czochralski developed a method for growth
of single semiconductor crystals which was then improved on by Teal,
Little and others \cite{Czochralski1918,TealLittle1950,Demianetsetal1990,
Scheel2000}. Crystal growth is a very active subject for mathematical
modelling and computer simulation~\cite{Prakashetal1987,
ElGanaouiBontoux1998,ElGanaouietal2002,Kolmychkovetal2005}. 

Anther complexity that arises is that of multiple disparate spatial or
temporal scales. For example, in a flow involving combustion, the
reaction rates can vary from femtoseconds to milliseconds while the
flow scale can be on the order of seconds. Air pollution in the
vicinity of a building is impacted by both local and meso-scale flows.
This requires resolution near a building at scales on the order of
meters while the atmospheric flow has a range of influence that extends
to hundreds or thousands of meters. If the Knudsen number is high, or a
problem involves interaction between nano and macro scales, the
governing equations may transition from Navier--Stokes to Boltzmann.
Though we may understand the processes involved, such problems, often
present computational challenges and require resources beyond those
currently available. Compromises have to be made to resolve these
multiple scales that impact the accuracy and reliability of
predictions.

An important class of phenomena that requires computational resources
beyond those currently available is that of the integro-differential
processes. Such processes imply influences that involve spatial or
temporal coupling between remote points and add spatial or temporal
integration terms to the governing equations. Thermal radiation is a
good example where each point in space radiates to the other. Other
prominent example is that of the aerosols where the individual droplets
can coalesce, break down and move at different speeds than that of the
fluid or other droplets. In general, such processes are poorly
understood and approximations have to be made.

Finally, we face some computational challenges. After 50 years of CFD,
we still have no satisfactory method to deal with the convection terms
in the equations. The flow regime changes from elliptic to parabolic to
hyperbolic because of interaction between convection and diffusion
terms. Diffusion is a second-order process that is linear in nature.
Convection, on the other hand, is a first-order process with
nonlinearities. No robust and accurate discretization scheme exists
that meets all the requirements of actual physics of flow. Once the
governing equations are discretized, we face another challenge. The
governing equations are non-linear and admit multiple solutions
(turbulence is a good example), but the available computational methods
generate a linearized matrix of equations. Further, in today's
engineering practice it is not unusual to have a matrix with millions
of algebraic equations. There is no satisfactory matrix solver that is
both robust and accurate. Most common practice is to use an iterative
solver that may or may not produce a solution. If a solution is
produced there is no guarantee that it is accurate and admissible. So,
there is an overwhelming need for a robust, accurate and economic
matrix solver. Another challenge on the horizon is that the computer
architecture is evolving towards parallel and cloud computing systems.
CFD methods are mostly based on implicit or partially implicit
algorithms that were primarily developed for computers that were
sequential or had a few parallel pipes. The CFD algorithms today do not
scale linearly with processors or make optimal use of cloud computing
with, say, thousands of processors. Most run into computational or
communicational limitations. 

\section{The future of CFD and CFD of the future}\label{sec6}

As mentioned earlier, there are many challenges that need to be
overcome. CFD and its affiliate subjects, are among the most active
fields of research. Considerable advancements are likely in the next
couple of decades that will put CFD on firmer grounds. These will
include developments in turbulence, physics, fluid and material
properties and processes, numerical algorithms, matrix theory, geometry
and mesh generation, visualization and post-processing, uncertainty
quantification, higher order schemes, and user interfaces. A more
complete discussion of these is given by Runchal and Rao
\cite{RunchalRao2020}.

The CFD codes as we know them today will disappear; they are time
consuming to learn and difficult to use. Devising a grid that produces
a good solution can take weeks if not months. The CFD codes are
resource intensive and it is not unusual for a simulation to take weeks
or even months. More important, it is difficult to quantify the overall
reliability of CFD predictions for critical applications. Even with
modern computing power and recent advances in physics-based models, it
is still very challenging to represent the physics of the flow
accurately. One can design a gas turbine engine by using CFD but it
would not be wise to put that in an aircraft without extensive tests
and modifications of the engine. Such fundamental limitations still
continue to plague us today. A more complete discussion of these
limitations and the future prospects of CFD is available in
Patterson~\cite{Patterson1978} and Runchal~\cite{Runchal2012}. Instead of the general
purpose CFD codes, it is most likely that most CFD applications will be
through ``Apps'' that relate to dedicated niche applications. These CFD
Apps will embed a CFD engine and extensive data libraries in an
intuitive user interface with Virtual Reality software and design
tools.

A fundamental driving force for change, in the form of Machine Learning
and Artificial Intelligence (ML/AI) has made its impact in the last
decade. Of course, ML has been explored since the 1980's 
\cite{Rumelhartetal1986}. But primarily due to limitations of
computational resources, its use was neither extensive not pervasive.
It is only in the last few years that it has evoked considered
interest. 

The ML technology is already transforming our lives and we believe that
it will have a substantial impact in almost every aspect of life.
Today, the profound impact of this technology is evident in
self-driving cars, natural language interfaces, speech recognition and
medical imaging---to name a few. Now ML is beginning to have a
considerable impact in the field of CFD
\cite{Raissietal2017,KaniElsheikh2017}.  Machine learning is being used
to improve the accuracy as well as the computing efficiency of CFD
solvers \cite{VinuesaBrunton2021}.  Several researchers have trained
deep neural networks to replace or enhance the usual spatio-temporal
discretization schemes
\cite{BarSinaietal2019,Asem2020,StevensColonius2020,Ranadeetal2021,Madduetal2021}. 
Learned discretizations have also been applied to shock capturing
\cite{StevensColonius2020}. Kossaczka \etal~\cite{Kossaczkaetal2021}
have developed a fifth-order WENO shock capturing scheme which uses
optimized coefficient output by a trained deep neural network. Partial
differential equation solvers assisted by trained deep neural networks
have demonstrated accuracy enhancements and performance gains
\cite{BelbutePeresetal2020}. Deep learning approaches have also been
applied to speed up the large linear systems that arise during the
numerical solution of partial differential equations 
\cite{Kanedaetal2022}.

Currently, the Artificial Neural Network (ANN) appears to be the most
powerful machine learning technique. It can be used for the
(approximate) solution of ordinary or partial differential equations as
well as differential algebraic systems. The first known use of a neural
network for the approximate solution of differential equations was
by Lagaris \etal~\cite{Lagarisetal1998}. 
A typical ANN consists of an input layer, one or
more hidden layers and an output layer. Each layer takes its input from
the preceding layer, performs a linear transformation with a set of
weights followed by a non-linear transformation via an activation
function and delivers its output to the next layer. ANNs with a large
number of hidden layers are typically referred to as Deep Neural
Networks [DNN]. When a DNN is constrained so that it obeys the
governing principles of physics, it is called a Physics Informed Neural
Network (PINN) or Physics Informed Machine Leaning (PIML)
\cite{Raissietal2017,Zhuetal2019}. \looseness=-1

For a typical Machine Learning CFD application, the known initial and
boundary conditions are input into a deep neural network. The output
obtained is then substituted into the governing set of equations that
may consist of Navier--Stokes and other multi-physics equations. If the
residue criteria for the governing equations are satisfied, it is
considered an admissible output. If not, then the input is
appropriately modified to re-train the neural network and update the
neural network parameters. Any experimental or real-life data available
for the system is used as an additional constraint to train the neural
network. ML, augmented with real-life data, can also be trained to deal
with highly nonlinear problems and augment knowledge of phenomena that
are not well understood. For example, if the RANS turbulence model
doesn't work for some particular application, a PINN can assist in
adding a term that allows the model to better simulate the flow
\cite{Duraisamyetal2019,Lingetal2016}. The PINN has been used to
provide data driven subgrid scale closures for large eddy simulations
\cite{Wangetal2018a} and for near wall modelling for LES 
\cite{Yangetal2019a}. Another recent innovation is the development of
a PINN for modelling of multiphase flows in porous media
\cite{Hannaetal2022}. Modifications have been proposed to PINNs to
increase their accuracy for problems that exhibit multi-scale behavior,
such as in turbulence, by enforcing temporal causality through a
weighted loss function that includes the residual at all times steps up
to the current step~\cite{Wangetal2022d}.

One use of ML is to develop a surrogate or reduced order model (ROM)
that captures some essential components of the behavior of a complex
system with substantially reduced storage and CPU requirements 
\cite{Lassilaetal2014}. These model order reduction methods (MOR) are
typically based on Proper Orthogonal Decomposition (POD), Principal
Component Analysis (PCA) or Dynamic Mode Decomposition (DMD) or 
Reduced Basis Methods
(RBM)~\cite{Berkoozetal1993,BuiThanhetal2004,Narasimha2011,
Tuetal2013,Proctoretal2016,Liangetal2018}.

An alternative, and perhaps more promising, ML technique for dealing
with the resource limitations is to develop a Digital Twin for a
complex system which is a virtual replica of a physical product,
process or service. Digital Twin was named as one of the top ten
strategic technology trends by Gartner~\cite{Kerremansetal2019}. In
the CFD context, the basic process consists of defining an objective
function that is of primary interest for engineering design. Then, a
sufficient spectrum of input variables for that system are fed into a
PIML that is constrained by the governing equations, or a ROM if
available, and any known empirical data. A matrix of the output
objective function is then generated as a function of the input state
of the system~\cite{GrievesVickers2017}. Once defined, a Digital Twin
replaces the target physical component in the CFD simulations, or it
can be used directly for purposes of engineering design.

Other parallel developments that are likely to impact CFD are the
Application-Specific-Integrated Chips (ASICs) and System-on-a-Chip
(SOC's). As the names imply, these essentially encode a specific
computational algorithm in hardware and usually provide an order of
magnitude, or more, speed up over the equivalent implementation in
software. A typical example of ASIC is that for Fourier-Transform
\cite{Reddyetal2014,Nooretal2018}. With these developments, CFD will find
application in almost every sphere of human activity. The current
generation of CFD codes will be overshadowed by niche and
application-specific software tools augmented by Digital Twins, ROMs,
ASICs or SoCs. This will make it possible to embed CFD engines in
engineering design systems. Beyond traditional applications, we expect
CFD engines to be embedded in common appliances and devices such as
refrigerators and HVAC systems. CFD engines may also be built into
analysis tools such as MS Excel or OpenOffice to assist in design
calculations. CFD analysis will be a standard tool for a practicing
engineer rather than the domain of a CFD expert.

\section{Conclusion}\label{sec7}

Since its inception in the mid-20{th} century, CFD has become an
invaluable tool for almost every sphere of human activity. However, a
number of challenges remain. In many instances we lack both the physics
and the mathematical tools to adequately understand the behavior of
fluids. As regards physics of fluids, there is a knowledge gap in our
basic understanding of many of the properties, processes and phenomena
that govern fluid flow. Additionally, the governing equations are
non-linear, but due to lack of adequate mathematical tools, we resort
to their linearized equivalents for numerical solution. Also, many
practical problems involve stochastic phenomena which are poorly
understood and hard to describe mathematically. The current generation
of CFD tools are resource intensive and require extensive training and
experience. They are neither robust nor provide a high level of
confidence in the prediction. One has to resort to testing, empirical
data or sensitivity studies to ensure that the predictions are reliable
for practical purpose. The on-going research and development activities
will lead to a better understanding of the physics of flow and
improvements in CFD tools. However, we expect that compromises will
have to be made for the foreseeable future.

The CFD of the future will be revolutionized by machine learning and
advances in computer hardware, big data and virtual reality tools. PIML
will enable development of reduced order models (ROM) and Digital
Twins. Advances in hardware will lead to ASICs and SoC's for niche
applications. Big data tools will make it possible to build libraries
of fluid properties and application-specific simulations that can be
used to train neural networks and enhance confidence in predictions.
Advances in hardware and software and big data libraries, and Virtual
Reality tools will enable a practicing engineer to effectively
visualize a fluid system, interact with its components, conduct CFD
simulation, and control its behavior through dynamic intervention.

The current generation of CFD codes will be overtaken by CFD software
tools that are application-specific, easy to use and limited in scope.
CFD engines will be embedded in engineering design systems such as
CAD/PLM, software or other design and analysis software tools. Built-in
CFD engines may become an integral part of appliances and devices to
enhance their performance and control. 

The engine in an automobile has evolved from being the focus of
attention for the user to being an inconspicuous component that the
user is seldom concerned about. The focus of the user is on the
interface with the vehicle and performance of the vehicle. CFD of the
future will be more like the engine that drives a car---it is essential
for the functioning of the device but the user does not have to be
concerned with its intricacies.

\section*{Declaration of interests}

The authors do not work for, advise, own shares in, or receive funds
from any organization that could benefit from this article, and have
declared no affiliations other than their research organizations.

\section*{Acknowledgements}

I am greatly indebted to Dr.~Madhukar Rao of Analytic \& Computational
Research, Inc, for his invaluable assistance in helping me improve this
paper with his extensive knowledge and contribution. He particularly
brought me up to date on the range of CFD applications today and the
promising developments in ML/AI as it applies to CFD. Mr.~Pachalla
Rajagopal and Mr.~Ayush Karmwar helped with their extensive literature
search and checking the manuscript for completeness. I am thankful and
gratefully acknowledge their contribution.

\back{}

\bibliographystyle{crunsrt}
\bibliography{crmeca20220900}
\refinput{crmeca20220900-reference.tex}

\end{document}