# Evolution of CFD numerical methods and physical models towards a full discrete approach

Jean-Paul Caltagirone

► To cite this version:

Jean-Paul Caltagirone. Evolution of CFD numerical methods and physical models towards a full discrete approach. Comptes Rendus. Mécanique, 2022, 350 (S1), pp.1-10. 10.5802/crmeca.131. hal-04037123

HAL Id: hal-04037123

https://hal.science/hal-04037123

Submitted on 20 Mar 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Comptes Rendus
## Mécanique

Jean-Paul Caltagirone

**Evolution of CFD numerical methods and physical models towards a full
discrete approach**

Published online: 11 October 2022

https://doi.org/10.5802/crmeca.131

**Part of Special Issue:** More than a half century of Computational Fluid Dynamics

**Guest editor:** Mohammed El Ganaoui (Professeur des Universités, Spécialiste :
Mécanique des Fluides et Transferts de Chaleur et de Masse, Université de
Lorraine)

<figure><img src="image_2.png" /></figure>

©

<figure><img src="image_3.png" /></figure>

This article is licensed under the
CREATIVE COMMONS ATTRIBUTION 4.0 INTERNATIONAL LICENSE.
http://creativecommons.org/licenses/by/4.0/

<figure><img src="image_4.png" /></figure>
<figure><img src="image_5.png" /></figure>

More than a half century of Computational Fluid Dynamics / Plus d'un
demi-siècle de mécanique des fluides numérique

Evolution of CFD numerical methods and
physical models towards a full discrete
approach

Jean-Paul Caltagirone® a

a Bordeaux INP, University of Bordeaux, CNRS UMR-5295, Arts et Métiers Institute of Technology, INRAE, I2M Bordeaux, 33400 Talence, France

E-mail: calta@ipb.fr

**Abstract.** The physical models and numerical methodologies of Computational Fluid Dynamics (CFD) are historically linked to the concept of continuous medium and to analysis where continuity, derivation and integration are defined as limits at a point. The first consequence is the need to extend these notions in a multidimensional space by establishing a global inertial frame of reference in order to project the variables there. In recent decades, the emergence of methodologies based on differential geometry or exterior calculus has changed the point of view by starting with the creation of entangled polygonal and polyhedral structures where the variables are located. Mimetic methods and Discrete Exterior Calculus, notably, have intrinsic conservation properties which make them very efficient for solving fluid dynamics equations. The natural extension of this discrete vision relates to the derivation of the equations of mechanics by abandoning the notion of continuous medium. The Galilean frame of reference is replaced by a local frame of reference composed of an oriented segment where the acceleration of the material medium or of a particle is defined. The extension to a higher dimensional space is carried from cause to effect, from one local structure to another. The conservation of acceleration over a segment and the Helmholtz-Hodge decomposition are two essential principles adopted for the derivation of a discrete law of motion. As the fields covered by CFD are increasingly broad, it is natural to return to the deeper meaning of physical phenomena to try a new research or new path which would preserve the properties of current formulations.

**Keywords.** Conservation of acceleration, Navier-Stokes equations, Mimetic methods, Discrete exterior calculus, Discrete mechanics, Helmholtz-Hodge decomposition.

Published online: 11 October 2022

**1. Introduction**

The last five decades of the evolution of mechanics from the theoretical point of view [1-4] have seen the emergence of numerical simulation as a discipline in its own right [5, 6]. Although unknowns remain, in particular on the problem of the existence and the uniqueness of solutions to the Navier-Stokes equations according to the initial conditions given in three dimensions of space, mathematics has opened up prospects for the development of new numerical methods and schemes. Computational Fluid Dynamics (CFD) has been marked by the extraordinary
progress of numerical methodologies and computer resources; they have led to in-depth im-
provements to the knowledge of the physical behavior of fluids and solids subjected to ever more
extensive stresses. However, numerical resources are not the only elements that have allowed the
acceleration of CFD. The development of schemes adapted to the resolution of boundary prob-
lems for incompressible or compressible flows have, for example, facilitated the consideration of
shock waves by limiting numerical diffusion.

As the diversity of the problems treated has increased with the mechanical equations taken
into account—Navier–Lamé, Stokes, Euler, Navier–Stokes—the formulations used needed to be
adapted, for example to approach in monolithic form problems of Fluid–Structure Interaction.
The formulation in velocity–pressure primitive variables (v, p) is often associated with a conser-
vation equation of mass, for example $\nabla \cdot v = 0$ for an incompressible movement, and a state law
to ensure the link between pressure and density $\rho$. The methods for solving partial differential
equations, in particular Navier–Stokes equations, are numerous, with finite differences, finite vol-
umes, finite elements spectral methods, etc. For almost a century, they have constantly evolved,
sometimes towards hybrid approaches with new properties. The adoption of staggered grids with
the Marker And Cell method [7] probably marked a turning point in the correspondence between
the locations of the variables and those of the geometrical structure considered; this has led many
authors to reflect on this duality. The progress of mathematical analysis applied to the resolution
of mechanical equations has also led to the development of specific formulations, such as projec-
tion methods [8, 9] which decouple the pressure and the velocity when the Navier–Stokes equa-
tion is solved. Strong coupling techniques like that of the Augmented Lagrangian [10] make it
possible to ensure incompressibility without introducing artifacts linked to decoupling. Analy-
ses on the conservative properties of these methods have shown their remarkable suitability with
balance equations [11–14].

