Quantum computing for simulation of fluid dynamics

We present a pedagogical introduction to a series of quantum computing algorithms for the simulation of classical fluids, with special emphasis on the Carleman-Lattice Boltzmann method.


I. INTRODUCTION
Quantum computing has been vigorously pursued in the recent years, mostly in view of the spectacular speedup it may provide as compared to the performance of electronic computers, for a number of scientific applications [1].Quantum computing can be traced back to Richard Feynman's epoch-making paper, in which he observed that physics "ain't classical", hence it ought to be simulated on quantum computers [2].Following Feynman's observations (with due credit to previous investigators such as T. Toffoli and E. Fredkin), early theoretical work on quantum computing was performed in the 1980s, e.g.Deutsch's work on the link between quantum theory, universal quantum computers and the Church-Turing principle [3].Then, with the publication of Shor's algorithm for integer factoring and Grover's search algorithm in the middle of the 1990s, the research area gathered significant momentum in terms of theoretical work and quantum computing hardware as well.The research area of quantum computing has continued to grow ever since [4,7,8].
In terms of applications, the simulation of quantum many-body systems has received special attention, due to its scientific and industrial relevance, as well as due to its close link with quantum hardware, meaning by this the possibility of mapping quantum hamiltonian directly into native quantum gates.In this paper, we shall focus on a much less beaten track, namely the use of quantum computers for the simulation of classical fluids [50].To this purpose, it is convenient to refer to the physics-computing plane defined by the following four-quadrants: CC: Classical computing for Classical physics; CQ: Classical computing for Quantum physics; QC: Quantum computing for Classical physics; QQ: Quantum computing for Quantum physics.as shown in Figure 1 (taken from [5]).
Feynman's observation pertains to the CQ sector shown in Figure 1, where one is often confronted with exponential complexity barriers associated with the phase-space of quantum many-body problems [11,12].The basic idea is that such exponential barriers can be handled by the corresponding exponential capacity of qubit representations offered by the QQ quadrant.In this Chapter, we shall focus on the opposite off-diagonal QC quadrant, where the power of quantum computing might be brought to the fruition of hard computational problems in classical physics.In particular, we shall focus on computational fluid dynamics, which presents a ceaseless quest for better and more efficient algorithms and computers, mostly on account of the problem of fluid turbulence [5,13].The physics of fluids presents two general features which are not shared by quantum physics: nonlinearity and dissipation.Both set a major challenge to quantum computing, which draws much of its power from the superposition principle and on the unitary and hamiltonian structure of quantum mechanics.In the following we shall present a few potential strategies to deal with both issues above.

Nature ain't classical…
II. CHALLENGES FACING QUANTUM COMPUTATIONAL FLUID DYNAMICS (QCFD) As mentioned in the introduction, realising the potential of quantum computing means leveraging distinctive features of quantum mechanics that, by definition, are not available on classical computers.However, it also follows that it is precisely these specific features that expose major challenges in realising simulations with a quantum advantage.
The main quantum mechanical concepts spawning the potential benefit of quantum algorithms are quantum superposition and quantum parallelism.The quantum state in an Q-qubit coherent register can be described by the Schrödinger wave function, defined by 2 Q complex amplitudes for 2 Q states in superposition.The square of each of these amplitudes defines the probability of finding the system in the corresponding state after quantum measurement.By encoding classical data in terms of these amplitudes an exponential saving in storage can be achieved when the number of qubits is compared to required number of classical bits.
Let us illustrate the idea for the specific case of turbulent flow simulation.Turbulence features a Re 3 complexity, where the Reynolds number Re represents the relative strength of convection (nonlinearity) versus dissipation.Most real life problems feature Reynolds numbers in the many-millions; for instance an airplane features Re ∼ 10 8 , implying O(10 24 ) floating-point operations per simulation.This is basically the best one can afford on a nearly-ideal Exascale classical computer.The simulation of regional atmospheric circulation flows takes us at least another two decades above in the Reynolds number, hence totally out of reach for any foreseeable classical computer [16,17].On the other hand, the minimum number of qubits Q required to represent Re 3 complexity can be roughly estimated as 2 Q = Re 3 , namely: This simple estimate shows that 80 qubits match the requirement of full-scale airplane design, while O(100) qubits would enable regional atmospheric models [16] or nuclear fusion applications [18].However, several key challenges stand in the way of this task.First, Quantum measurement: Extracting classical information collapses the quantum wave vector into a single eigenstate (a house of cards).Hence, to get classical values for each of the amplitudes, multiple realisations of the quantum state vector are needed with an associated set of measurements.
Second, Conditional operations: In the quantum circuit model, the 'classical' information is not available to the quantum gate operations performed in the circuit.Specifically, gate operations can be conditional on the state of one or more control qubits, while specifying gate operations conditional on one or more of the complex amplitudes defining the wave function is impossible.Therefore, when classical data are encoded in terms of amplitudes, a rotation angle in a quantum gate operation, this information is required at the time the circuit is compiled.This angle cannot be changed at run-time as a function of 'classical' data encoded in quantum amplitudes.
Third, Time marching: In an algorithm involving multiple iterations or time-steps, the overhead associated with quantum measurements used to extract classical data and re-initialization of the quantum state for a next iteration, scale quadratically with the grid size [19] and consequently they severely tax the quantum CFD efficiency.It should be observed that the well-known HHL [20] algorithm for linear system solution, assumes that the input and output data are encoded in terms of quantum amplitudes without including the cost of setting up the quantum state and extracting the classical solution.
Fourth, Circuit depth: With the exception of quantum measurement operations, quantum mechanical operators are unitary, linear and reversible, hence they can be directly implemented on a series of unitary quantum gates (quantum circuit model).Usually, in order to quantify the computational cost, producers and specialists take track of the number of two-qubit gates employed in the circuit.In fact, two-qubit gates, such as the CNOT gate, come in the hardware with an error that is several orders of magnitude (depending on the platform) larger than single-qubit gates errors.We know from the DiVincenzo theorem [49] that the number of two-qubits gates needed to implement a generic Q qubits unitary operation scales as 4 Q , thus voiding the exponential advantage of the embedding process.In order to sidestep this issue, the unitary operation that we want to perform in our circuit has to be combinatorially block diagonal [48].Hence, this provides a stringent constraint on the possible algorithms that can actually be implemented to simulate classical fluids.
If we further assume that the velocity field is represented in terms of amplitude encoding [21], i.e. the components of velocity vector at all grid points are represented in terms of the amplitudes defining the wavefunction, then the nocloning theorem prevents the use of (temporary) copies of any of these amplitudes.So, evaluation of u 2 or u ∂u ∂x cannot be performed by storing a temporary copy temp = u, to perform the computation of the value of u 2 as u × temp.Also, for data encoded in terms of the complex amplitudes of the Schrödinger wavefunction, there is a need for this state vector to have a unit norm, since these amplitudes represent probabilities of states.This means that for an operator attempting to compute the squares of these complex amplitudes, the resulting state vector loses unitarity.So, even without the no-cloning theorem complicating such a step, this points to a further problem with computing non-linear terms.This loss of unit norm of the quantum state vector represents an example of a non-unitary operation, typically associated with a corresponding loss of information.A similar argument runs for orthogonality of the eigenstates, since a nonlinear propagator rotates the state by an angle which depends on the state itself.
Fifth, Dissipation: Dissipation breaks hermicity of the quantum propagator.This can be recovered by augmenting the system with its hermitian conjugate, so that the doubled system is hermitian by construction.One can dispense with doubling degrees of freedom by representing the non-unitary update as a weighted sum of unitaries.This introduces a failure probability which needs to be handled probabilistically, i.e. by executing the update conditional on the successful measurement of an ancilla qubit.
To summarize, dealing with non-linear and dissipative terms raises major extra challenges for Quantum Computational Fluid Dynamics (QCFD) besides the ones commonly encountered for the simulation of quantum systems.