In addition to solution methods based on partial differential equations, many other tech-
niques are available to access solutions of boundary problems in fluid mechanics. A number of
them mimic physical phenomena at small scales up to the macroscopic level, Molecular Dynam-
ics, Lattice Gas Cellular Automata, etc. The most promising one seems to be the Lattice Boltz-
mann method which is particularly flexible to represent many physical phenomena with a re-
duced modeling [15, 16].

Recently, resolution methods resulting from differential geometry and exterior calculus have given new impetus to CFD, in particular, based on the Discrete Exterior Calculus (DEC) and mimetic methods. Under the impulse of Shashkov [17], mimetic methods were developed to solve the classical equations of mechanics [18, 19], and the Maxwell equations [20–22]. Since then, methods preserving certain quantities, for example a mass, energy, enstrophy and vorticity conserving (MEEVC) [23, 24] in spectral or by summation-by-parts operators [25], have shown that high precision levels could also be obtained. Schemes for solving the Stokes or Navier–Stokes equations have also been developed on these bases, in particular the Compatible Discrete Operator (CDO) [26, 27] method. At the same time, differential geometry and exterior calculus allowed the emergence of new methods of resolution associated with simplicial meshes, successfully applied to incompressible Navier–Stokes equations [28].

The numerical methods and the various formulations of the equations of mechanics pre-
sented succinctly above are all based on the notion of continuous medium and apply to equa-
tions formulated in this context. It is possible to step outside of this historical framework to ap-
proach the derivation of the equations of mechanics from a geometric point of view, which is
the object of discrete mechanics as conceptualized recently [29]. This discrete physical model ex-
presses that the intrinsic acceleration of a material medium or of a particle is equal to the sum
of the accelerations which are applied to it within the geometrical structure. The two main accel-
erations corresponding to the effects of compression and shearing are immediately dissociated
<figure><img src="image_6.png" /><figcaption>**Figure 1.** Geometric structures of discrete mechanics, primal in blue and dual in red. The primal surface S is defined by the edges $\sigma$, whose collection forms the boundary $\Gamma$ of the surface; each edge oriented by the unit vector $\mathbf{t}$ is delimited by two points $a$ and $b$. Similarly, the dual surface $D$ is composed by segments $\delta$, which are associated with the barycenters $c$ of each of the facets S and oriented by the unit vector $\mathbf{n}$ and whose collection forms the boundary $\Delta$. The unit vectors are orthogonal by definition, $\mathbf{n} \cdot \mathbf{t} = 0$. Besides, the direct orientation of circulation along $\Gamma$ or $\Delta$ is on agreement with the Maxwell's corkscrew rule.</figcaption></figure>

within the equation of discrete motion; the corresponding terms within this equation are orthogonal. However, it is possible to revert to the concept of continuous medium if necessary; all terms of the discrete equation of motion can indeed be transformed to find the continuous form of the Navier-Stokes equations [30]. The discrete equation of motion resulting from this formalism leads directly to a numerical methodology of the mimetic or DEC type, insofar as the differential operators used for the physical model correspond to the discrete operators.

## 2. Discrete mechanics framework

### 2.1. *Geometrical structures*

The discrete mechanics developed in other articles [31–34] is condensed here to present the essential points. Geometric construction is the foundation of the physical model; it is formed by two entangled structures, one primal and the other dual, presented in Figure 1.

The acceleration and the velocity are located on oriented segment $\sigma$ of the primal structure of ends $a$ and $b$ and of unit vector $\mathbf{t}$. The primal planar face $S$ is bounded by a collection of segments $\sigma$ whose unit normal is $\mathbf{n}$. The dual face $D$, delimited by segments $\delta$, is also assumed to be planar. The global physical domain is formed by a tessellation of elementary structures similar to those in Figure 1 using the meshes generator.