III. HYBRID QUANTUM/CLASSICAL APPROACHES
The challenges sketched above have resulted in the fact that most of the existing work in QCFD is based on a hybrid quantum/classical approach, with the quantum processor performing computations for which efficient quantum algorithms exist, while the output is then passed on to classical hardware to perform further computational tasks not (yet) amenable to quantum algorithms.
Figure 2 illustrates the idea.The quantum state ψ 0 is advanced to the next quantum state ψ 1 via a QQ algorithm.The quantum state ψ 1 is then measured to generate classical observables C 1 which are advanced to C 2 by a CC algorithm.The classical observables C 2 are then used to reconstruct the quantum state ψ 2 , ready for the next QQ step.The Q2C conversion shown in Figure 2 involves quantum measurements and requires repeated measurements over a statistical sample of quantum states, since none of them can be reused.The C2Q reconstruction requires the preparation of all the quantum eigenstates, hence a full reset of the quantum circuits from scratch.Both operations impose a substantial computational burden on the hybrid algorithm, typically scaling exponentially with the number of qubits.
Examples of previous works using the hybrid quantum/classical approach include the works of Steijl [22], Gaitan [24] and Budinski [25].The algorithm presented by Gaitan uses Kacewiz's quantum amplitude estimation ODE algorithm [26] as applied to the set of nonlinear ODE's resulting from standard discretization of the Navier-Stokes equations.As shown, for certain 'non-smooth' problems (illustrated using the quasi-1D flow in converging-diverging duct with normal shock wave), the complexity analysis shows potential for exponential speed-up, so that the challenges associated with hybrid quantum/classical computing can potentially be overcome.For the linear advection-diffusion equation, Budinski [27] presented a quantum algorithm based on the Lattice Boltzmann method [28,29].The algorithm can perform multiple successive time steps with no need for quantum measurements and re-initialization of qubit register between the time steps if suitable re-scaling of solution vector is applied to deal with non-unitarity.In the quantum circuit model, the fact that velocity field is unchanged between successive time steps implies that the operator implementing u ∂u ∂x can be re-used in multiple time steps.Budinski [25] then extended the work to Navier-Stokes equations in streamfunction-vorticity formulation.Then, velocity-field updates during each time step means that quantum-circuit implementation of convection terms cannot be re-used during multiple time steps and that the classical value of u (as well as other flow field data) is needed to define quantum circuit implementation for the next time step.This shows that it is the non-linearity that forces the use of a hybrid quantum/classical approach, similar to previous work where Navier-Stokes equations were discretized, e.g. the hybrid quantum/classical algorithms for fluid simulations based on quantum-Poisson solvers [22].
In summary, key challenges for the hybrid quantum/classical approach are: i) Cost and complexity of (repeated) measurements; ii) Statistical noise due to sampling of the quantum solution; iii) Cost and complexity of (repeated) re-initialization.
The development of more efficient (re-)initialization techniques appears to be a key factor for the future success of quantum computing for fluids 2: Sketch of a hybrid quantum-classical algorithm.The illustration shows the steps involved in a single time step or single iteration, including preparation of the subsequent step/iteration.

IV. QUANTUM FLUID DYNAMICS STRATEGIES
For the sake of concreteness, let us refer to the Navier-Stokes equations for time-dependent incompressible flows: where u is the velocity vector, P is the pressure, defined for location x as function of time t.The kinematic viscosity is defined by ν and density has been conventionally set to a unit constant value.Equation (3) enforces mass conservation, while Equation ( 2) is based on momentum conservation in each coordinate direction.Equations ( 2) and (3) highlight that it is the convection term that represents the non-linearity, i.e. the second term on the left-hand side of Equation (2).Writing the Navier-Stokes equations in non-dimensional form, such that x and u are scaled by reference length L ref and U ref , respectively, it follows that in Equation ( 2) the term ν is replaced by 1/Re, where Reynolds number Re = U ref L ref /ν.For Stokes flow, i.e. with Reynolds approaching 0, nonlinear terms are vanishingly small, but still not to be neglected since they are responsible for non trivial long-range correlations especially important in biological flows [23].For high Reynolds number (turbulent) flows, obviously the nonlinear terms play the leading role.Different strategies have been adopted to simulate fluid dynamics on quantum computer, mostly involving hybrid algorithms.Since they have been reviewed in the recent Perspective [5], here we simply list them, directing the reader to the original literature:

A. Nonlinear quantum ODE solvers
This approach consists in treating the discretized Navier-Stokes equations as a set of nonlinear ODE's and advance them in time by developing bespoken quantum nonlinear ODE solvers.This approach has been pioneered by Gaitan [24] and demonstrated for the case of a Laval nozzle on grids with O(30 − 60) grid points over 1400 time steps.The critical issue, though, is that the concrete implementation of the quantum oracle is left to a classical computer.

B. Nonlinear variational quantum eigenvalue solvers
it has recently been proposed that variational quantum eigenvalue (VQE) solvers, a major tool of the QQ sector, might be extended to the fluid equations as well [30].The basic idea is to use quantum computers for the efficient generation of variational eigenfunctions by suitable circuit parametrization.The associated energy functional is then minimized through a classical procedure.The appeal of these methods is the possibility to import algorithmic knowhow from quantum mechanical applications, especially quantum chemistry.To date, however, there is no evidence that the procedure can meet the standards of accuracy required by computational fluid dynamics.

C. Functional approach
It is long known that any nonlinear problem can be mapped into a linear one in a space of higher dimensions, a technique also known as Carleman embedding or Carleman linearization [31].The Carleman embedding trades nonlinearity for infinite dimensions and nonlocality.An alternative procedure to recover linearity is to take a probabilistic approach and formulate a functional kinetic (Liouville) equation for the the probability distribution function (PDF) of the fluid velocity field.Formally [39]: where P = P [u] is the functional PDF and u = f (u) is the equation of motion (hence f (u) is a nonlinear and nonlocal operator in function space).This is linear by construction, regardless of the nonlinearity of the dynamics.Once the fluid equations are discretized on a grid with G lattice nodes, the Liouville distribution P G becomes a 3G-variate PDF P G (u 1 . . .u G ) which lives in a O(3G)-dimensional space.For large-scale flow applications, G can reach values up to multi-billions, hence fully under the "curse of dimensionality".Fortunately, for most practical purposes, low order marginals prove often capable to retain the essential physical information.Marginalization is known to introduce coupling terms between different marginals, thus raising the need for appropriate closures expressing the high order marginals as functionals of the low order ones.This is a classical topic in non-equilibrium statistical physics, which may draw significant benefits from modern developments in tensor networks theory [40,41].
The idea can be taken one level further by noting that Eq. ( 4) is linear, but still not unitary.In [6] the authors presented a method to relate the classical evolution and the unitary evolution of a quantum system, making use of the Carleman-Koopman embedding, which transforms the Liouville equation into an equivalent many-body Schrödinger equation.
This is formally very appealing, but still subject to the dimensional curse; even in quantum chemistry nobody knows how to handle molecules with billions atoms with ab-initio methods.In addition, [6] shows that for chaotic systems the spectrum of the hamiltonian contains continuous modes that are pretty hard to tame.
Always in the spirit of functional methods, in [42], the authors develop another form of Koopman-von Neumann approach based on the level-set method as applied to classical nonlinear field theories [42].The formalism is applied to hyperbolic PDE's and Hamilton-Jacobi equations, but its extension to the Navier-Stokes equations is addressed only marginally, making it difficult to draw firm conclusions.

D. Carleman embedding
Introducing indices m and n to denote interacting neighbors, and using Einstein summation, the fluid equations can be recast in terms of a linear operator L and quadratic operator Q, as follows, The pressure term was removed for simplicity, although this all but a minor item, since pressure-flow-stress coupling is going to affect the structure of the Carleman matrix.
At this stage, it should be noted that in the discretization of the velocity derivatives at lattice point l, the values of the velocity at one or more neighboring lattice points is used.This means that when the Carleman linearization is used, it introduces a new variable U lm = u l u m to formally generate a linear system.Marching the system of equations forward in time produces an ever growing hierarchy in which the Carleman variables at level k involve the Carleman variables at the next level k + 1.Note that at each level we are faced with a tensor of rank k, with O(Gκ k−1 ) independent components, κ being the sparsity of the Q matrix.
Furthermore the tensors occupy a neighborhood of the original field whose diameter grows linearly with the Carleman level.This shows that there is a major prize to pay to uplift the nonlinearity of the fluid equations, both in terms of increasing dimensionality and loss of locality.Nevertheless, in [32] the authors present an algorithm based on Carleman linearization along with the use of a quantum linear system solution for the solution of the one-dimensional Burgers equation.The presented algorithm shows a polylog scaling with the number of grid points, i.e. an exponential improvement over classical approach.However, for a given time span T , the algorithm shows a complexity T 2 P oly(logT ), i.e. a significantly increased time-complexity compared to classical case.The authors simulate shock formation in a one-dimensional Burgers flow, with 16 grid points over 4000 time-steps, with a fourth order Carleman truncation.They reach Reynolds numbers up to 40, an order of magnitude larger than predicted by the theoretical no-go analysis, a welcome discrepancy which begs for further analysis.The key question then is: is there a maximum Reynolds number above which the Carleman procedure fails to converge?Also, when it does converge, how many Carleman levels are required to reach the desired accuracy at a given Reynolds number?
The answer may well depend on the chosen representation of the fluid equations and in the following we turn to a promising candidate in this respect, namely the Lattice Boltzmann method [28,46].

V. CARLEMAN LATTICE BOLTZMANN
In a recent paper, Cheung and coworkers argued that lattice Boltzmann might be more convenient for Carleman linearization, since the fluid nonlinearity is formally encoded in the Mach number instead of the Reynolds number, which makes of course a huge difference.They performed a classical Taylor-Green vortex simulation based on a Carleman-Lattice Boltzmann scheme, showing excellent agreement at moderate Mach number, with just three Carleman levels [38].More recently, in Ref. [15], the authors showed the convergence of the Carleman method to the analytical solution with the Carleman system truncated at second order for fluids at moderate-low Reynolds number, up to 100 (see Figure 3).
Let us remind that the Lattice Boltzmann (LB) method is based on the following scheme: In equations: where f i (x, t), is the probability to find a representative fluid parcel with discrete velocity v i at position x and time t.
In the above, f eq i is the corresponding discrete local equilibrium and τ is a local relaxation time, fixing the viscosity of the lattice fluid.Importantly, the local equilibrium is a quadratic function of the Mach number M a = u/c s , c s being the speed of sound.A defining feature of the LB method is that it allows to use a fixed, generally small, number of discrete velocities v i .For instance, it turns out that just nine discrete velocities are sufficient to fully recover the physical properties of a two-dimensional flow defined on a square lattice.The nine velocities are depicted in Fig. 4 and described in Table I.The model is called D2Q9.A key feature of LB, not shared by the macroscopic representation of fluid dynamics, is that streaming proceeds along straight lines, defined by the discrete velocities.This operation is exact, in that it does not involve any floating-point calculation but just a memory shift from the lattice position x to another lattice position x + v i ∆t.Also to be observed that collisions are completely local, hence perfectly amenable to parallel computing.Because of this, dissipation is an emergent property which does not require any second order spatial derivative, as it is the case for the Navier-Stokes equations.Besides their mathematical elegance, these properties are at the roots of the computational efficiency of the method and its amenability to parallel computing.
A key point of using the discrete-velocity Boltzmann formalism instead of Navier-Stokes is that, owing to the double dimensional phase-space, in the latter non-locality (streaming) is linear while non-linearity (collision) is local, while in Navier-Stokes the two merge into a single u∇u convective term.
On classical computers this disentanglement proves extremely beneficial because of the local nature of the collision term, which makes the system highly amenable to parallel computation.On quantum computers instead it is not FIG.4: The D2Q9 model, defined on a two-dimensional square lattice, with the corresponding set of nine velocity vectors labeled by index i, cf.Table I.
Since the collision step is local and depends only on the velocities of the fluid in the lattice site, one may apply an operator Ω which performs the collision step on the quantum register adapted to the velocities and an other operator Ŝ encoding the streaming process to be applied to all quantum registers, being it both non-local and velocity dependent.This ideal circuit is shown in Fig. 5(a).It turns out that since Ω is a non-linear and non-unitary operator, this circuit cannot be realized, unless specific symmetry conditions hold, such as those used in Ref. [9].Unfortunately, these do not conform to the universal quantum computation paradigm.
We may then opt for hybrid strategies whereby only streaming (collisions) would be performed on a quantum computer, leaving collisions (streaming) to a classical one.Interestingly, this leads to different circuits, as sketched in Figure 5. Fig. 5(b) and Fig. 5(c) show the circuits for collision-free streaming and streaming-free collision processes, that we are going to discuss next.