The vectors or, rather, the components on $\sigma$ of the vectors of space, like acceleration $\mathbf{\gamma}$ and velocity $\mathbf{v}$, are assumed to be constant on this segment. Axial vectors or pseudo-vectors, or rather their components, like $\psi$, are orthogonal to the planes of primal space; the two unit vectors are orthogonal $\mathbf{n} \cdot \mathbf{t} = 0$ by construction. The scalars are assigned to the vertices of the primal structure.
In this context, the discrete differential operators take on meanings that are different from their classical forms defined in a continuous medium. The gradient $\nabla\phi = (\phi_b - \phi_a)/dht$ is the restriction of the classical gradient on the oriented segment $\sigma$. The primal curl $\nabla \times \mathbf{v}$ corresponds to the classical curl projected on the normal $\mathbf{n}$. The divergence of velocity $\nabla \cdot \mathbf{v}$ is calculated on the dual volume and assigned to point $a$. Finally, the dual curl $\nabla^d \times \mathbf{\psi}$ is calculated as the circulation of the vector along contour $\Delta$ formed by elementary segments $\delta$; this operator projects the result onto segment $\sigma$. This last operator does not exist within the framework of the continuous medium; which leads to a differentiation of the vectors to which they apply: axial or polar.

The overall structure exhibits remarkable properties of orthogonality. Some of them mimic those of the continuous medium, $\nabla \cdot \nabla \times \mathbf{\psi} = 0$ and $\nabla \times \nabla\phi = 0$, whether the geometric structures' element bases are polygons or polyhedra.

## 2.2. Physical model

Acceleration, the material derivative of velocity, is considered an absolute quantity in discrete mechanics, unlike other quantities such as velocity or displacement; indeed, the derivatives of higher order have no particular physical significance. The integration of $\gamma$ along its rectilinear trajectory on segment $\sigma$ represents the energy per unit mass necessary to change the velocity of the material medium or of a particle between $a$ and $b$ over a time interval $dt$. Thus the acceleration of the medium subjected to an external acceleration $\mathbf{h}$ is kept on the segment. The fundamental law of dynamics in terms of accelerations is then written as follows:

$$ \gamma = \mathbf{h} \qquad (1) $$

where $\gamma$ is the intrinsic acceleration of the material medium or of the particle with or without mass on a segment, and $\mathbf{h}$ the sum of the accelerations applied to it.

The time lapse $dt$ is a physical parameter which fixes the time scale in which a phenomenon is observed. The current time $t = t^o + dt$ is defined with respect to a reference time $t^o$ corresponding to a previous moment of mechanical equilibrium. Likewise, the length scale must be in line with it. A characteristic length, called discrete horizon $dh = [a, b]$ corresponding to the length of segment $\sigma$, is then defined by $dh = c_l dt$ where $c_l$ is the celerity of the longitudinal waves. Velocity $\mathbf{v}$ and displacement $\mathbf{u}$ at time $t$ are related to time $t^o$ by:

$$ \begin{cases} \mathbf{v} = \mathbf{v}^o + \gamma dt \\ \mathbf{u} = \mathbf{u}^o + \mathbf{v} dt. \end{cases} \qquad (2) $$

The principle of velocity relativity is satisfied when velocity $\mathbf{v}^o$ or displacement $\mathbf{u}^o$ at instant $t^o$ no longer appear in the equation of motion when these are uniform, translated or rotated at a constant velocity.

Like any vector, acceleration can be decomposed into a component with curl-free and another with divergence-free in the form of a decomposition of Helmholtz-Hodge [35-37]:

$$ \gamma = -\nabla\phi + \nabla^d \times \mathbf{\psi} \qquad (3) $$

where $\phi$ and $\psi$ are the scalar and vector potentials of the acceleration respectively.

From a physical point of view this observation is essential; in discrete mechanics the acceleration of the material medium on a segment is due to a direct action represented by the scalar potential linked to the propagation of longitudinal waves and an induced action fixed by the vector potential related to transverse waves. This is the object of Maxwell's remarkable discovery [38], which introduced the entanglement of magnetic and electric fields through their respective dynamics. In mechanics this translates the fact that the effects of compression can in no case be transformed into shearing or vice versa without an evolution in time. This parallel between electromagnetic and mechanical effects is not just a simple analogy, it is the same concept.
Therefore, the variation in total energy per unit mass between the ends of segment $\sigma$ is equal to:

$$ \Phi_b - \Phi_a = \int_a^b \gamma \cdot t dl = - \int_a^b \nabla \phi \cdot t dl + \int_a^b \nabla^d \times \psi \cdot t dl. \quad (4) $$

The physical modeling of current potentials $\phi$ and $\psi$ according to their values at instant $t^o$, $\phi^o$ and $\psi^o$ is carried out based on simple compression and shear experiments [29]. The energy balance on the segment is read as:

$$ \int_{\sigma} \gamma \cdot t dl = - \int_{\sigma} \nabla (\phi^o - dt c_l^2 \nabla \cdot \mathbf{v}) \cdot t dl + \int_{\sigma} \nabla^d \times (\psi^o - dt c_t^2 \nabla \times \mathbf{v}) \cdot t dl \quad (5) $$