Collision-free streaming
The quantum computing algorithm for streaming is based on the approach used by Steijl and co-workers [34].The streaming is a non-local process shifting the particles in the next neighboring site, depending on the corresponding velocity, without affecting their value.Thus, this step can be executed by moving the probability function f i (x, t) to f i (x + v i , t + 1) at the next time step.
The corresponding operation on the quantum state reads as follows: i.e.only the positions register is affected by the streaming process, whereas the velocity register acts as control.The streaming operator is obtained as a combination of the streaming in the different directions, for the D2Q9 model with 8 velocities (plus the null velocity) is Ŝ = 8 i=1 Ŝi , that are therefore applied depending on the value of the qubit in the velocity register.The streaming can be obtained by applying a series of multiqubit controlled operations on the position register.We show in Figure 6 an example for the Ŝ1 streaming operator when applied on a 8 × 8 lattice.I, for a D2Q9 model define on a 8 × 8 lattice.The first three qubits are used to encode the position in the x direction and the latter three for the position in the y direction.They are applied to the quantum state only if the velocity register has value 1 or 2 respectively.

Streaming-free collision
A completely different approach is employed to perform the collision step: Carleman linearization.Based on the Lattice Boltzmann method, Itani et al. [33] developed an approach termed Carleman for second-quantized Lattice Boltzmann.This terminology stems from the fact that the Boltzmann operator, defined by: can be expressed in terms of the second-quantization annihilation and generation operators via the relation ∇ = (â − â+ ).Collisions follow the bosonic encoding first proposed by Mezzacapo et al. [35], whereby dissipative effects are represented as a weighted sum of two unitary operators.The scheme is as "Feynmansque" as it can possibly get, since it builds on a one-to-one analogy between LB and the Dirac equation, as first proposed in [36].As a result, the formal solution f t = e Bt f 0 can be computed in analogy with quantum mechanics.However, it is subject to a number of additional questions: primarily truncation effects due to the finite number of bosonic excitations and the long-time behaviour of non-unitarity errors [37].On a similar but operationally different line, Sanavio et al. [15] developed a quantum circuit implementing the collision step of a Carleman-Lattice Boltzmann algorithm.Since the number of variables with the Carleman method increases exponentially with the truncation order k, the circuit needs more qubits to embed the relevant information.The collision step requires k + 1 additional quantum registers, embedding the truncation order and the products between local Boltzmann functions with different velocities, say f i (x), f i (x)f j (y), f i (x)f j (y)f k (z) and so on.With this setting, the collision quantum circuit is genuinely local, as it depends only on the velocity registers and it can be implemented using only a fixed number of two-qubit gates, regardless of the number of lattice sites.This meets the promise of quantum advantage but unfortunately it is inconsistent with the streaming process, which must therefore be directed to a classical machine.This may seem an ideal option for hybrid quantum computing, since streaming is a floating-point free operation.However, accessing data is by no means a cost-free operation and it is indeed known to take a significant fraction of the computer time spent on collisions [14].