where $\phi^o$ and $\psi^o$ are named the retarded potentials, similar to those of Liénard and Wichert [39]:

$$ \phi^o = - \int_0^{t^o} c_l^2 \nabla \cdot \mathbf{v} d\tau; \quad \psi^o = - \int_0^{t^o} c_t^2 \nabla \times \mathbf{v} d\tau. \quad (6) $$

The transition from this one-dimensional vision to two- or three-dimensional spaces is achieved without difficulty, considering that the interactions between the segments of all the entangled primal and dual structures are cause and effect; the same point $a$ of the primal structure is shared by other segments of space. Finally, the system of equations of discrete mechanics in terms of displacements is put in the form:

$$ \left\{ \begin{aligned} \frac{d^2 \mathbf{u}}{dt^2} &= -\nabla(\phi^o - c_l^2 \nabla \cdot \mathbf{u}) + \nabla^d \times (\psi^o - c_t^2 \nabla \times \mathbf{u}) \\ \phi &= (1 - \alpha_l) \phi^o - c_l^2 \nabla \cdot \mathbf{u} \\ \psi &= (1 - \alpha_t) \psi^o - c_t^2 \nabla \times \mathbf{u}. \end{aligned} \right. \quad (7) $$

The longitudinal $c_l$ and transverse $c_t$ celerities are known physical quantities which may depend on potentials or on other variables such as temperature. The quantities $\alpha_l$ and $\alpha_t$ are called the attenuation factors of the longitudinal and transverse waves. They express the dissipation of these waves in the form of heat; the energies per unit of mass $\phi$ and $\psi$ are then reduced by this degradation. As expected, the conservation of mass is not added to this equation of motion insofar as it expresses the conservation of total energy. The equivalence of mass and energy resulting from the theory of relativity rules out maintaining redundant physical quantities.

To express the acceleration as a Helmholtz-Hodge decomposition, it is necessary for all the other terms of (3) to be in this same form. This is not the case for the classical expressions of inertia, $\mathbf{v} \cdot \nabla \mathbf{v}$, $\nabla \cdot (\mathbf{v} \otimes \mathbf{v}) - \mathbf{v} \nabla \cdot \mathbf{v}$ or from $\nabla(|\mathbf{v}|^2/2) - \mathbf{v} \times \nabla \times \mathbf{v}$. Discrete mechanics establishes inertia in the form of two terms [33], the first to curl-free and the second to divergence-free, $\nabla(|\mathbf{v}|^2/2) - \nabla^d \times (|\mathbf{v}|^2/2\mathbf{n})$. The equation of motion in terms of velocity is then written as:

$$ \left\{ \begin{aligned} \frac{\partial \mathbf{v}}{\partial t} + \frac{1}{2} \nabla (|\mathbf{v}|^2) - \frac{1}{2} \nabla^d \times (|\mathbf{v}|^2 \mathbf{n}) &= -\nabla (\phi^o - dt c_l^2 \nabla \cdot \mathbf{v}) + \nabla^d \times (\psi^o - dt c_t^2 \nabla \times \mathbf{v}) \\ \phi &= (1 - \alpha_l) \phi^o - dt c_l^2 \nabla \cdot \mathbf{v} \\ \psi &= (1 - \alpha_t) \psi^o - dt c_t^2 \nabla \times \mathbf{v}. \end{aligned} \right. \quad (8) $$

For a viscous Newtonian fluid $\alpha_t = 1$, the shear energy is entirely dissipated in heat. Indeed, for time constants $\tau > 10^{-10}$ s corresponding to the molecular reorganization under the effect of a shear, quantity $dt c_t^2$ is replaced by the kinematic viscosity $\mathbf{v}$.

In the general case of movements of several fluids, for example for two-phase flows, the interfaces must be advected using adapted methods in Eulerian or Lagrangian description (TVD, Front-Tracking, Volume Of Fluid, Level-Set, ALE, etc.). This problem of advection of fields or interfaces, which also arises for a single fluid, does not interfere with the proposed formulation, it is described by the term $-\mathbf{v} \cdot \nabla \rho$ within of the mass conservation equation. The rest of the continuity equation, $\partial \rho / \partial t = - \rho \nabla \cdot \mathbf{v}$, can be integrated separately. Indeed, irrespective of the
reasons for variations in density, effects of compression, fluxes of mass or heat, this is only a function of the local divergence of the velocity. The density can then be integrated at a point explicitly by:

$$ \rho = \rho^0 e^{-\int \nabla \cdot \mathbf{v} \, dt} \qquad (9) $$

where $\rho^0$ is the density at time $t^0$. Even if this, like pressure, is no longer a variable retained in discrete mechanics, it can be updated like the scalar potential of the acceleration. It again becomes a physical characteristic of the environment, which may be dependent on the main variables $(\phi, \mathbf{v})$.

In discrete equation (8), it is possible to integrate source terms resulting from gravitational, capillary or other interactions, in the form of accelerations. Each of these terms will be written as a Helmholtz-Hodge decomposition, preserving the unity of the formulation of the equation.

In form (8) or that of (7), the law of motion is like Maxwell: an equation of waves. Only the celerity of light is replaced in mechanics by the celerity of acoustic waves. To be convinced of this, it is sufficient to consider a celerity $c_l = c_t = c$, neglect the inertia terms and use the formula of vector calculus $\nabla^2 \mathbf{u} = \nabla \nabla \cdot \mathbf{u} - \nabla \times \nabla \times \mathbf{u}$ to get:

$$ \frac{\partial^2 \mathbf{u}}{\partial t^2} - c^2 \nabla^2 \mathbf{u} = -\nabla \phi^0 + \nabla^d \times \boldsymbol{\psi}^0 \qquad (10) $$

where the first member of this equation is a d'Alembertian $\Box \mathbf{u}$:

$$ \Box \mathbf{u} \equiv \frac{1}{c^2} \frac{\partial^2 \mathbf{u}}{\partial t^2} - \nabla^2 \mathbf{u} \qquad (11) $$

where the second member corresponds to the Helmholtz-Hodge decomposition of the retarded acceleration.

The boundary conditions are not expressed as continuum mechanics by direction using partial derivatives on the velocity. In discrete mechanics they are the same operators as those used in the equation of motion, which translate the interactions from one local geometrical structure to another. If $\sigma_p$ represents a jump in scalar potential, for example due to capillary pressure, and $\tau_p$ a jump in vector potential, the general boundary conditions or those at the interfaces (fluid and solid) are written as:

$$ \begin{cases} dt c_{l1}^2 \nabla \cdot \mathbf{v}_1 = dt c_{l2}^2 \nabla \cdot \mathbf{v}_2 + \sigma_p \\ dt c_{t1}^2 \nabla \times \mathbf{v}_1 = dt c_{t2}^2 \nabla \times \mathbf{v}_2 + \tau_p. \end{cases} \qquad (12) $$

The equation of motion (8) has particular properties due to its Helmholtz-Hodge type configuration. In particular, the application of the divergence and curl operators directly to the equation of motion eliminates existing terms in classical formulations. They reflect the observation that a compression effect cannot generate shearing without assigning it a dynamic way. If the time derivative is zero, $\partial \mathbf{v} / \partial t = 0$, then it is possible to separate the equation of motion (8) into two components:

$$ \begin{cases} \nabla(\phi_B^o - dt c_l^2 \nabla \cdot \mathbf{v}) = \mathbf{h} \\ \nabla^d \times (\psi_B^o - dt c_t^2 \nabla \times \mathbf{v}) = \mathbf{h} \end{cases} \qquad (13) $$

where $\phi_B^o = \phi^o + |\mathbf{v}|^2/2$ and $\psi_B^o = \psi^o + |\mathbf{v}|^2/2 \mathbf{n}$ are the Bernoulli scalar and vector potentials. Vector $\mathbf{h}$ is a harmonic term with both divergence-free and curl-free.

It should be noted that (7) is composed of two Lagrangians, the first one relating to the effects of compression and the second associated with the effects of shearing. Each of the Lagrangians consists of an oscillator, $\phi^o$ or $\psi^o$, setting aside or releasing the energy per unit mass over time. These particular forms of Lagrangian or Hamiltonian equations give them properties of symmetry and conservation according to the theorems of Noether [40, 41]. Each invariance in time or space, in translation or rotation, corresponds to the conservation of a fundamental
quantity: energy, angular momentum, etc. The conservation of these quantities is of course essential for the laws of mechanics and the numerical methodologies which make it possible to solve them.

## 2.3. *Numerical methodology*

The distinction between the physical model and the associated numerical methodology is very minimal in discrete mechanics. The step of spatial discretization of continuous equations to form an algebraic system is not necessary; the discrete operators are those used to establish the algebraic system on the principal unknown, the component of velocity $v$ on each edge. The numerical methodology already described with precision [34] is similar to mimetic methods; applying the gradient, primal curl, dual-curl and divergence operators to the scalar and vector potentials reprojects the results as accelerations on segment $\sigma$. Integration in time makes it possible to determine the solution at the current time according to the solution at the previous times.

As shown by the results of the simulations in many situations with a simple fluid or in two-phase flows, whether unsteady, incompressible or not, in two or three-dimensional spaces, the solutions are of order two in space and time [31–34].

# 3. Application to few laws of mechanics

The mechanics of continuous media, supposed to unite solid mechanics and fluid mechanics, has not completely achieved its objective. This observation does not in any way affect the validity of each equation—Euler, Stokes, Navier-Lamé, Navier-Stokes—or their representativeness or the results of the simulations carried out with each of them. Solid and structural mechanics, fluid mechanics or the propagation of nonlinear waves still appear to be related fields. Even if the Navier-Stokes equation makes it possible to comprehend certain problems in linear elasticity or of wave propagation, essential differences remain. For example, Navier-Stokes equations are made up of a vector equation on the momentum and an adjoint equation, the conservation of mass, while the Navier-Lamé equation is autonomous. The links between moduli of elasticity and the Lamé coefficients for fluids have not been established satisfactorily and the controversies of the last century [42] on this point are not over. In fact, Stokes' law is incorrect [43, 44] and the bulk viscosity is not zero. The impact of the confusion between compressive viscosity and the compressibility of a fluid is greatly reduced thanks to the overlap due to the association of the continuity equation. Other paradoxes or biases present in these equations are due to the nature attributed to them, which is independent of the time frame in which the observation of the physical system is carried out.

Table 1 summarizes the expressions of the discrete equations for a number of equations of continuum mechanics.

The Stokes equation $-\nabla p + \mu \nabla^2 v = 0$ can be reformulated without the scalar potential $-\nabla(r\nabla \cdot v) - \nabla^d \times (v\nabla \times v) = 0$ where the numerical parameter $r$ ensures incompressibility. For the Euler equation, the conservation of mass is not required and the state law is integrated within the scalar potential. For example, for an isentropic transformation, the pressure is deduced *a posteriori* from $p = \rho\phi^o$. The incompressible Bernoulli law can be completed here for a compressible unidirectional flow. The generalized law of Darcy is obtained by associating with the term of compression a term of linear drag function of permeability $K$ of the porous medium. The velocity which appears in the diffusion equation represents a flux of mass or heat and its divergence serves to update the scalar potential associated with the concentration or with the temperature.
<table><caption>**Table 1.** Discrete expressions of some equations of continuum mechanics where $v$ is the kinematic viscosity, $\mu$ the dynamic viscosity, $\rho$ the density, $p$ the pressure, $\phi_g$ the gravitational potential, $K$ the permeability of porous medium, $a$ the heat diffusivity, $\kappa$ a dissipation coefficient and $r$ a numerical parameter imposing incompressibility</caption><thead><tr><th>Equation</th><th>Application</th></tr></thead><tbody><tr><td>Navier-Stokes</td><td>$\frac{d\mathbf{v}}{dt} = -\nabla(\phi^o - dt c_l^2 \nabla \cdot \mathbf{v}) - \nabla^d \times (v \nabla \times \mathbf{v})$</td></tr><tr><td>Navier-Lamé</td><td>$\frac{d^2\mathbf{u}}{dt^2} = -\nabla(\phi^o - c_l^2 \nabla \cdot \mathbf{u}) + \nabla^d \times (\psi^o - c_t^2 \nabla \times \mathbf{u})$</td></tr><tr><td>Stokes</td><td>$-\nabla(r \nabla \cdot \mathbf{v}) - \nabla^d \times (\frac{\mu}{\rho} \nabla \times \mathbf{v}) = 0$</td></tr><tr><td>Euler</td><td>$\frac{d\mathbf{v}}{dt} = -\nabla(\phi^o - dt c_l^2 \nabla \cdot \mathbf{v}), \phi^o = \frac{p}{\rho}$</td></tr><tr><td>Bernoulli</td><td>$-\nabla(\phi^o + \frac{|\mathbf{v}|^2}{2} - dt c_l^2 \nabla \cdot \mathbf{v} + \nabla \phi_g) = 0$</td></tr><tr><td>Darcy</td><td>$\frac{\partial \mathbf{v}}{\partial t} = -\nabla(\phi^o - dt c_l^2 \nabla \cdot \mathbf{v}) - \frac{v}{K} \mathbf{v}$</td></tr><tr><td>Diffusion</td><td>$\frac{\partial \mathbf{v}}{\partial t} = -\nabla(\phi^o - a \nabla \cdot \mathbf{v}) - \kappa \mathbf{v}$</td></tr><tr><td>Waves</td><td>$\frac{d^2\mathbf{u}}{dt^2} - c^2 \nabla^2 \mathbf{u} = -\nabla \phi^o + \nabla^d \times \psi^o$</td></tr></tbody></table>

The adaptation of discrete equation (8) to other areas of field theory is possible by abandoning the classical physical variables and parameters and defining only the physical parameters of this equation; for example, thermal diffusion at all space and time scales [32] can thus be represented using the discrete model. Given the origin of the discrete physical model, which is very similar to the concepts of electromagnetism, its application to this framework is natural.

## 4. Outlook

Contrary to the numerical methodologies of CFD which progress regularly, the laws of mechanics have evolved little over the centuries; this observation is probably the mark of the excellent representativeness of these equations. As always, the new laws and theories do not call into question the more traditional laws established previously, validated by innumerable observations; instead they prolong them by widening, for example, their field of validity. The differences then become much more subtle, of the second order. They can be completely negligible in simple cases and considered as simple biases or paradoxes of classical equations. The use of these equations to predict solutions in complex or extreme situations can then give rise to real failures. The equations of mechanics and the numerical methodologies to solve them, in particular for CFD, become more and more inseparable in the perspective of solving problems of interactions of fluids, solids and waves.

This evolution towards discrete equations in no way alters their capacities to represent physical phenomena at all spatial scales greater than the mean free path of molecules. The reduction towards a point for continuum mechanics is replaced in discrete mechanics by a homothetic reduction of the primal and dual structures by decreasing discrete horizon $dh$. Discrete equations thus preserve the continuous character of classical mechanics. In addition, the absence of
density within the equation of motion explicitly makes it possible to define each of its terms as an acceleration, and this applies equally well to a material medium as to particles with or without mass. Finally, all the quantities, variables or physical properties can only be expressed with two fundamental units, a length and a time.

## Conflicts of interest

The authors declare no competing financial interest.

## References

[1] L. Landau, E. Lifshitz, *Fluid Mechanics*, Pergamon Press, London, 1959.

[2] C. Truesdell, W. Noll, *The Non-Linear Field Theories of Mechanics*, Encyclopedia of Physics, vol. III/3, Springer-Verlag, New York, 1992.

[3] P. Germain, P. Muller, *Introduction à la Mécanique des Milieux Continus*, 2ᵉ ed., Masson, Paris, 1995.

[4] J. Salencon, *Mécanique des milieux continus*, Editions de l'Ecole Polytechnique, Palaiseau, 2002.

[5] A. Chorin, "The numerical solution of the Navier-Stokes equations for an incompressible fluid", *Bull. Am. Math. Soc.* **73** (1967), p. 928-931.

[6] R. Temam, *Navier-Stokes Equations, Theory and Numerical Analysis*, North-Holland, Amsterdam, New-York, 1984.

[7] E. Harlow, J. Welch, "Numerical calculation of time-dependent viscous incompressible flow of fluid with a free surface", *Phys. Fluids* **8** (1965), p. 2182-2189.

[8] A. Chorin, "A numerical method for solving incompressible viscous flow problems", *J. Comput. Phys.* **135** (1967), p. 118-125.

[9] J. Guermond, P. Minev, J. Shen, "An overview of projection methods for incompressible flows", *Comput. Methods Appl. Mech. Eng.* **195** (2006), p. 6011-6045.

[10] M. Fortin, R. Glowinski, *Méthodes de lagrangien augmenté, Application à la résolution numérique de problèmes aux limites*, Dunod, Paris, 1982.

[11] R. Nicolaides, "Flow discretization by complementary volume techniques", in *9th AIAA Computational Fluid Dynamics Conference*, vol. 1, 1989, p. 464-470.

[12] B. Perot, "Conservation properties of unstructured staggered mesh schemes", *J. Comput. Phys.* **159** (2000), p. 58-89.

[13] B. Perot, V. Subramanian, "Discrete calculus methods for diffusion", *J. Comput. Phys.* **224** (2007), p. 59-81.

[14] B. Perot, "Discrete conservation properties of unstructured mesh schemes", *Annu. Rev. Fluid Mech.* **43** (2011), p. 299-318.

[15] G. R. McNamara, G. Zanetti, "Use of the Boltzmann equation to simulate lattice-gas automata", *Phys. Rev. Lett.* **61** (1988), p. 2332-2335.

[16] S. Succi, *The Lattice Boltzmann Equation for Fluid Dynamics and Beyond*, Clarendon Press, Oxford, 2001.

[17] M. Shashkov, *Conservative Finite-Difference Methods on General Grids*, CRC Press, Boca Raton, 1996.

[18] J. M. Hyman, M. Shashkov, "Natural discretizations for the divergence, gradient ans curl on logically rectangular grids", *SIAM J. Numer. Anal.* **36** (1999), p. 788-818.

[19] J. M. Hyman, M. Shashkov, "The orthogonal decomposition theorems for mimetic finite difference methods", *SIAM J. Numer. Anal.* **36** (1999), p. 788-818.

[20] J. M. Hyman, M. Shashkov, "Mimetic discretizations for Maxwell's equations", *J. Comput. Phys.* **151** (1999), p. 881-909.

[21] J. M. Hyman, M. Shashkov, "Mimetic discretizations for Maxwell's equations and the equations of magnetic diffusion", *Prog. Electromagn. Res.* **32** (2001), p. 89-121.

[22] K. Lipnikov, G. Manzini, M. Shashkov, "Mimetic finite difference method", *J. Comput. Phys.* **257** (2014), p. 1163-1227.

[23] A. Palha, M. Gerritsma, "A mass, energy, enstrophy and vorticity conserving (MEEVC) mimetic spectral element discretization for the 2D incompressible Navier-Stokes equations", *J. Comput. Phys.* **328** (2017), p. 200-220.

[24] G. de Diegoa, A. Palha, M. Gerritsma, "Inclusion of no-slip boundary conditions in the MEEVC scheme", *J. Comput. Phys.* **378** (2019), p. 615-633.

[25] H. Ranocha, K. Ostaszewski, P. Heinisch, "Discrete vector calculus and Helmholtz Hodge decomposition for classical finite difference summation by parts operators", *Commun. Appl. Math. Comput.* **2** (2020), p. 581-611.

[26] J. Bonelle, A. Ern, "Analysis of compatible discrete operator schemes for the stokes equations on polyhedral meshes", *IMA J. Numer. Anal.* **35** (2015), p. 1672-1697.

[27] J. Bonelle, D. Di Pietro, A. Ern, "Low-order reconstruction operators on polyhedral meshes: application to compatible discrete operator schemes", *Comput. Aided Geom. Des.* **35** (2015), p. 27-41.
[28] M. Mohamed, A. Hirani, R. Samtaney, "Discrete exterior calculus discretization of incompressible Navier-Stokes equations over surface simplicial meshes", *J. Comput. Phys.* **312** (2016), p. 175-191.

[29] J. P. Caltagirone, *Discrete Mechanics, Concepts and Applications*, ISTE, John Wiley & Sons, London, 2019.

[30] J. P. Caltagirone, "On a reformulation of Navier-Stokes equations based on Helmholtz-Hodge decomposition", *Phys. Fluids* **33** (2021), article no. 063605.

[31] J. P. Caltagirone, S. Vincent, "On primitive formulation in fluid mechanics and fluid-structure interaction with constant piecewise properties in velocity-potentials of acceleration", *Acta Mech.* **231** (2020), p. 2155-2171.

[32] J. P. Caltagirone, "Non-Fourier heat transfer at small scales of time and space", *Int. J. Heat Mass Transfer* **160** (2020), article no. 120145.

[33] J. P. Caltagirone, "On Helmholtz-Hodge decomposition of inertia on a discrete local frame of reference", *Phys. Fluids* **32** (2020), article no. 083604.

[34] J. P. Caltagirone, "Application of discrete mechanics model to jump conditions in two-phase flows", *J. Comput. Phys.* **432** (2021), article no. 110151.

[35] P. Angot, J. P. Caltagirone, P. Fabrie, "Fast discrete Helmholtz-Hodge decomposition in bounded domains", *Appl. Math. Lett.* **26** (2013), p. 445-451.

[36] H. Bhatia, G. Norgard, V. Pascucci, P. Bremer, "The Helmholtz-Hodge decomposition—a survey", *IEEE Trans. Visual. Comput. Graph.* **19** (2013), p. 1386-1404.

[37] H. Ranocha, P. Öffner, T. Sonai, "Summation-by-parts operators for correction procedure via reconstruction", *J. Comput. Phys.* **311** (2016), p. 299-328.

[38] J. Maxwell, "A dynamical theory of the electromagnetic field", *Phil. Trans. R. Soc. Lond.* **155** (1865), p. 459-512.

[39] A. Liénard, "Champ électrique et magnétique produit par une charge électrique concentrée en un point et animée d'un mouvement quelconque", *L'Éclairage électrique XVI* (1898), p. 5-112.

[40] Y. Kosmann-Schwarzbach, *Noether Theorems. Invariance and Conservations Laws in the Twentieth Century*, Springer-Verlag, New York, 2011.

[41] G. Sardanashvily, *Noether's Theorems. Applications in Mechanics and Field Theory*, Lomonosov Moscow State University, Moscow, 2016.

[42] L. Rosenhead, "Introduction. The second coefficient of viscosity: a brief review of fundamentals", *Proc. R. Soc. Lond. A: Math. Phys. Sci.* **226** (1954), p. 1-6.

[43] M. Gad-El-Hak, "Stokes hypothesis for a newtonian, isotropic fluid", *J. Fluids Eng.* **117** (1995), p. 3-5.

[44] K. Rajagopal, "A new development and interpretation of the Navier-Stokes fluid which reveals why the "Stokes assumption" is inapt", *Int. J. Non-Linear Mech.* **50** (2013), p. 141-151.