Fully quantum algorithm: Collision and Streaming
In Ref. [15], the authors managed to generalize the circuit for the streaming step to Carleman variables of higher orders, thus allowing the subsequent application of streaming and collision step on Carleman variables with a single embedding and a single readout of the results after evolution.This results in a fully quantum algorithm.In order to embed different non-linear terms one needs 2(k + 1) additional quantum registers, one encoding the information about the degree of non linearity up to order k, and the other 2k to embed position and velocity of the multiproduct functions.
This method was found to meet with an exponential depth of the quantum circuit as a function of the number of qubits.The reason is that the collision operator cannot be written in the form of a sparse and combinatorially block-diagonal matrix, nor in terms of a compact combination of Pauli matrices.Hence, to date, it cannot be quantum compiled.
Possible ways out are currently under investigation: one option is the application of the Carleman linearization method directly to the Navier-Stokes equations.The advantage of the Carleman-Navier-Stokes procedure rests a much lesser number of Carleman variables, with a corresponding reduction of the circuit depth.The flip side, however, is a more complex structure of the Carleman matrix, which, combined with the Reynolds versus Mach number as a measure of nonlinearity, may impair Carleman convergence.A detailed analysis along this direction is currently underway.

VI. CONCLUSIONS
In summary, we have presented a survey of the a few leading approaches to the quantum simulation of classical fluids.Various obstacles stand in the way of the efficient simulation of fluid flows on quantum computers, especially at high-Reynolds numbers, mostly on account of the strong nonlinear effects.A few ways potential ways out have been illustrated, but their practical implementation commands major advances both on the algorithmic side, let alone quantum technology [47].
Before closing, it should be observed that while, everybody would expect quantum computers to handle turbulence, the physics of fluids is littered with interesting problems at low-Reynolds, especially in microfluidics, soft matter and biology [44][45][46], but also in high-energy physics, typically quark-gluon plasma experiments.Hence, even in case turbulent flows would prove beyond reach to quantum computing, there would nevertheless be plenty of interesting fluid applications in the low-moderate Reynolds regime.
For instance, it could be of great interest to devise a quantum multi-scale application, coupling quantum algorithms for biomolecules in a water solvent described by a quantum algorithm for low-Reynolds fluid flow.Given that low-Reynolds flows are nonlocal, perhaps the inherent nonlocality of quantum mechanics could prove helpful in representing the classical nonlocality of low-Reynolds flows.
On a more philosophical ground, should the physics of fluids show a case for "classical advantage", there would still be interesting lessons to learn.Indeed, as observed earlier on, while it is true that Nature isn't classical, it is equally true that Nature has a very strong tendency to become classical at macroscopic scales/high temperatures.Withstanding such tendency, or perhaps even take of advantage of it, is indeed the prime struggle of quantum computers.Looking into ways to achieve this very hard goal is not only of practical but also of major foundational interest, since it is likely to shed new light into the still open problem of "classicalization", i.e. the emergence of classical behaviour out of underlying quantum mechanical microscopic physics.This is only one (outstanding) example of the fact that quantum computers offer the unique chance to ask questions that the founding fathers of quantum mechanics could only formulate as "Gedanken Experiments".

FIG. 1 :
FIG. 1: The four quadrants in the Physics-Computing plane.The CC and CQ quadrants are the mainstay of current scientific computing.QQ is the quadrant invoked by Feynman, and QC is the quadrant relevant to this paper. FIG.

2 FIG. 3 :
FIG. 3: (a) The deviation ϵ in time of the second order Carleman solution from the original LB solution for the case of a Kolmogorov-like flow at Re ∼ 10, on a 32 × 32 grid.(b) The mean deviation ⟨ϵ⟩ as a function of the Reynolds number for the same flow.
The circuit for the collision-free streaming in the D2Q9 model.

FIG. 5 :FIG. 6 :
FIG.5:The sketch of three quantum circuits for different Lattice-Boltzmann implementations for a grid with G lattice sites for the D2Q9 model.In (a) the collision is applied on the velocity register whereas the streaming is a global operation applied on both the position and the velocity registers.In (b) The evolution is obtained as a series of streaming in the 8 velocities of the D2Q9 model.In (c) the quantum register |k⟩ holds the information about the truncation order, whereas |i⟩ and |j⟩ are the velocities associated to the Carleman variables.

TABLE I :
The velocity set of the D2Q9 model.simple to relate together the streaming and the collision steps, as we are going to show next.When writing the Lattice Boltzmann functions f i (x, t) in a quantum state we can employ amplitude encoding, using a quantum register |•⟩ v for embedding the discrete velocities v i and a second quantum register for encoding the lattice sites |•⟩ x .The fluid is then described by the following quantum state