A Generalized Model for Predicting the Drag Coefficient of Arbitrary Bluff Shaped Bodies at High Reynolds Numbers
1,2]Yousef M. F. El Hasadi [Email address for correspondence<EMAIL_ADDRESS>1] Johan T. Padding
[1]Process and Energy department, Delft University of Technology,Leeghwaterstraat 39,2628 CB Delft, the Netherlands
[2]Eindhoven University of Technology, Mechanical Engineering Department, The Netherlands

We propose an accurate model for the drag coefficient of arbitrary bluff bodies that is valid for high Reynolds numbers ($Re$). The model is based on the drag coefficient model derived for the case of a sphere:, $C_D = a_1 +{\frac{K a_2}{Re}} +{a_3\log(Re)+ a_4\log^2(Re) + a_5\log^4(Re)}$ (El Hasadi and Padding, Chemical Engineering Science, Vol. 265, 2023). The coefficients $a_2$, $a_3$, $a_4$, and $a_5$ do not depend on the object's shape or its orientation with respect to the flow, and $K$ is the Stokes drag correction factor, which for the case of the sphere, is equal to 1.0.
The shape and orientation effects are included in the value of $a_1$ for the high Reynolds number flow regime. Interestingly, we found a strong correlation between the value of the $a_1$ coefficient and the frictional drag derived from boundary layer theory. One of the main findings of this investigation is that the rate of change of the drag coefficient with respect to the Reynolds number in the inertial flow regime is independent of the shape of the body or its orientation. Our model successfully predicts, with acceptable accuracy, the historical data of Wieselsberger (Technical Report, 1922) for the case of an infinite cylinder. Additionally, the model predicts the drag coefficient of other bluff body geometries such as oblate and prolate spheroids, spherocylinders, cubes, normal flat plates and irregular non-spherical particles.

§ INTRODUCTION

The drag force, which is the resistance force that occurs when a body moves through a fluid, is a critical parameter in various fields of science and engineering. Organisms from tiny bacteria to massive whales have evolved adaptations to minimize the drag force they experience when swimming through fluids [1, 2, 3, 4, 5]. In the transportation industry, drag force is a crucial optimization parameter for reducing fuel consumption in vehicles, aeroplanes, and ships [6, 7]. In the military sector, drag force is vital in designing projectiles, submarines, and unmanned aerial vehicles [8, 9]. Sustainable energy applications such as wind turbines and fluidized beds also rely on understanding and controlling drag force [10, 11]. Fluid drag is one of nature's most vital forces, yet we have a limited ability to predict it, as we will describe in the following paragraphs.

One of the earliest systems of partial differential equations derived were those of Euler [12] in 1752, which described the motion of an ideal fluid (inviscid incompressible flow). However, shortly after their introduction, d' Alembert [12] found that fluid resistance on an object moving in a static fluid is zero after applying Euler's equations, which contradicts the behaviour of real fluids. This statement made by d'Alembert is now widely known as d'Alembert's paradox, which was a great surprise, since Euler's equations are based on a rigorous mathematical framework of  Newton's laws.
With the advancement of the molecular theory of matter, mainly that of gases, Navier and Stokes independently derived a new set of equations that better resembled fluid behaviour, known as the Navier-Stokes equations [12]. These equations introduced a shearing stress term that depends on the fluid's viscosity. In the meantime experiments in the nineteenth century showed that fluids with finite but small viscosity that flow over bluff bodies such as spheres and cylinders are unstable, especially at the rear of the body where eddies are formed. Also, the experiments  verified that the drag force depends on the quadratically on the flow velocity, as Newton suggested [12]. The instability of the fluid at the rear side of the bluff bodies were very difficult to  explain by either Euler's, or the  Navier-Stokes equations. The  non-linearity  of the
Navier-Stokes makes it impossible to find a generic analytical  solution. This  inability of the theory to match the experiments led to a division between theoreticians and empiricists, especially when it comes to  predicting the drag force of an object. This division continued for the whole period of the nineteenth century.

One of the first triumphs of the Navier-Stokes equations was when Stokes [13] found an analytical expression for the drag force experienced by a sphere in slow (creeping) flow conditions and showed that the drag force scales linearly  with the velocity. However, for high-flow velocities, theoreticians attempted to develop new approaches to tackle the discrepancy in the prediction of the drag force. Examples of those are Rayleigh's dead water theory and Poncelet and Saint-Venant's eddy resistance theories, but with no success[12].

A partial resolution of d'Alembert's paradox came at the beginning of the twentieth century [14] by Prandtl [15] using the idea of the boundary layer, which states that the effects of viscosity are confined to a boundary layer near the wall. In contrast, beyond the  the boundary layer, the flow is governed by inviscid flow theory. Prandtl [15] suggested that suction pressure generated at the rear part of the bluff body contributes significantly to the drag force, and it is formed  due to the separation of the boundary layer from the wall.
Prandtl's ideas are based on the boundary layer equations, a simplified  version of the Navier-Stokes equations. Prandtl used the non-slip boundary condition to describe the fluid velocity at the wall. The  mathematical correctness of this boundary condition  has yet to be proven [16]. However, substantial evidence from experiments suggests that the velocity of the fluid attached to a static wall is zero [16], which supports the concept of the non-slip boundary condition. Blasius [17] solved the boundary layer equations for the case of a thin flat plate placed horizontally to the flow direction and came up with an estimate of the drag force for high-velocity laminar flow conditions.
Despite the usefulness of the boundary layer theory, it has its limitations. For example, it is not valid for separated flows, such as those that occur for flow over bluff bodies, and it consists of several paradoxes, as mentioned by Birkhoff [18]. However, the boundary layer theory led to significant progress in understanding fluid mechanics in general. Also, it helped significantly to close the gap between experimental observations and theoretical predictions, and above all, it agglomerated the theoretical and experimental societies together.

Wieselsberger [19] performed experiments with cylinders that were oriented at right angles to the fluid flow direction. He was one of the first researchers to propose using the drag coefficient to represent the drag force. The drag coefficient is defined as $C_D = F_D/\frac{1}{2}\rho v^2_{\infty} A_p$, where $F_D$ is the drag force, $\rho$ is the density of the fluid, $v_{\infty}$ is the undisturbed velocity of the flow, and $A_p$ is the projected area. Also, he showed that the drag coefficient $C_D$ is strongly dependent on   the Reynolds number that can be defined as $Re = \dfrac{\rho v_{\infty} d}{\mu}$, where $d$ is a length scale related to the object, for the case of the cylinder is the diameter, and $\mu$ is the viscosity of the fluid. Wieselsberger [19] found that the cylinder's drag coefficient depends in a complex non-expected way on the Reynolds number  $Re$, and the flow behaviour over a cylinder can take on a range of different morphologies. Churchill's [20] further elucidates this spectrum of flow patterns, which can range from symmetrical and non-separated at low Reynolds numbers to separated and recirculating at moderate Reynolds numbers to shedding vortexes at higher Reynolds numbers. When the Reynolds number is around $3\times 10^5$, the boundary layer on the surface of the cylinder changes from laminar to turbulent. The cascade of the different  flow morphologies significantly effects the evolution   of the drag coefficient with the Reynolds number. Several other experiments for different bluff body geometries, such as sphere and the cube [21, 22] showed the same dependence of the drag coefficient on the Reynolds number.

Only a few analytical solutions are available for the drag coefficient, such as for spherical and cylindrical geometries [23, 24]. However, these solutions are only valid for Reynolds numbers $<$ 1.0. At higher Reynolds numbers, analytical solutions of the Navier-Stokes equations no longer exist due to their non-linearity. In these cases, the drag coefficient is typically evaluated using  boundary layer theory [25, 26] or by solving the Navier-Stokes equations numerically [27, 28, 29, 30]. The drag coefficient results from the numerical simulations are typically expressed in curve-fitted equations, which we will discuss in more detail later on.

In real-life applications, the correlations used for the drag coefficient are predominantly based on experimental data rather than theory. This is because obtaining theoretical results, even from numerical simulations, can be prohibitively computationally  expensive, particularly at high Reynolds numbers. Despite several decades of using boundary-layer theory and CFD simulations, our understanding of the drag force experienced by a bluff body has remained the same since Proudman [31] expressed his concerns about the field's lack of progress and made the following statement  nearly 60 years ago:  “It is a sobering thought that, despite the fifty or so years that the techniques of boundary-layer theory have been at our disposal, so little progress has been made with the problem that even such a bulk characteristic as the dependence of drag coefficient on (large) Reynolds number is still a matter of pure conjecture.".
Several correlations exist for the drag coefficient of arbitrarily shaped particles, but they are only valid for moderate Reynolds numbers of up to 100. One of the earliest correlations proposed for multi-shaped particles is that of Haider and Levenspiel [32]. They introduced the effects of the shape by introducing the drag coefficient dependence on the sphericity. The mathematical form of their correlation is similar to that of the Clift [33] correlation, developed for the case of sphere geometry at high Reynolds numbers.

Recently, we employed  a physics-informed symbolic regression machine learning algorithm that uses sparse volumes of data and incorporates knowledge such as analytical solutions in the training process [34]. We applied this algorithm to obtain mathematical models that describe the drag coefficient of a prolate spheroid and sphero-cylinder particles for different particle orientations and aspect ratios for $Re$ less than 2000, using only 72 data points for training [34]. For the geometries we considered, the drag coefficient is best described by logarithmic terms of the Reynolds number.

We used our machine learning algorithm to verify the conjecture made by Proudman and Pearson [23], which states that the drag coefficient of a sphere can be expressed as an expansion of logarithmic terms and powers of the Reynolds number. Our findings showed that the drag coefficient for a sphere from the creeping flow regime up to the drag crisis ($Re \simeq 10^5$) can indeed be represented by an expansion of logarithmic terms of the Reynolds number . We trained our algorithm using data from the Brown and Lawler [35] correlation for Reynolds numbers ranging from 0.1 to $10^5$. Our logarithmic-based equation was able to extrapolate beyond its training range. For example, it converged to the Oseen constant as $Re \to$ 0 and accurately predicted the onset of the drag crisis for high $Re$. Furthermore, we demonstrated that our algorithm could obtain analytical solutions from data based on experimental correlations. For instance, the equation we derived for the Nusselt number in the case of forced convection over a sphere converged to the analytical solution of Acrivos and Goddard [36] in the limit as $Pe$ (Peclet number) $\to \infty$.

We also found a strong connection between the equations we derived for the drag coefficient of both spherical and non-spherical particles [37, 34]. Regardless of the particle shape, they all contained logarithmic terms, and the coefficients of these terms had very close values. This implies that a general model governing the variation of the drag coefficient could possess logarithmic terms of the Reynolds number.

The main focus of this paper is to develop a drag coefficient model for non-spherical geometries  based on the derived drag coefficient for a sphere, incorporating logarithmic terms of the Reynolds number. The paper will be divided into several sections. Firstly, we will discuss and analyze various drag coefficient formulas available in the literature, providing readers with an understanding of the historical evolution of these formulas. Subsequently, we will delve into the theory underlying our drag coefficient model and extensively validate it against different geometries documented in the literature.

§ CORRELATIONS FOR THE DRAG COEFFICIENT

This section will summarize and describe the main correlations available for the drag coefficient for different non-spherical shapes. We will emphasize the general correlations that can be applied to different particle shapes. Also, we will cover correlations that are specific to the geometries, such as that of oblate and prolate spheroids, since for those geometries there are vast numbers of correlations that are based on data from CFD simulations, which will allow us to compare our results directly with the approximate solution of the Navier-Stokes equations.

§.§ Correlations for arbitrary particle shapes

There are few equations in the literature for the drag variation of an arbitrary particle shape; usually, those correlations are derived using either experimental or numerical available data. For all the correlations from the different sources from the literature that we will present for different particle shapes, the Reynolds number is defined based on the volume-equivalent sphere diameter, and the frontal area is based on a frontal area of a sphere with the same diameter. This type of definition of the Reynolds number and the drag coefficient is widely used in chemical engineering applications. It will not be a problem comparing our results with the available correlations, as we will show later. For this reason, we did not change the symbols of the drag coefficient and the Reynolds number in the upcoming equations.

As we mentioned in the introduction, Haider and Levenspiel [32] were among the first to suggest a general correlation for the drag for non-spherical particles. By fitting a comprehensive variety of data from different experiments for different particle geometries, including spheres, cube octahedrons, octahedrons, cubes, tetrahedrons, and disks. The parameter they used to describe the non-spherical particle geometry variation is the sphericity  ($\phi$). Sphericity is defined as the ratio between the surface area of the sphere having the same volume as the particle, divided by the actual surface area of the particle. By definition, the value of  $\phi$ is  $\leq$ 1.0,  the Haider and Levenspiel [32] correlation is given by the following formula:

\begin{equation}\label{eqpowg1}
C_D = \frac{24}{Re}\left(1+ARe^B\right)+\frac{C}{1+\frac{D}{Re}}
\end{equation}

The constants in Eq.(<ref>) are given as the follows:

\begin{equation}\label{eqpowg2}
\begin{aligned}
A=& \exp \left(2.3288-6.4581 \phi+2.4486 \phi^{2}\right) \\
B=& 0.0964+0.5565 \phi \\
C=& \exp \left(4.905-13.8944 \phi+18.4222 \phi^{2}
-10.2599 \phi^{3}\right) \\
D=& \exp \left(1.4681+12.2584 \phi-20.7322 \phi^{2}+15.8855 \phi^{3}\right)
\end{aligned}
\end{equation}

The correlation of Haider and Levenspiel [32] shows that the particle geometry influences the drag coefficient  from the  low to high Reynolds numbers since all constants from $A$ to $D$ are functions of $\phi$. The most accurate  general correlation for predicting the drag coefficient of non-spherical particles is that of Holzer and Sommerfeld [38], it is based on the same philosophy as that of  Haider and Levenspiel [32] by fitting the experimental drag coefficient values for different particle shapes. Compared to the Haider and Levenspiel [32], they added crosswise sphericity($\psi$), which is defined  as the ratio between the cross-sectional area of the volume equivalent sphere and the projected cross-sectional area of the considered particle perpendicular to the flow. Crosswise sphericity accounts for the effects of the orientation of the flow field on the drag coefficient. At the same time, it significantly increases the complexity of the resulted correlation. The Holzer and Sommerfeld correlation in its simple form is given by the following :

\begin{equation}\label{eqpowg3}
C_D = \frac{8}{Re}\frac{1}{\sqrt{\psi}}+\frac{16}{Re} \frac{1}{\sqrt{\phi}}+\frac{3}{\sqrt{Re}}\frac{1}{\phi^{\frac{3}{4}}}+0.42\times 10^{0.4\left(-\log(\phi)\right)^{0.2}}\frac{1}{\psi}
\end{equation}

§.§ Oblate and prolate spheroids

A few correlations deal with predicting the drag coefficient for oblate spheroids, especially at  high Reynolds numbers. The correlation provided by Ouchene[30]  is applicable for a wide range of aspect ratios that spans from 0.2 to 1.0 and Reynolds numbers from the Stokes regime up to 100,  for case that the flow is parallel to the minor axes of the oblate spheroid  the orientation $\alpha$ = $0^{\circ}$, while if the flow is perpendicular to the minor axes of the oblate spheroid  $\alpha$ = ${90^{\circ}}$. The correlation takes the following form for the case of $\alpha = 90^{\circ}$:

\begin{equation}\label{eqpowo1}
C_{D_{90}} = \frac{24}{Re}\left( K_{o_{90}} +0.15p_a^{100.7}Re^{0.687}+0.14\left(1-p_a^{24.75}\right)\right) Re^{0.7143}
\end{equation}

Where  $K_{o_{90}}$ is the correction of the Stokes drag for the oblate spheroid derived by Happel [39] for $\alpha$ = $90^{\circ}$ given by the following relation:

\begin{equation}\label{eq2olog}
K_{o_{90}} = \frac{8}{3} p_a^{-\frac{1}{3}}\left[-\frac{p_a}{1-p_a^2}-\frac{2p_a^2+3}{(1-pa^2)^{\frac{3}{2}}}\sin^{-1}(1-p_a^2)\right]^{-1}
\end{equation}

Eq.(<ref>) shows that $C_D$ is dependent on the $Re$, and $p_a$. Ouchene[30] also provided a correlation for the oblate spheroids for $\alpha = 0^{\circ}$. However, In order to keep the paper concise and short, we did not include it.

There is  a correlation that is based on experimental data provided by Clift [33] for oblate spheroids with aspect ratio $p_a$ =0.5, for $\alpha$ = $90^{\circ}$ valid for a wide range of Reynolds numbers ranging  between 40  and  $10^4$.  The mathematical  formula of the  Clift[33] correlation is the following:

\begin{equation}\label{eqpowo3}
\log_{10}(C_{D_{90}}) = 2.035-1.660w+0.3958w^2-0.0306w^3
\end{equation}

where $w = \log_{10}(Re_s)$, and $Re_s = Re/(2p^{1/3}_a) $.

Example of drag coefficient correlation for prolate spheroids is given in our previous work   [34] based on the numerical  data of Sanjeevi et al.[40] for the case of $p_a$  = 2.5 and $\alpha$ = $0^{\circ}$, it takes the following form:

\begin{equation}\label{eq1plog}
C_{D_0} = a_1+ \frac{a_2 K_{p{0}} + a_3+ a_4\log^2(Re)}{Re}+a_5\log(Re)+a_6\log^2(Re)
\end{equation}

The values of the coefficients of Eq.(<ref>) are listed in Table <ref>. Since the drag coefficient and the Reynolds number of the Sanjeevi et al.[40] data  are based on the volume equivalent of a sphere, the coefficients of the logarithmic equation (Eq.(<ref>)) are also based on the volume equivalent of a sphere as well. The orientation angle $\alpha$ for the case of prolate spheroids is the following $\alpha$ = $0^{\circ}$ when the flow direction of the flow is in a parallel direction to the major axis of the spheroid, and  $\alpha$ = $90^{\circ}$, when the flow direction is perpendicular to the major axis of the spheroid.

to 0.8 | X[c] | X[c] |
Coefficients      Eq.(<ref>)

$a_1$      3.151

$a_2$      25.873

$a_3$      -2.258

$a_4$       0.230

$a_5$       -0.841

$a_6$     0.057

Coefficients for Eq.(<ref>)

Where  $K_{p_{0}}$ is a Stokes correction factor for the prolate spheroids  given by the following equation:

\begin{equation}\label{eq3plog}
K_{p_0}=\frac{8}{3} p_a^{-\frac{1}{3}}\left(\frac{-2 p_a}{p_a^2+1}+\frac{2 p_a^2-1}{(p_a^2-1)^\frac{2}{3}} \log \frac{p_a+\sqrt{p_a^2-1}}{p_a-\sqrt{p_a^2-1}}\right)^{-1}
\end{equation}

Several power-based drag coefficient correlations exist for prolate spheroids across a wide range of Reynolds numbers and aspect ratios. Ouchene et al.[41]  developed a specific correlation for Reynolds numbers in the range $1\leq Re\leq 240$ and aspect ratios in the range $1.25\leq p_a \leq 32$. This correlation is given by:

\begin{equation}\label{eqpowp1}
C_{D_{0}} = \frac{24}{Re}[K_p +0.15p^{-0.8}_aRe^{0.687}+\frac{(p_a-1)^{0.63}}{24} Re^{0.41}]
\end{equation}

Frohlich et al. [42]  developed correlations for prolate spheroids using Lattice Boltzmann method and a large number of simulations (around 4400 runs) [42]. The resulting correlation is valid for Reynolds numbers in the range of 1 to 100 and aspect ratio values $p_a$ between 1 and 8. The correlation is given by:

\begin{equation}\label{eqpowp2}
C_{D{_0}} = \frac{24 K_p}{Re}f_{d_0}
\end{equation}

Where $f_{d0}$ is given by the following relation:

\begin{equation} \label{eqpowp3}
f_{d_0} = 1+0.15Re^{0.687}+c_{d1}\log(p_a)^{c_{d2}}Re^{c_{d3}+c_{d4}\log(p_a)}
\end{equation}

$c_{d1}$, $c_{d2}$,  $c_{d3}$, and  $c_{d4}$ are constants. Saneevi et al. [29] recently developed correlations for prolate spheroids with aspect ratios $p_a$ in the range of 1 to 16 and Reynolds numbers in the range of 0.1 to 2000. Their correlations are given by:

\begin{equation} \label{eqpowp4}
C_{D{_0}}=\left(\frac{a_1}{R e}+\frac{a_2}{R e^{a_3}}\right) \mathcal{R}+a_5(1-\mathcal{R})
\end{equation}

Where $\mathcal{R}$ is given by:

\begin{equation}\label{eqpowp5}
\mathcal{R}=e^{-a_4 R e}
\end{equation}

The coefficients in Eq.(<ref>) are functions of the aspect ratio. Saneevi et al.'s [29] correlation will be used as the main data source for accurate comparison due to its applicability for a wide range of Reynolds numbers and aspect ratios for prolate spheroids. The above correlations for the prolate spheroids have analogous correlations for the case of $\alpha$ = $90^{\circ}$. The interested reader can refer to the original papers, the same implies to the Stokes correction factor $K$.

Recently Chen et al.[43] developed a correlation for the drag coefficient specifically   aimed at the oblate and prolate spheroids. The correlation was obtained by fusing the numerical data obtained by the authors for prolate spheroids geometries with drag coefficient data from other  investigations for prolate and oblate spheroids.  For example, they used the oblate, and prolate, data  of Sanjeevi et al.[40] that we used to derive  equations for prolate spheroids [34], simultaneously  they used the data of prolate  spheroids  from Ouchene et al.
<cit.>. The correlation is valid for the case of oblate  spheroids  for aspect ratio ($p_a$) and $Re$  in the range of   $0.2\leq p_a <1$, $ 1\leq Re\leq 100$, and for prolate spheroids $1\leq p_a\leq 5$, $1\leq Re\leq 2000$.  The mathematical form  of the correlations that they provide, for the horizontal orientation ($\alpha = 0^\circ$) and for the vertical orientation ($\alpha = 90^\circ$) are the following:

\begin{equation}\label{eqpowg5}
C_{D_0} = \frac{26.47}{Re}p^{-0.120}_a +\frac{4.576}{Re^{0.343}}p^{-0.525}_a+\frac{0.413}{1+16300Re^{-1.09}}
\end{equation}

\begin{equation}\label{eqpowg4}
C_{D_{90}} = \frac{26.58}{Re}p^{0.066}_a +\frac{4.607}{Re^{0.343}}p^{-0.166}_a+\frac{0.413}{1+16300Re^{-1.09}}+p_a^{0.352}-1.0
\end{equation}

The drag coefficient depends on the aspect ratio for high and low Reynolds number regimes regardless of the flow orientation, as shown by both equations Eq.(<ref>) and Eq.(<ref>). It will also be interesting to evaluate the generalization behaviour of Eq.(<ref>) and Eq.(<ref>) since they are derived by using a substantial amount of data.

The correlations that we presented for the case of prolate spheroids show a significant  diversity in their mathematical forms, even though they describe the same geometry. This shows that, in general, we lack a good understanding of the theory that governs the flow over obstacles.

§ THEORY

Our predictions for the drag coefficient will be made using a single logarithmic equation that we derived for the case of a sphere, where the Reynolds number is based on the diameter of the sphere  [37].The drag coefficient is defined as $C_D = \dfrac{F_D}{\frac{1}{2} \rho v_{\infty}^2 A_p}$, where $A_p$ is the frontal area of the sphere, which is defined as $\dfrac{\pi}{4} d_p^2$, and $d_p$ is the diameter of the sphere. The equation is valid for large Reynolds numbers up to $10^6$. The equation takes the following form:

\begin{equation}\label{eqloggeneral}
C_D = a_1 +{\frac{K a_2}{Re}} +{a_3\log(Re)+ a_4\log^2(Re) + a_5\log^4(Re)}
\end{equation}

to 1.0 | X[c]  | X[c]|
Coefficients      Eq.(<ref>)

$a_1$      shape and orientation dependent

$a_2$       24.205

$a_3$      -0.818

$a_4$       0.064

$a_5$        -0.000107

Coefficients for Eq.(<ref>)

The values of the coefficients are tabulated in Table <ref>. Comparing the prefactors of the logarithmic terms in Equation (<ref>) with the corresponding prefactors of the logarithmic terms in Equation (<ref>), we found that they are alike. This similarity suggests that the prefactor values are universal and applicable in cases where the drag coefficient is defined based on the volume equivalent sphere. The value of the $a_1$ coefficient  for the case of the sphere is 3.286. For Eq.(<ref>) to predict the drag coefficient for a geometry other than a sphere, the following assumptions are made:

* The rate of change of the drag coefficient with respect to $Re$ at high Reynolds number will not depend on the object's geometry or its orientation with respect to the flow.

*  The coefficient $K$ and $a_1$ are the only coefficients in Eq.(<ref>) that change with the object's geometry and orientation Throughout this paper, if we mention the $a_1$ coefficient, we mean the $a_1$ coefficient related to Eq.(<ref>), unless we say otherwise.  Several analytical solutions are available to evaluate the correction factor $K$ for different bluff body geometries. In cases where we lack sufficient information about the value of $K$, we assign a value of zero to the correction factor $K$ to cancel out the Stokes term from Eq.(<ref>). The cancellation of the Stokes term will not significantly impact the value of $C_D$ at high Reynolds numbers since the Stokes term does not play an essential role in that flow regime.

* The value of the $a_1$ coefficient by the following two way:

* The value of the coefficient $a_1$ can be deduced from a single measurement of the drag coefficient at high Reynolds numbers, whether obtained through experimental or numerical means. This requirement for only one drag coefficient reading is rooted in one of the outcomes of the mean value theorem [44], which states that if two functions share the same derivative over a specific interval, they differ only by a constant. This constant difference corresponds to the coefficient $a_1$. Once the value of coefficient \(a_1\) is determined, the logarithmic Equation (Eq.(<ref>)) incorporates all definitions of the drag coefficient and Reynolds numbers from the data source used for $a_1$ calculation. This unification of frameworks occurs because of the logarithmic Equation (Eq.(<ref>)) and the data source becomes identical functions.

* In the forthcoming sections of this paper, we will demonstrate that the coefficient $a_1$ is intricately connected to the shear stress at the wall. This connection is established through the prefactor of the $Re^{-1/2}$ term, which governs the skin friction coefficient as derived from the boundary layer theory. By adopting this method to ascertain the value of coefficient $a_1$, Equation (<ref>) becomes independent from any requirement for numerical or experimental data pertaining to the drag coefficient.

We do not have mathematical proof for the above two assumptions, but we will provide extensive comparison with experimental and numerical results to demonstrate their validity.

In our previous work [37], we demonstrated a potential link between the $a_1$ coefficient and the boundary layer separation point. In this work, we aim to provide an additional explanation for the nature of the $a_1$ coefficient by showing that it is strongly correlated with the prefactor of the $Re^{-1/2}$ term in the frictional  drag (shear stress at the wall) experienced by the object and derived using the boundary layer theory. To illustrate this correlation, we will compare the coefficient values $a_1$ with the available prefactor values of $Re^{-1/2}$ obtained from  solutions for the skin drag coefficient in the literature for various object geometries.

For relatively high $Re$, the Stokes drag $\dfrac{Ka_2}{Re}$ becomes insignificant. The drag coefficient for high Reynolds numbers is given by:

\begin{equation}\label{eqhighRe}
C_{D_{h}} = a_1  +{a_3\log(Re)+ a_4\log^2(Re) + a_5\log^4(Re)}
\end{equation}

Where the total drag  at  high Reynolds number is $C_{D_{h}}$. The total drag at high Reynolds numbers  is the summation of the friction drag and pressure drag:

\begin{equation}\label{eqfphighRe}
C_{D_{h}} = C_{D_{hf}}+C_{D_{hp}}
\end{equation}

From the boundary layer theory [45]  we can assume that the skin friction drag in the inertial flow regime  can be approximated by the  following expression:

\begin{equation}\label{frictiondrag}
C_{D_{hf}} = \frac{\delta}{\sqrt{Re}}
\end{equation}

We will demonstrate that there is a strong correlation between the value of $a_1$ and the shear stress at the wall, as represented by the constant $\delta$. For the case of the sphere, for example,  $\delta = a_1 = 3.286$. We know that the sphere is bluff body in which the boundary layer will separate and the above approximation for the frictional drag  may not be valid. However,  there are  plenty of drag correlations  in which  the   $Re^{-1/2}$ term exists, such that of Abraham [46].  For this reason we will keep the current approximation for the frictional drag Eq.(<ref>), and test its validity by comparing our results for the frictional drag  with numerical, and experimental results available in the literature. Churchill [20] integrated further the  results of Frossling [26] for the boundary layer theory for the case of  sphere, and found that the frictional drag is given is by $\dfrac{3.464}{\sqrt{Re}}$. The value of $\delta$ obtained from the logarithmic equation, which is 3.286, is very close to the value of $\delta$ obtained from Churchill's formula [20], which is 3.464. The difference between the two values is less than 5%.

The overall friction  and pressure drag are obtained by adding the contribution from the Stokes drag, as given by the following equations for the case of the sphere:

\begin{equation}\label{eqfrictiondraag}
C_{D_{f}} = \frac{2}{3} \frac{24}{Re}+ \frac{3.286}{\sqrt{Re}}
\end{equation}

\begin{equation}\label{eqformdrag}
C_{D_{p}} = \frac{1}{3} \frac{24}{Re}+(C_{D_{h}}- \frac{3.286}{\sqrt{Re}})
\end{equation}

to 1.25 | X[c] | X[c] | X[c] |X[c] | X[c] |
$Re$   $C_{D_{f}}$ Ref[27]    $C_{D_{f}}/2$ Eq.(<ref>)(%)    $C_{D_{p}}$ Ref[27]   $C_{D_{p}}/2$ Eq.(<ref>)(%)

0.5    17.23    18.46 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 17.23, CD2 = 18.46][2]$%)   8.62   7.68($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 8.62, CD2 = 7.68][2]$%)

1.0    9.13    9.71 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 9.13, CD2 = 9.7][2]$%)   4.58   4.03($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 4.58, CD2 = 4.03][2]$%)

5.0    2.36    2.34 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 2.36, CD2 =2.34][2]$%)   1.23   1.13($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.23, CD2 = 1.13][2]$%)

10.0    1.42    1.32 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.42, CD2 =1.32][2]$%)   0.78   0.75($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.78, CD2 = 0.75][2]$%)

20.0    0.85    0.77 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.85, CD2 =0.77][2]$%)   0.51   0.53($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.51, CD2 = 0.53][2]$%)

40.0    0.53    0.46 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.53, CD2 =0.46][2]$%)   0.36   0.40($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.36, CD2 = 0.40][2]$%)

Comparison of the predictions of the friction drag Eq.(<ref>), and pressure drag Eq.(<ref>) with the results of Dennis and Walker [27].

We will compare the predictions of the friction drag and the pressure drag from Eq.(<ref>), and Eq.(<ref>) with the numerical results of Dennis and Walker [27], the results of the comparison is shown in Table <ref>. Our predictions are close to those of Dennis and Walker[27] for low and moderate $Re$, the largest deviation is about 13%, which is acceptable for an approximation. We will further compare our results with the results of Achenbach [21] for extremely high Reynolds numbers, as shown in Table <ref>. We also included comparison with the results of Abraham [46] equation:

\begin{equation}\label{eqfaridspher}
C_D = 0.29+\frac{24}{Re}+\frac{5.291}{\sqrt{Re}}
\end{equation}

We specifically choose the specific semi-analytical equation because, in its core derivation, it uses the boundary layer theory. We will assume that for the Abraham [46] semi-analytical equation the $C_{D_{f}} = \dfrac{5.291}{\sqrt{Re}}$. The comparison with the experimental results shows the great sensitivity of the results on the value of the constant  ($\delta$). Our results show a maximum difference of about 20% near the drag criss.   At the same Reynolds number, the results of Abraham [46] show a deviation of about 200%, which is substantial as shown in Table <ref>.  From the results of the compression in
Table <ref> we conclude that the value of 5.291  from the Abraham [46] correlation does not correspond to a $\delta$ value related to the frictional drag.The comparison that we made for the frictional drag with different results from the literature is shown in Table <ref>, and Table <ref>, plus that the value of the $a_1$ is close to the prefactor value of the $Re^{-1/2}$ term in Churchill's [20] frictional drag formula for the sphere solidifies the idea that the coefficient $a_1$ may also have strong links to the shear stress around the sphere. We will  show in the coming sections, that we found similar relation for $a_1$ coefficient of the logarithmic equation with frictional drag for the geometries of flat plate and cylinder, and normal flat plate. We emphasise that $\dfrac{1}{\sqrt{Re}}$ terms also appeared in the power-based drag coefficient equations that we derived in our previous work [37]. For example, the most accurate power-based equation that we derived (Eq.(8) Ref[37]) contains $Re^{-1/2}$ term, similar to the frictional drag, and its $\delta$ value is close to the value of the $a_1$  coefficient of the logarithmic equation for the case of  a sphere. The inverse square root Reynolds number term also appears in the power-based equation that we derived directly [37] from the experimental data of Brown and Lawler [35]:

\begin{equation}\label{eqexpspherebrown}
C_D =0.5+\frac{23.22}{Re}+\frac{2.76}{\sqrt{Re}}
\end{equation}

Here the value of the coefficient of the square root term is 2.76, which is still close to the value of 3.28.  The reappearance of the constant $a_1$ coefficient value in different drag correlations of various mathematical structures derived from different data sources shows that the value of the $a_1$ coefficient plays an essential role in the evolution of the drag coefficient for the sphere.

to 1.25 | X[c] | X[c] | X[c]| X[c]|
$Re$    Ref[21]$\dfrac{C_{D_{f}}}{{C_D}}$ (%)    $\dfrac{C_{D_{f}}}{{C_D}}$ (%)    Ref [46] $\dfrac{C_{D_{f}}}{{C_D}}$(%)

$8.9\times 10^{4}$    2.46    2.36 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 2.46, CD2 = 2.36][2]$%)   5.79 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 2.46, CD2 =5.79][2]$%)

$1.1\times 10^{5}$   1.92   2.08 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.92, CD2 = 2.08][2]$%)   5.15 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.92, CD2 =5.15][2]$%)

$1.5\times 10^{5}$    1.52   1.78 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.52, CD2 = 1.78][2]$%)   4.40 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.52, CD2 =4.4][2]$%)

$1.9\times 10^{5}$    1.34   1.61 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.34, CD2 =1.61][2]$%)   3.95 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.34, CD2 =3.95][2]$%)

Comparison of the predictions of the ratio of between the friction drag and the total drag from our predictions, and the predictions of Abraham [46], with the experimental results of Achenbach[21].

Another interesting property of the evolution of the drag coefficient that we believe that plays an essential role evolution of the drag coefficient is the  partial derivative of drag coefficient with respect to the Reynolds number for high $Re$. If we differentiate Eq.(<ref>)  we get the following equation :

\begin{equation}\label{eqhighRe2}
\Delta =\frac{dC_{D_{h}}}{dRe} =\frac{a_3+2a_4\log(Re)+4a_5\log^3(Re)}{Re}
\end{equation}

Where $C_{D{s}}$ is the Stokes drag. We mention that when we will take the derivative of the drag coefficient, we will omit any terms containing the Stokes term ($\sim Re^{-1}$). We will plot the value of $\Delta$   for a number of correlations and models available for a sphere as shown in Figure <ref>. The values of logarithmic-based Eq.(<ref>) are close to the results of Abraham [47], even though the two equations have different mathematical structures.  The trend of the two previous equations shows that the derivative values obtained from Eq.(<ref>) do not belong to a specific source of the training data. On the contrary, it shows that derivative values obtained from Eq.(<ref>) represent the actual general behaviour of the change of the drag coefficient with the Reynolds number for a sphere. The derivative values obtained from  Eq.(<ref>), and Holzer and Sommerfeld [38] correlation, differ significantly from those obtained by the logarithmic equation or that of Abraham [46]. Despite all of these equations having similar accuracy in predicting the drag coefficient at the specific $Re$ regime that are tested, the deviation in the $\Delta$ values indicates that the coefficients in the drag equations must be carefully adjusted to produce the $\Delta$ values obtained from Eq.(<ref>). For example,  the correlations  by Abraham [46], Holzer and Sommerfeld [38], and Eq.(<ref>) have the same mathematical structure. Still, the difference in the coefficient values produces a significant deviation in $\Delta$ values of the different correlations, and models. In the coming sections, we will demonstrate that the values of $\Delta$ obtained from Eq.(<ref>) are universal and apply  not only spheres but also to other non-spherical geometries. The values of $\Delta$ are indicators of the accuracy of drag coefficient equations, particularly in the high Reynolds number regime.

Comparison  between the value of the derivative of the drag coefficient (exclusion the Stokes term $\sim Re^{-1}$) with respect to the $Re$ for different correlations,  for the case of a sphere.

§ RESULTS

We will first confirm our assumption about the universality of the drag coefficient's rate of change with the Reynolds number. Then, we will demonstrate how our theory can be used to determine the drag coefficient for different shapes, including flat plates and irregularly shaped particles.

§.§ The variation of $\Delta$ for different geometries

We will demonstrate in the current subsection that at moderate and high Reynolds numbers, regardless of the shape of the object and their orientation to the flow direction, the value of $\Delta$ is the same. The value of  $\Delta$ indicates how the drag force changes with the changes in the flow field around the object.   To the best of  our knowledge, no one investigated a similar parameter such $\Delta$ before. We used an automatic differentiation Python-based library to obtain the derivative of different correlations and predictive models. In the process of estimating the derivatives of the different correlations, we omitted the Stokes term $Re^{-1}$ when we calculated the derivatives of the different correlations .  Only for the case of the experimental results of Wieselsbergerg [19] for the case of the cylinder, we could not omit the Stokes drag from the data, because the analytical solution of Lamb's for the case of Stokes flow around the cylinder diverges for $Re > 10$. However, this will not affect the accuracy of the results since the Reynolds numbers considered are at the high end.

To assess our results of $\Delta$, we will examine the variation of $\Delta$ from Eq.(<ref>) with other correlations and experimental data found in the literature
Figure <ref>. Our selection includes the oblate spheroid geometry from Sanjeevi et al.'s numerical results [40], as well as two prolate spheroid geometries from Sanjeevi et al. [29]. These geometries were chosen due to their capacity to provide information on the $C_D$ variation with Reynolds numbers up to 3000, which would enable us to test the $\Delta$ value variation for a broad Reynolds number range. Furthermore, experimental-based results will be utilized to compare our findings. We have opted to use Wieselsberger's results [19] for the scenario of an infinite cylinder with right angles to the flow.

Comparison of the $\Delta$ value obtained from Eq.(<ref>) with $\Delta$ values  from different sources available in the literature.

From Figure <ref>, it appears that the values of $\Delta$ are converging, regardless of the source of the data, the type, size, and orientation of the bluff body. This observation is consistent with our theory, and leads us to conclude that a single value for $C_D$ at high $Re$ can be used to determine the value of $a_1$ in Eq.(<ref>) for flows beyond that of a sphere. The deviations we observe in the $\Delta$ values from Wieselsberger's [19] experimental results can be attributed to the fact that the digitally extracted data is not uniformly distributed, resulting in inaccurate gradient calculations. However, overall, the values of $\Delta$ for Wieselsberger's [19] experiments are close to the $\Delta$ values calculated by Eq.(<ref>), even for high $Re$ where $\Delta$ values are on the order of $10^{-4}$.

§.§ Drag coefficient

In this subsection we will deploy Eq.(<ref>) to different bluff body geometries to test our theory.

§.§.§ Horizontal flat plate

Flat plate geometry plays a vital role in aerodynamics as it serves as a benchmark for testing fluid mechanics theories. The boundary layer theory solution by Blasius and Prandtl[48, 15] applies to a flat plate with flow velocity parallel to its main side.  It is worth noting that for flat plate geometry, only the skin friction drag contributes to the total drag force. Consequently, we can evaluate the validity of our theory at moderate Reynolds numbers where skin friction drag prevails over pressure drag (form drag). To further bolster our findings, we will compare our results to those obtained by Melnik and Chow [49] using triple deck boundary layer theory. Specifically, their analytical solution for a single wetted side finite flat plate is the following:

\begin{equation}\label{eqflat1}
C_{D_{f}} = \frac{1.328}{\sqrt{Re}}+\frac{2.67}{Re^{7/8}}
\end{equation}

While, for the double wetted side finite flat, the drag coefficient is given by the following equation:

\begin{equation}\label{eqflat2}
C_{D_{f}} = \frac{2.656}{\sqrt{Re}}+\frac{5.34}{Re^{7/8}}
\end{equation}

In Melnik and Chow [49] model, the Reynolds number is defined based on the plate length, and the drag coefficient is based on the plate's surface area.  We will compare the $\Delta$ values from the above equations with the logarithmic based equation Eq.(<ref>) to gauge the validity of our assumptions, as illustrated in Figure <ref>. A close examination of Figure <ref> reveals two intriguing insights. Firstly, the $\Delta$ values for the single-wetted flat plate deviate significantly from those of the logarithmic-based Eq.(<ref>). In contrast, the $\Delta$ values for the double-wetted flat plate are relatively closer to the values of Eq.(<ref>). Secondly, the $\Delta$ values for the double wetted flat plate tend to converge with the results of Eq.(<ref>) at moderate Reynolds numbers ($Re$ = 20 to 100) as we had hypothesized. In the range where the $\Delta$ values of the logarithmic and double wetted flat plate equations converge, Eq.(<ref>) and Eq.(<ref>) will be identical functions. The logarithmic-based equation (Eq.(<ref>)) will inherit the definitions of $Re$ and $C_D$ from Eq.(<ref>) as our theory describe.

In order to use the logarithmic-based equation (Eq.(<ref>)) for a flat plate, it is necessary to determine the value of the coefficient $a_1$. Since the logarithmic equation (Eq.(<ref>)) and Eq.(<ref>) have the same derivative within a narrow Reynolds number range of 20 to 100, we can evolute the value of the $a_1$ by plunging the value of $C_D$ from Eq.(<ref>) at $Re$ = 50 to Eq.(<ref>). The value of the coefficient $K$ will be assigned a value of zero since we are not working in the creeping flow regime the coefficient $a_1$ is found to be 2.795.

The $a_1$ value obtained from the logarithmic equation  is  about $\eval*{-100\frac{2.656-2.795}{2.656}}[2]$% different from the actual value of  Blasius'[17] 2.656, indicating a close relationship between the $a_1$ coefficient  and the shear stress at the wall obtained from the  boundary layer theory. The value of $a_1$ emerged without using any elements of the boundary layer theory, suggesting that it could be part of a more general theory of fluid mechanics yet to be discovered. For the sake of comparison, the widely used commercial CFD software Ansys Fluent [50] requires more than half a million grid points to converge to the Blasius'[17] solution, indicating that the computational process is  expensive.
The logarithmic equation Eq.(<ref>) is compared to the analytical solution from the boundary layer theory and the experimental data of Janour[51] as shown in Figure  <ref>. The comparison shows that the logarithmic equation closely follows the analytical solution in the Reynolds number range of 20 to 100. This is expected as the logarithmic equation was derived from the bluff body data of a sphere, where at high Reynolds numbers, form drag dominates over skin friction drag.

Suppose we use the Brown and Lawler [35] correlation the source of training data that we used for the derivation of the basic version of the logarithmic equation(Eq(.<ref>)) and try to find the constant that makes the difference between its value of the drag coefficient and that from Eq.(<ref>) equal to zero (similar to the $a_1$ coefficient) at $Re$ = 50. We will find that the value of this constant is -0.486, which has no physical significance. Only logarithmic-based equations produce values close to the ones of the boundary layer theory. This makes the logarithmic-based equations a strong candidate to represent the solution of the Navier-Stokes equations more closely than the available power-based correlations, as we showed in our previous publication [37].

Comparison of $\Delta$ values from  Melnik and  Chow [49] analytical solution for the case of the flat plate for the single, and double wetted cases  with the logarithmic based equation Eq.(<ref>) for different Reynolds numbers.

Comparison of $C_D$ values predicted by the logarithmic equation, and the boundary layer theory, and experiments for the case of horizontal flat plate.

§.§.§ Normal flat plate

The flow over a flat plate normal to its direction is the simplest case of a bluff body. At the edges of the plate, the flow separates. In the middle of the nineteenth century, Kirchhoff and Rayleigh [52] calculated the drag force for the normal flat plate. They assumed that the shear layers leave the flat plate at the edges and that the dynamic pressure equals the free stream pressure, and that the total drag force is due to the dynamic pressure. They derived the following formula for the drag coefficient:

\begin{equation}\label{eqnormal_flat_plate}
C_D = \frac{4\pi \sin{\alpha}}{4+\pi \sin{\alpha} }
\end{equation}

The drag coefficient is based on the frontal area of the plate. The angle of attack ($\alpha$) is 90 degrees for a normal flat plate case. Eq. (<ref>) gives a drag coefficient of 0.88, which differs by
*$(2.0-0.88)/(2.0)\times 100 $% from the  measured value of 2.0 from experiments. For the case of the normal flat plate, the drag coefficient data are sparse, and even finding correlations is difficult.  We will not use any experimental data or data from correlations to calculate the value of the $a_1$ coefficient,  we will rely on the boundary layer theory instated .

The normal flat plate case is a special case for the family of the Falkner-Skan boundary layer flows, which their numerical solution is known.  The skin friction coefficient for the case of the normal flat plate is given by [53], and its formula is
the following:

\begin{equation}\label{eqskin_friction_normal_plate}
C_{D_f} = \frac{4.932}{\sqrt{Re}}
\end{equation}

following our theory will assign the value of 4.932 for the $a_1$ coefficient of Eq.(<ref>) for the case of the normal flat plate.

to 1.25 | X[c] | X[c] | X[c]| X[c]|
$Re$   $C_D$   $C_D$ Eq.(<ref>) $a_1$ = 4.932, $K$ =0

Ref [54](DNS)   250    2.36    2.26 ($\eval*{$\frac{CD1-CD2}{CD1}100$}[CD1 = 2.36, CD2 = 2.266][2]$%)

Ref[55](DNS)   $1000$   2.31     2.09 ($\eval*{$\frac{CD1-CD2}{CD1}100$}[CD1 = 2.31, CD2 =2.09][2]$%)

Ref [52](Exp)    $5000$   2.0    2.04 ($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 2.0, CD2 =2.04][2]$%)

Ref[56](Exp)   $1.5\times 10^{5}$    2.13    2.11($\eval*{$\frac{CD1-CD2}{CD1}100$}[CD1 = 2.13, CD2 =2.11][2]$%)

Ref[57](LES)   $1.5\times 10^{5}$    2.2    2.11($\eval*{$\frac{CD1-CD2}{CD1}100$}[CD1 = 2.2, CD2 =2.11][2]$%)

Relative error between the predictive $C_D$ coefficient from the
logarithmic equations, and those form different sources from the literature
for different Reynolds numbers for the geometry of normal flat plate.

We will compare our logarithmic equation with the existing results in the literature, as shown in Table <ref>. Since the drag coefficient for the normal flat plate case is predominantly influenced by form drag (pressure drag), we expect the logarithmic equation to hold at high Reynolds numbers. The results presented in Table <ref> confirm this expectation, with the relative error of the logarithmic equation being below 5% in the majority of cases tested across a wide range of Reynolds numbers, ranging from $250$ to $10^5$. Moreover, the logarithmic equation
Eq.(<ref>) with $a_1 = 4.932$ and $K = 0$ outperforms the only analytical model available for the normal flat plate case
Eq.(<ref>). One intriguing aspect of the normal flat plate case is the strong relationship between the $a_1$ coefficient value and the shear stress value derived from boundary layer theory for the same geometry.

§.§.§ Infinite cylinder

The study of fluid flow over infinite cylinders has been a subject of immense  investigation, primarily due to the proliferation of complex, unsteady flow patterns associated with this particular bluff body geometry. Simulations using large eddy turbulence models have shown strong dependence of the flow parameters on factors such as grid resolution, time steps, and aspect ratio. The system's dimensionality also plays a crucial role, with two-dimensional simulations only reliable up to a Reynolds number of 250 [58]
. Beyond that, three-dimensional flow simulations are necessary to adequately capture complex flow structures and obtain accurate results for parameters such as the drag coefficient. Machine learning methods, primarily using different types of neural networks, have been applied to the flow over cylinders and trained on extensive flow field data. However, these studies have yet to accurately predict flow parameters, such as the drag coefficient at high Reynolds numbers [59, 60]
. For instance, despite including $Re$ = 1000 flow data in their training set, Vlachas et al. [61] only achieved a 15% relative error in their prediction of the drag coefficient at that specific Reynolds number.

In our current investigation, we will examine the geometry of an infinite cylinder. The configuration is similar to the experiments conducted by Wieselsberger [19], where the flow velocity flows perpendicular to the cylinder's axis. In the experimental data of Wieselsberger [19], the $Re$ is defined on the bases of the diameter of the cylinder, and the $C_D$ is defined based on the frontal area. With this specific geometry in mind, Thom [25] utilized boundary layer theory to calculate the skin friction coefficient for the front part of the cylinder up to $60^{\circ}$ degrees from the forward stagnation point. Through his calculations, he derived the following relation:

\begin{equation}\label{eqcylinederThom1}
C_{D_{f}} = \frac{3.84}{\sqrt{Re}}
\end{equation}

For the rear of the cylinder
Thom [25] calculated the skin friction coefficient by slightly modifying the coefficient presented in Eq.(<ref>) based on experimental results. As a result, he proposed the following equation as a comprehensive representation of the skin friction for the entire surface of the cylinder:

\begin{equation}\label{eqcylinederThom2}
C_{D_{f}} = \frac{4.0}{\sqrt{Re}}
\end{equation}

Takami and Keller [62] derived the necessary constants for the Imai [63] semi-theoretical equation, which is based on the Kirchhoff type of flow, from their numerical findings. The Imai [63] semi-theoretical equation takes the following form:

\begin{equation}
C_D  = \sqrt{0.5+\frac{3.547}{\sqrt{Re}}}
\end{equation}

Munson et al. [64] utilized data from experiments to construct a semi-empirical equation specifically designed for high Reynolds numbers. The resulting equation is as follows:

\begin{equation}
C_D = 1.17+ \frac{5.93}{\sqrt{Re}}
\end{equation}

To apply the logarithmic-based equation Eq.(<ref>) to the cylinder case, it is necessary to determine the value of the coefficient $a_1$. We will assign a value of zero to the coefficient of $K$ for the Stokes term, since we our only interested for high $Re$ regime. We will calculate $a_1$  utilizing  two methods.
The first approach is to calculate the $a_1$ coefficient based on the value of the $C_D$ at $Re$ = 114.58 from the Wieselsberger experiments [19], $C_D$ is equal to 1.37. Using this information, we obtain $a_1$ = 3.865.
The second method based on our theory which links the value of the $a_1$ with the frictional
drag coefficient. For this reason we will assign a value of 4.0 to $a_1$ coefficient following the estimations of Thom [25] for frictional drag coefficient for the cylinder Eq.(<ref>).

The value  of  the $a_1$ coefficient calculated using the  Wieselsberger's[19] data
differs only by  $\eval*{-100\frac{3.84-3.86}{3.84}}[2]$%  from the value calculated by Thom [25] for the front part of the cylinder using boundary layer theory, affirming our hypothesis that the $a_1$ coefficient is linked to the cylinder's skin friction. The evolution of the drag coefficient with Reynolds number is shown in Figure <ref>. The predictive models of Takami and Keller [62] and Munson et al. [64] were unable to produce accurate predictions, but the logarithmic equations were found to be more accurate predicting the  Wieselsberger's[19]
experimental results. The logarithmic equations were the first mathematical models of the drag coefficient that were able to accurately predict the results of Wieselsberger's [19]
experiments, which are over a century old. Eq.(<ref>) with $a_1$ = 3.865, and $K$ = 0 is capable of accurately predicting the experimental results  for Reynolds numbers up to 3000, while Eq.(<ref>) with $a_1$ = 4.0, and $K$ = 0 is effective for Reynolds numbers higher than
$10^{4}$. An intriguing observation is that Eq.(<ref>) with $a_1$ = 3.865, and $K$ = 0 can predict the minimum point in the drag coefficient curve, which occurs around Reynolds number of 3000, and how close its predictions are to the LES simulations conducted by Jiang and Cheng[65]. The logarithmic equation Eq.(<ref>) with $a_1$ = 3.865, and $K$ = 0 accurately predicts the numerical results of Suker and Bruuer [66] for low Reynolds numbers, and Eq.(<ref>) with $a_1$ = 4.0, and $K$ = 0 predicts the LES results of Cheng et al. [28]. The overall performance of the logarithmic equations supports our belief that the rate of change of the drag coefficient with the  Reynolds number is given by Eq.(<ref>), which underscores its universality. We aim to develop our theory further so that we can arrive at a single equation that can describe the entirety of the cylinder data. Table  <ref> reveals that the accuracy of the logarithmic equation Eq.(<ref>) with $a_1$ = 3.865, and $K$ = 0 is satisfactory in the lower and medium Reynolds number regimes. Even at a high Reynolds number of $10^5$, the deviation from accuracy is only 14% which is still acceptable. Conversely, the logarithmic equation Eq.(<ref>) with $a_1$ = 4.0, and $K$ = 0 performs exceptionally well in the high Reynolds number range.

Comparison of $C_D$ values predicted by the logarithmic equations, other predictive models, experiments, and numerical simulations for the case of infinite cylinder.

to 1.0 | X[c] | X[c] | X[c] |X[c] |
$Re$   Wieselsberger[19]     Eq.(<ref>)$a_1$ = 3.865, $K$ = 0(%)    Eq.(<ref>)$a_1$ = 4.0, $K$ = 0(%)

1.59     1.63($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.59, CD2 = 1.63][2]$%)    1.77($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.59, CD2 = 1.77][2]$%)

$\eval*{$\frac{573.28 }{10^2}$}[2]$ $\times 10^2$

1.07($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.16, CD2 = 1.07][2]$%)    1.21($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 =  1.07, CD2 =1.21][2]$%)

$\eval*{$\frac{2302.77 }{10^3}$}[2]$ $\times 10^3$

0.98($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.92, CD2 = 0.98][2]$%)    1.11($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.92, CD2 =1.11][2]$%)

$\eval*{$\frac{48888 }{10^4}$}[2]$ $\times 10^4$
1.04($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.16, CD2 = 1.04][2]$%)    1.17($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.16, CD2 =1.17][2]$%)

$\eval*{$\frac{179207 }{10^5}$}[2]$ $\times 10^5$

1.04($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.21, CD2 = 1.04][2]$%)    1.17($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.21, CD2 =1.17][2]$%)

Relative error between the predictive $C_D$ coefficient from the logarithmic equations , and the experiments of Wieselsberger[19] for different Reynolds numbers  for the geometry of infinite  cylinder.

Figure <ref> shows the variation of the skin friction drag coefficient for infinite cylinder. The skin friction for logarithmic equations are the following $\dfrac{3.865}{\sqrt{Re}}$, and $\dfrac{4.0}{\sqrt{Re}}$ receptively.
It can be seen from Figure <ref> that both logarithmic equations provide close predictions of the skin friction drag coefficient compared to numerical results obtained by Sucker and Bruuer [66] at low Reynolds numbers. The predictions also agree with the experimental results by Goldstein [67] for both low and high Reynolds numbers.

In comparing the form drag coefficient $C_{D_{p}}$, the results from logarithmic equations as depicted in Figure <ref> are analyzed against the Large-Eddy Simulations (LES) carried out by Jiang and Cheng[65]. Interestingly both equations capture  the increasing trend of the form drag coefficient with the Reynolds number as described by the results of Jiang and Cheng[65]. It is noted that the logarithmic equation (Eq.(<ref>) with $a_1$ = 3.865, and $K$ = 0) exhibits an excellent agreement with the simulation results of Jiang and Cheng[65] across a broad range of Reynolds numbers, with an average relative difference of merely 1.5% as depicted in Table <ref>. However, a higher deviation is observed in the logarithmic equation (Eq.(<ref>) with $a_1$ = 4.0, and $K$ = 0), with an average relative difference of approximately 16%.

to 1.0 | X[c] | X[c] | X[c] |X[c] |
$Re$   Jiang and Cheng[65]      Eq.(<ref>)$a_1$ = 3.865, $K$ = 0(%)    Eq.(<ref>)$a_1$ =4.0, $K$ = 0(%)

$\eval*{$\frac{1.516557829568893112 *10^3 }{10^3}$}[2]$ $\times 10^3$

0.873     0.899($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.873, CD2 = 0.899][2]$%)    1.031($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.873, CD2 = 1.031][2]$%)

$\eval*{$\frac{2.020521197191781994 *10^3 }{10^3}$}[2]$ $\times 10^3$

0.885     0.901($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.885, CD2 = 0.901][2]$%)    1.033($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 =  0.885, CD2 =1.033][2]$%)

$\eval*{$\frac{2.505542495503466398 *10^3 }{10^3}$}[2]$ $\times 10^3$

0.904($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.892, CD2 = 0.904][2]$%)    1.036($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.892, CD2 =1.036][2]$%)

$\eval*{$\frac{3.994673144433314064 *10^3 }{10^3}$}[2]$ $\times 10^3$
0.915($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 =  0.922, CD2 = 0.915][2]$%)    1.031($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 =  0.922, CD2 =1.031][2]$%)

Relative error between the predictive $C_D$ coefficient from the logarithmic equations , and the LES simulations of
Jiang and Cheng[65] for different Reynolds numbers for the geometry of infinite  cylinder.

Comparison of $C_{D_{f}}$ (skin friction drag coefficient)  values predicted by the logarithmic equations,  experiments, and numerical simulations for the case of infinite cylinder.

Comparison of $C_{D_{p}}$ (form drag coefficient)  values predicted by the logarithmic equations, and numerical simulations for the case of infinite cylinder.

§.§.§ Cubes

A Cube is another bluff geometry,  for which we have experimental data for the drag coefficient at high Reynolds numbers, so we can test further the validity of our logarithmic equation. The only correlation that we know  of for the drag coefficient of a cube is that of Saha [68], which is based on numerical simulations. It has the following form:

\begin{equation}\label{eqcubepower}
C_D = \frac{24}{Re}(1+0.232Re^{0.628})
\end{equation}

We will determine the coefficient $a_1$ in Equation of the logarithmic equation (Eq.(<ref>))  by utilizing the drag coefficient obtained from the Saha correlation [68] at a Reynolds number ($Re$) of 100, the value of $K$ = 1.0 following from the correlation of Saha[68]. The resulting value of $a_1$ is found to be 3.488.

We will compare our results to the experimental findings of Khan et al. [22] for the case of a cube placed normal to the flow. In addition, we will also include in our comparison the correlations of Saha [68] and Haider and Levenspiel [39], as illustrated in Figure <ref>. Khan et al. [22] estimated the drag coefficient from their experimental readings by employing both the Wake survey and the modified Wake survey methods.

For low Reynolds numbers, the drag coefficient values predicted by the logarithmic equation closely matches those from the correlations of Haider and Levenspiel [39] and Saha [68]. However, at high Reynolds numbers, the logarithmic equation
(Eq.(<ref>), $a_1$ = 3.488, $K$ = 1.0) predicts the experimental results of Khan et al. [22] with acceptable accuracy, capturing the increasing and decreasing trends in their data. The logarithmic equation predicts the results of the modified Wake survey method more closely.

On the other hand, the correlations of Saha [68] and Haider and Levenspiel [39] under-predict and over-predict the experimental results of Khan et al. [22], respectively.

Comparison of $C_D$ values predicted by the logarithmic equation, correlations, and  experiments, for the case  of a cube.

The drag coefficient for the incoming geometries relies on the concept of a volume-equivalent sphere, while the Reynolds number is determined based on the diameter of this volume-equivalent sphere.

§.§.§ Oblate spheroids

For oblate spheroids, we can compare our theory with a significant number of spheroid geometries, as well as different flow orientations.
For the sake of brevity we will focus mainly for case of flow that is parallel to the equatorial diameter of the spheroid ($\alpha = 90^{\circ}$), and for different aspect ratios $p_a$. To determine the $a_1$ coefficient in the logarithmic equation, we will utilize the $C_D$ value from  the Ouchene [30] correlation and the numerical data provided by Sanjeevi et al. [40]. We will be using the values of $C_D$ from the Ouchene correlation [30] for aspect ratios of $p_a = 0.8$, $0.5$, and $0.2$ at $Re = 100$ and $\alpha = 90^{\circ}$. The corresponding values of $a_1$ that we calculated are $3.361$, $3.365$, and $3.365$, respectively. It is interesting to note that the value of the $a_1$ coefficient remains almost constant for all aspect ratio values used. Therefore, we have decided to keep the value of the $a_1$ coefficient at $3.365$.

We will start by comparing the predictions of our  logarithmic equation  with those of Ouchene [30], and with other predictive models available from the  literature as   shown in Figure <ref> for values of $p_a$ $\geq$ 0.1. The drag coefficient predicted by the Ouchene [30] correlation follow a pattern  that is independent  of the aspect ratio at high Reynolds numbers, similar to the observation by El Hasadi and Padding [34], and Frohlich et al. [42] for the prolate spheroids. Beyond the Ouchene [30] correlation only the logarithmic predictive model of  Eq.(<ref>),$a_1 = 3.365$, and $K = K_o$ predicted the aspect ratio independency of the drag  coefficient independently, here what we mean by $K = K_o$ is that the $K$ value in the logarithmic equation is calculated form Eq.(<ref>) for each value of the aspect ratio $p_a$.

Comparison  between the drag coefficient $C_D$ predicted by different predicted models for the case of oblate shaped particle, for different Reynolds numbers, and high aspect ratio values, and $\alpha$ = $90^{\circ}$.

The power based correlations such as Chen et al.[43], Haider and Levenspiel [32], and Holzer and Sommerfeld [38] show different degree of accuracy with respect to the Ouchene [30] results. For the case of the Chen et al. [43] correlation  its divergence from the predictions of Ouchene [30] is only due to its power based  mathematical structure that they used as we showed in our previous work [37] in which we showed that power based drag coefficient correlations are poorly extrapolate beyond their training data.
For the case of the Haider and Levenspiel [32] correlation its basic concept that the geometry of the particle (i. e. $\phi$) control the evolution of the drag coefficient in low and high Reynolds number regimes contradicts  the numerical predictions of Ouchene et al.[30]. The Holzer and Sommerfeld [38] is the only correlation beyond the predictive logarithmic equations that approximately predicted the aspect ratio independence of the drag coefficient at high $Re$. However, if we zoom further  we will conclude that Holzer and Sommerfeld [38] correlation  predicts a weakly non-linear relation between the drag coefficient with the aspect ratio of $p_a \geq$ 0.1, and the Reynolds number.

to 0.8 | X[c] | X[c] |
Average Relative Error (%)

Eq.(<ref>) $a_1$ = 3.365,$K= K_o$     4.36 %

Chen et al.[43]     19.09%

Haider and Levenspiel [32]     43.38%

Holzer and Sommerfeld [38]     18.74%

Relative average error of the predictions of  different predictive models with the respect to the results of Ouchene [30].

The average relative error for different predictive models and correlations with respect to the data of
Ouchene [30]   for the different $Re$ considered in Figure <ref>  is given in Table <ref>. The worst performing correlation is that of Haider and Levenspiel [32] with average error of about 43%, which is not suitable for making accurate  predictions for oblate geometries in general especially for the ones with low  value of the aspect ratio. The logarithmic equation Eq.(<ref>), out performs all the other predictive models by average error of 4.36%.   One issue in fluid flow investigations, particularly for multiphase problems, is the lack of benchmark cases for validating numerical simulations or experimental results. This makes it difficult to assess the accuracy of Ouchne's [30] results. Therefore, further comparisons with available results in the literature must be made. Thus, will compare our predictive models with further results from the literature either from numerical simulations or from experiments, for different Reynolds numbers, and different  $p_a$ values. Finding relevant experimental results to compare with our predictive models is difficult due to the lack of empirical investigations. For this reason, we selected only one experiential investigation to compare our results with.

The only available experimental results for oblate spheroids at high Reynolds numbers are found in the empirical correlation provided by Clift [33] for the case of an oblate with $p_a = 0.5$. Clift's experimental correlation (Eq.(<ref>)) provides us with the opportunity to compare our results for Reynolds number ranges beyond those covered by Ouchene's [30] numerical results. The comparison is shown in Figure <ref>. The logarithmic equation (Eq.(<ref>), $a_1 = 3.365$, and $K = 0.998$) predictions are very close to the empirical results of Clift [33] for the entire range of Reynolds numbers considered. These results solidify our assumption that $a_1$ is constant for aspect ratios in the range $0.2\leq p_a\leq 0.8$.
The results of Ouchene [30] are close to those of Clift [33] for Reynolds numbers up to 500. However, for higher Reynolds numbers, they start to deviate. On the other hand, the predictions of Chen et al. [43] and Holzer and Sommerfeld [38] significantly under-predict and over-predict, respectively, the empirical results of Clift [33].
The logarithmic equation predictions are in excellent agreement with the empirical results of Clift [33], with an average difference of 5.69%, as shown in Table <ref>. In contrast, the correlation of Chen et al.[43] significantly under-predicts the results of Clift [33] by 27.72%.

Comparison  between the drag coefficient $C_D$ predicted by different predictive  models for the case of oblate shaped particle with the experimental  results of Clift [33] for $p_a$ = 0.5, and for  $\alpha$ = $90^{\circ}$.

to 1.0 | X[c] | X[c] |
Average Relative Error (%)

Eq.(<ref>),$a_1 = 3,365, K = 0.998$     5.69%

Chen et al.[43]     27.72%

Ouchene [30]     9.34%

Holzer and Sommerfeld [38]     17.91%

Relative average error of the predictions of  different predictive models with the respect to the results of  Clift [33] experimental correlation for the case of oblate of $p_a$ =0.5, and $\alpha = 90^{\circ}$.

The only numerical results that we have for high Reynolds number flow over an oblate geometry is that of Sanjeevi et al.[40] for the case of $p_a$ = 0.4, for different flow orientations  that includes flow parallel to the equalutral axis ($\alpha = 90^{\circ}$), and polar axis ($\alpha = 0^{\circ}$). In addition, we have calculated the $a_1$ coefficient using the numerical data of Sanjeevi et al. [40] for $Re$ = 100, $\alpha$ = $90^{\circ}$, and $p_a$ = $0.4$. The value we obtained for $a_1$ is 3.226. The difference between the values of the $a_1$ coefficient obtained from the data of Ouchene [30] and that of Sanjeevi et al. [40] is due to the different numerical schemes that they used, which produce results with varying levels of accuracy. The Reynolds number for both Ouchene [30], and  Sanjeevi et al. [40] data is based on the volume equivalent sphere. We need now to calculate the $a_1$ for the case of flow orientation of $\alpha = 0^{\circ}$, and we will do it for both the  Ouchene [30] and that of Sanjeevi et al. [40] data. For the case of Ouchene [30], we will evaluate the $a_1$ coefficient using the data for the $C_D$ for case of $Re$ = 100, $p_a$ =0.4, and for this case we obtain a value of $a_1$ = 4.242. While for the   Sanjeevi et al. [40] data we will use the data of the $C_D$ at $Re$ = 300, for better accuracy, for this case we obtained a value for the $a_1$ equal to 4.186. We want to emphasize that Sanjeevi et al. [40] data are discrete points directly obtained from numerical simulations, and not influenced by any correlations.

This time we will test our logarithmic equations for two different flow orientations, as shown in Figure <ref>. For the case of $\alpha = 90^{\circ}$, the logarithmic equation with $a_1 = 3.365$ and $K = 1.016$, which is derived from the Ouchene [30] data (Eq.(<ref>)), follows the Sanjeevi et al. [40] numerical data from low Reynolds numbers up to a value of 100. However, for higher Reynolds numbers, there is a slight deviation between the drag coefficient predicted by the logarithmic equation and the numerical data of Sanjeevi et al. [40]. On the other hand, the logarithmic equation with $a_1 = 3.226$ and $K = 1.016$, which is derived from the Sanjeevi et al. [40] data (Eq.(<ref>)), follows the numerical results of Sanjeevi et al. [40] from low to high Reynolds numbers.
For the case of $\alpha = 0^{\circ}$, where the oblate particle experiences higher  drag force due to its bluff shape compared to $\alpha = 90^{\circ}$, the logarithmic equations with $a_1 = 4.242$ and $K = 1.206$, and $a_1 = 4.186$ and $K = 1.206$, which are derived from the Ouchene [30] and Sanjeevi et al. [40] data (Eq.(<ref>)), respectively, follow the numerical data of Sanjeevi et al. [40] from low to high Reynolds numbers. However, the correlation of Ouchene [30] significantly underpredicts the values of the numerical data of Sanjeevi et al. [40], which shows the limitation of its applicability. We have shown that our theory holds for the oblate spheroid with different flow geometries.

Comparison  between the drag coefficient $C_D$ predicted by different predictive  models for the case of oblate shaped particle with the numerical results of Sanjeevi et al.[40] for $p_a$ = 0.4, $\alpha = 0^{\circ}$, and $90^{\circ}$.

§.§.§ Prolate spheroids

The flow behaviour of prolate spheroids differs from that of oblate spheroids, and it is more sensitive to the geometry aspect ratio. For example, Sanjeevi et al. [29] found that at low Reynolds numbers ($Re$), the drag coefficient reaches a minimum at $p_a \simeq 1.95$ and then increases as the aspect ratio increases. Therefore, we will separately test our equations in the low [1-2] and high [2-15] aspect ratio ranges. We will investigate extensively, the case when the flow is parallel to the major axis of the spheroid  $\alpha$ = $0^{\circ}$.

To calculate the coefficient $a_1$, we will use the data from the correlation provided by Sanjeevi et al. [29] at a Reynolds number of $Re$ = 100, considering different aspect ratios of the prolate spheroid geometry.
For the low range of aspect ratios [1-2], we will use the data for the case of $p_a$ = 1.2. From the data, we found that $a_1$ = 3.1416, which is approximately equal to $\pi$.
For the higher range of aspect ratios [2-15], we will calculate the $a_1$ coefficient for three specific aspect ratio values: $p_a$ = 5, 10, and 15. The calculated values of $a_1$ are approximately 2.989, 2.994, and 2.997, respectively, which shows that they are equal. To simplify, we will approximate the value of $a_1$ equal to 3.00.

For aspect ratios between 1.1 and 1.7, we plot the variation of the drag coefficient with $p_a$ for different $Re$ as shown in Figure <ref>. We have minimal data points for the drag coefficient in this aspect ratio range, even from CFD investigations. Therefore, the correlations we will use to compare our results with, usually only contain a single point in this aspect ratio range in their fitting data. For the case of $Re$  = 10, all the logarithmic-based equations produce results that are close to the results of the majority of the correlations that are based on CFD simulations, while the correlations of Chen et al.[43], and Holzer and  Sommerfeld [38], they over, and under predicts respectively the results of the CFD based correlations. Increasing the Reynolds further to values  50 and 500 in the medium Reynolds number regime, we observe that the values of logarithmic equation Eq.(<ref>) start to deviate and under predicts the results of the correlations of the CFD simulations.  We mean by CFD-based correlations the correlations of Sanjeevi et al. [29], Ouchene et al.[41], and  Frohlich et al. [42]. Interestingly,the logarithmic equation (Eq.(<ref>),$a_1$ = $\pi$, and $K=K_p$) predictions are close to the those of Sanjeevi et al. [29], and  Frohlich et al. [42]. If we increase further the Reynold number to 1000, the logarithmic-based equations keep similar predictive behaviour as in the medium Reynolds number flow regime, with the predictions of the logarithmic equation (Eq.(<ref>), $a_1$ = $\pi$, and $K=K_p$) getting closer to the predictions of Sanjeevi et al. [29], and Ouchene et al.[41], it is expected that the deviation will be increased from the predictions of  Frohlich et al. [42] since their correlations is valid only up to Reynolds number 100.  Still, there needs to be more information on how to drag coefficient varies in this specific aspect ratio range to have a better conclusion of functional dependence of the drag coefficient on the geometry .

to 0.8 | X[c] | X[c] |
Average Relative Error (%)

Eq.(<ref>)    23.58 %

Eq.(<ref>)$a_1$ = $\pi$, $K$ =$K_p$      8.18 %

Chen et al.[41]    23.68%

Ouchene et al. [30]     8.69%

Holzer and Sommerfeld [38]     48.27%

Frohlich et al.[42]     8.046%

Relative average error of the predictions of  different predictive models with the respect to the results of Sanjeevi et al. [29] for prolate spheroids with low aspect ratios and  $\alpha$ = $0^{\circ}$.

Comparison  between the drag coefficient $C_D$ predicted by different predicted models for the case of prolate spheroid, for different Reynolds numbers, and low aspect ratio values, and $\alpha$ = $0^{\circ}$.

Table <ref> shows the relative error of the logarithmic-based models, and the available correlations with respect to the Sanjeevi et al. [29] correlation. The logarithmic-based equation Eq.(<ref>), and Chen et al. [43] correlation show a substantial difference from the results of Sanjeevi et al. [29]. While the logarithmic-based equation (Eq.(<ref>)$a_1$ = $\pi$, $K$ =$K_p$ ) shows a relative error of 8.18%, which is similar to the relative errors of Ouchene [41] and Frohlich et al. [42].  The predictions of Holzer and Sommerfeld [38] have the worst relative error of 48.27%. This over prediction may be attributed to the divergence of the $\Delta$ for the Holzer and Sommerfeld [38] correlation from those of Eq.(<ref>) as shown in Figure <ref>.

For values of $p_a$ greater than or equal to 2, the drag coefficient increases as the aspect ratio increases for low Reynolds numbers (10 and 50), as shown in Figure <ref>. Most predictive models follow this trend, except for the Chen et al.[43] correlation, which exhibits non-physical behavior. This is because the Chen et al.[43] correlation does not accurately account for the Stokes term. At higher Reynolds numbers,the logarithmic-based equation Eq.(<ref>) and the logarithmic equation
(Eq.(<ref>)$a_1 = 3.0,$ and $K = K_p$)  predictions are close to those of the Sanjeevi et al.[29] correlation. The predictions of the Sanjeevi et al.[29] and Ouchene[41] correlations show that the drag coefficient is independent of the aspect ratio for high Reynolds numbers, which is consistent with our  predictions  in our previous work [34] for prolate geometry. The Holzer and Sommerfeld [38] correlation fails to predict the interdependence of the drag coefficient on the aspect ratio at high Reynolds numbers because two of its terms, $\frac{3}{\sqrt{Re}}\frac{1}{\phi^{\frac{3}{4}}}$ and $0.42\times 10^{0.4\left(-\log(\phi)\right)^{0.2}}\frac{1}{\psi}$, show  aspect ratio dependency at high Reynolds numbers, indicating that the theoretical framework they use in their predictive model does not correctly reflect the physics of the problem.

Table <ref> show average relative error for  the different models, and correlations for the case   of Figure <ref>.  The logarithmic based Eq.(<ref>) that used the training data of a single aspect ratio performs significantly better then the Frohlich[42], and Holzer and Sommerfeld[38], which they use multiple  aspect ratio prolates to build their correlations. The logarithmic equation (Eq.(<ref>), $a_1$ = 3.0, $K$ = $K_p$) performs better than the other correlations. Also, we noticed from the different geometries we used so far in the paper that the more streamlined the body is, the smaller the value of the $a_1$ is  .

to 0.8 | X[c] | X[c] |
Average Relative Error (%)

Eq.(<ref>)    11.94 %

Eq.(<ref>) $a_1$ = 3.0, $K$ = $K_p$      3.15 %

Chen et al.[41]    25.07%

Ouchene et al. [30]     11.92%

Holzer and Sommerfeld [38]     46.71%

Frohlich et al.[42]     32.07%

Relative average error of the predictions of  different predictive models with the respect to the results of Sanjeevi et al. [29] for prolate spheroids with high aspect ratios, and $\alpha$ = $0^{\circ}$ .

Comparison  between the drag coefficient $C_D$ predicted by different predicted models for the case of prolate  shaped particle, for different Reynolds numbers, and high aspect ratio values, and $\alpha$ = $0^{\circ}$.

To gain a more comprehensive understanding of the accuracy of our predictions, we will compare them to the numerical results of Khair and Chisholm [69] for two aspect ratios of $p_a$ = 10, 100, over a range of Reynolds numbers up to 100. We will examine the variation of $\frac{C_D}{C_{D_s}}$ with Reynolds number, as shown in Figure <ref>. Both logarithmic equation Eq.(<ref>), and the logarithmic equation (Eq.(<ref>),$a_1$ = 3.0, $K$ = 1.22, 2.99), and the Frohlich correlation [42] closely match the numerical results of Khair and Chisholm for both aspect ratios. We did not compare our drag coefficient results to the analytical solution provided by Khair and Chisholm [69] as their solution is only valid for Reynolds numbers up to 20.

Figure <ref> compares $\Delta$ values obtained from the numerical results of Khair and Chisholm [69]. We calculated the $\Delta$ for the results of Khair and Chisholm [69] by digitally extracting discrete points from Figure 4 of their paper and then calculating the gradient. The $\Delta$ values from Khair and Chisholm's [69] numerical investigation are close to our results from the logarithmic equations. This indicates that the $\Delta$ values from the logarithmic equations Eq.(<ref>) and Eq.(<ref>) accurately reflect the behaviour of Navier-Stokes equations solutions and are not limited to any particular drag coefficient correlation. The fluctuations in the  $\Delta$  results of Khair and Chisholm [69] is due to the fact we calculated the gradient from discrete points which were not uniform, adding to that the values of $\Delta$ are small of the order less than $O(10^{-2})$.
The $\Delta$ values for the case of $p_a$ = 100, from Sanjeevi et al. [29] correlation  are significantly deviating from those of  Khair, and Chisholm [69], and this is the main reason that they over-predicted the values of drag coefficient as shown in Figure <ref> previsuly. Frohlich [42] correlation results show a slight dependency on the aspect ratio. However, this is inconsistent with the results of the direct numerical simulations of Khair and Chisholm [69].

Comparison of $\frac{C_D}{C_{D_s}}$ values from Khair and Chisholm's [69] numerical results and predictive models for prolate spheroids with aspect ratios $p_a$ = 10 and 100 at various Reynolds numbers for $\alpha$ = $0^{\circ}$.

Comparison of $\Delta$ values from Khair and Chisholm's [69] numerical results and from predictive models for prolate spheroids with aspect ratios $p_a$ = 10 and 100 at various Reynolds numbers, and for $\alpha$ = $0^{\circ}$.

We will test our theory for a flow perpendicular to the major axis $\alpha$ = $90^{\circ}$ for two prolate spheroids with $p_a$ values of 5 and 10. For this purpose, we will use the drag coefficient from the correlation of Sanjeevi et al. [29] at $Re$ = 300. The obtained values for $a_1$ are as follows: $a_1$ = 4.063 for $p_a$ = 5.0 and $a_1$ = 4.7 for $p_a$ = 15.0. We used a higher Reynolds number to evaluate the $a_1$ coefficient because Sanjeevi et al. [29] reported that their predictions overestimate the drag coefficient values in the Stokes regime compared to the analytical solutions.

The comparison of the drag coefficients for two selected aspect ratios is illustrated in Figure <ref>. In the low Reynolds number regime, the logarithmic equations are in close agreement with the correlation of  Frohlich et al. [42] that is based on numerical simulations and the general correlation of Holzer and Sommerfeld [38]. At high Reynolds numbers, the drag coefficients obtained from the logarithmic equations match perfectly with the results of Sanjeevi et al. [29], suggesting that the rate of change of the drag coefficient follows Eq.(<ref>) even under extreme flow conditions as at $\alpha = 90^{\circ}$. For the case of $\alpha = 0^{\circ}$, the $a_1$ coefficient was found to be constant, independent of the aspect ratio. However, for $\alpha = 90^{\circ}$, the $a_1$ coefficient strongly depends on the aspect ratio. In addition, increasing the aspect ratio increased the skin drag for $\alpha = 90^{\circ}$.

We will compare the results obtained from the logarithmic-based equation with the high-fidelity CFD simulations for a prolate spheroid with an aspect ratio of $p_a = 6.0$ and an orientation angle of $\alpha = 45^{\circ}$ with respect to the flow velocity, as described in the work of Jiang et al. [70]. To perform this comparison, we first calculated the parameter $a_1$ using the predicted value of the drag coefficient, $C_D$, obtained from the correlation developed by Sanjeevi et al. [29] at a Reynolds number of 300. Our calculations revealed that $a_1$ = 3.607.

We have presented our comparison results in Table <ref>. Our findings further confirm our theory, with only a slightly significant difference observed at the lowest Reynolds number. The advantage of using our approach is that it significantly reduces the number of simulations required for each geometry and orientation. The reduction in the number of simulations is because we only need a single result for the drag coefficient at a moderate Reynolds number to predict the evolution of $C_D$ for the specific geometry and orientation.

to 1.2 | X[c] | X[c] | X[c] |X[c] |
$Re$    Jiang et al.[70]    Sanjeevi et al.[29]  (%)    Eq.(<ref>), $a_1$ = 3.607,$K$ = 1.276(%)

1.691     1.662($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.691, CD2 = 1.662][2]$%)    1.514($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.691, CD2 =1.514][2]$%)

1.040    0.954($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.040, CD2 = 0.954][2]$%)    0.964($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.040, CD2 =0.964][2]$%)

0.770($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.805, CD2 =0.770][2]$%)    0.750($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.805, CD2 =0.750][2]$%)

Relative error between the predictive $C_D$ coefficient from the logarithmic equations , and other correlations with the numerical results of
Jiang  et al.[65] for different Reynolds numbers for the geometry of prolate with $p_a$ = 6.0, and $\alpha$ = $45^{\circ}$.

Comparison of $C_D$ values predicted by the logarithmic equations, and different correlations for the case of prolate spheroids with $p_a$ = 5, 10 and for $\alpha$ = $90^{\circ}$.

§.§.§ Spherocylinders

For the case of spherocylinders, we will test our theory against the data of Sanjeevi et al. [40] for a spherocylinder with an aspect ratio of 4.0 at two angles of attack, $\alpha = 0^\circ$ and $\alpha = 90^\circ$. We will also compare our results with the correlation of Feng and Michaelides [71] for $pa = 2.0$ at the same angles of attack as used by Sanjeevi et al. <cit.>.
The discrete data provided by Sanjeevi et al. [40] allows us to directly compare our theory with the approximate solution of the Navier-Stokes equations. To evaluate the coefficient $a_1$ for the two orientations, we will use the drag coefficient at $Re = 100$  from Sanjeevi et al. [40] data, and assign a zero value to the coefficient $K$. The obtained values for $a_1$ are 3.017 for $\alpha = 0^\circ$ and 4.045 for $\alpha = 90^\circ$.
The comparison between the logarithmic equations (Eq. (<ref>), $a_1 = 3.017$, $K = 0$) and (Eq. (<ref>), $a_1 = 4.045$, $K = 0$) and the numerical results of Sanjeevi et al. [40] is depicted in Figure <ref>. The logarithmic equations accurately predict the numerical results of Sanjeevi et al. [40] for both angles of attack. The value of $a_1$ indicates that when $\alpha = 90^\circ$, the spherocylinder particle behaves approximately like an infinite cylinder.

Recently Feng and Michaelides[71] provided a correlation for the drag coefficient for different flow orientations, from their numerical simulations. The correlation is complicated function in the Reynolds number, and the aspect ratio. We will evaluate $a_1$, at $Re$ = 100 from Feng and Michaelides[71] correlation. The comparison   between the logarithmic equations and the drag coefficients from the  Feng and Michaelides[71] are shown in Table <ref>. Our results compare well for the two orientation with the results Feng and Michaelides[71], only for the case of $Re$ = 300, where there is a slight deviation.

Comparison of $C_D$ values predicted by the logarithmic equations, with the numerical results of Sanjeevi et al.[40] for the case of a spherocylinder with aspect ratio of 4.0 for $\alpha$ = $0^{\circ}$, and  $90^{\circ}$.

to 1.2 | X[c] | X[c] | X[c] |X[c] |X[c]|
$Re$    Feng and Michaelides[71] $\alpha$ = $0^{\circ}$    Eq.(<ref>), $a_1$ = 3.142, $K$ =0, $\alpha$ = $0^{\circ}$  (%)   Feng and Michaelides[71] $\alpha$ = $90^{\circ}$   Eq.(<ref>), $a_1$ = 3.566,$K$ = 0,$\alpha$ = $90^{\circ}$ (%)

0.907     0.896($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.907, CD2 = 0.896][2]$%)    1.342   1.320($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.342, CD2 =1.320][2]$%)

0.683     0.683($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.683, CD2 = 0.683][2]$%)    1.108
1.108($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.108, CD2 =1.108][2]$%)

0.504     0.520($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.504, CD2 = 0.520][2]$%)    0.906
0.944($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.944, CD2 =0.906][2]$%)

0.445     0.418($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.445, CD2 = 0.418][2]$%)    0.796
0.869($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 0.796, CD2 =0.869][2]$%)

Relative error between the predictive $C_D$ coefficient from the logarithmic equations , with the correlation of Feng and Michaelides[71]
for different Reynolds numbers for the geometry of spherocylinder with $p_a$ = 2.0, and $\alpha$ = $0^{\circ}$, $90^{\circ}$.

§.§.§ Irregular non-spherical particles

Finally, we will implement our theory for irregular non-spherical particles. These types of particles are complicated to characterize geometrically. Recently, Lain et al. [72] developed a correlation for the drag coefficient of irregular non-spherical particles using sphericity as the main geometrical characterization parameter. They performed thousands of direct numerical simulations for random orientations of particles, considering Reynolds numbers up to 200 and a narrow range of sphericity between 0.7 and 0.95. They derived a correlation similar to Haider and Levenspiel [32], as shown in Eq. (<ref>), only the coefficients in Eq.(<ref>) are different .
The correlation proposed by Lain et al. [72] provides an opportunity to test the universality of the $\Delta$ parameter beyond regular-shaped particles. We will evaluate the coefficient $a_1$ at $Re = 100$ for two values of sphericity, $\phi = 0.7$ and $\phi = 0.8$, with $K = 1.0$. The results of the logarithmic equations are listed in Table <ref>, and our predictions for the drag coefficient align closely with the results of Lain et al. [72], confirming the universality of our theory.

The Reynolds number used in all correlations presented from oblate subsections until now is based on the volume-equivalent sphere diameter.

to 1.2 | X[c] | X[c] | X[c] |X[c] |X[c]|
$Re$    Lain et al.[72],$\phi$ = 0.8     Eq.(<ref>), $a_1$ = 3.630, $K$ =1.0, $\phi = 0.8$  (%)   Lain et al.[72],$\phi$ = 0.7   Eq.(<ref>), $a_1$ = 3.870,$K$ = 1,$\phi$ =0.7 (%)

1.893     1.867($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.893, CD2 = 1.867][2]$%)    2.181   2.108($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 2.181, CD2 =2.108][2]$%)

1.414     1.413($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.414, CD2 = 1.413][2]$%)    1.654
1.653($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.654, CD2 =1.653][2]$%)

1.236     1.231($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.236, CD2 = 1.231][2]$%)    1.458
1.471($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.458, CD2 =1.471][2]$%)

1.140    1.128($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.140, CD2 =1.128][2]$%)    1.354
1.369($\eval*{$-\frac{CD1-CD2}{CD1}100$}[CD1 = 1.354, CD2 =1.369][2]$%)

Relative error between the predictive $C_D$ coefficient from the logarithmic equations , with the correlation of Lain et al.[72]
for different Reynolds numbers for different irregular particles with $\phi$ = 0.8, and 0.7.

§.§ Conclusions

We showed that after a rigorous comparison with available results in the literature, the drag coefficient can be evaluated with high accuracy by Eq.(<ref>)  which was  initially derived for the case of a sphere but can be applied to any bluff body geometry by linearly adding a constant ($a_1$ coefficient)  which  depends on the shape and   orientation of the  body. The logarithmic Eq.(<ref>) is applicable for the case of incompressible uniform flow. We draw the following conclusions:

* One of the key contributions of this investigation is the novel finding that, in the inertial flow regime, the rate of change of the drag coefficient with respect to the Reynolds number is independent of the geometry or orientation of the object. This property has profound implications for drag coefficient prediction of different shapes, and we believe it is directly linked to the vorticity generated in the fluid by the solid boundary.

* We have discovered a strong correlation between the $a_1$ coefficient of the logarithmic equation and the frictional drag predicted by the boundary layer theory. Moreover, we have found that only using logarithmic equations we can get values of the $a_1$ coefficient that is linked to  boundary layer theory. The  link between the $a_1$ coefficient  and  boundary layer theory may be a result of a yet undiscovered theory in fluid mechanics.

* By utilizing Eq.(<ref>), we can significantly reduce the computational time required to determine the drag coefficient for a single body at a specific flow orientation over a wide range of Reynolds numbers. This efficiency arises from the fact that we only need to compute the drag coefficient at a single Reynolds number, and predict the rest by implementing the logarithmic equation Eq.(<ref>) .

* We showed that by  using the physics informed symbolic regression algorithm [34, 37] there is a possibility to obtain accurate solutions that can be compared well with  high fidelity numerical solutions by using  empirical data, and previous knowledge of analytical solutions.

§ ACKNOWLEDGEMENTS

The first author  thanks Dimitra Damianidou, and Mohamed Fakhri El Hasadi for the enlightening discussions about the subject.  Finally, the authors thank the European Research Council for its financial support under its consolidator grant scheme, contract no. 615096 (NonSphereFlow).

[1]
Vogel, S.
Life in moving fluids: the physical biology of flow-revised and
expanded second edition.
Princeton university press,  (2020).

[2]
Cohen, I., Whitehead, S. C., and Beatus, T.
Nature Reviews Physics 1(11), 638–639 (2019).

[3]
El Hasadi, Y. M. and Crapper, M.
Biomicrofluidics 10(6) (2016).

[4]
El Hasadi, Y. M. and Crapper, M.
Applied Thermal Engineering 127, 857–869 (2017).

[5]
El Hasadi, Y. M. and Crapper, M.
Journal of Molecular Liquids 313, 113548 (2020).

[6]
Hucho, W. and Sovran, G.
Annual review of fluid mechanics 25(1), 485–537 (1993).

[7]
Anderson, J.
EBOOK: Fundamentals of Aerodynamics (SI units).
McGraw hill,  (2011).

[8]
Eckert, M.
Ludwig Prandtl: A Life for Fluid Mechanics and Aeronautical
Springer,  (2019).

[9]
Winslow, J., Otsuka, H., Govindarajan, B., and Chopra, I.
Journal of Aircraft 55(3), 1050–1061 (2018).

[10]
Hansen, M. O.
Aerodynamics of wind turbines.
Routledge,  (2015).

[11]
Mahajan, V. V., Mehmood, J., El Hasadi, Y. M., and Padding, J. T.
Chemical Engineering Science: X 3, 100030 (2019).

[12]
Darrigol, O.
Worlds of flow: A history of hydrodynamics from the Bernoullis
to Prandtl.
Oxford University Press,  (2005).

[13]
Stokes, G. G.
On the effect of the internal friction of fluids on the motion
of pendulums, volume 9.
Pitt Press Cambridge, Cambridge, UK,  (1851).

[14]
Stewartson, K.
SIAM review 23(3), 308–343 (1981).

[15]
Prandtl, L.
In Verhandlg. III. Intern. Math. Kongr. Heidelberg,  484–491.
Math Congress,  (1904).

[16]
Day, M. A.
Erkenntnis 33(3), 285–296 (1990).

[17]
Blasius, H.
Grenzschichten in Flüssigkeiten mit kleiner Reibung.
University of Gottingen, Germany,  (1908).

[18]
Birkhoff, G. and Hydrodynamics, A.
study in Logic, Fact and Similitude.
Princeton Univ. Press,  (1960).

[19]
Wieselsberger, C.
Technical report,  (1922).

[20]
Brenner, H.
Viscous flows: the practical use of theory.
Butterworth-Heinemann,  (2013).

[21]
Achenbach, E.
Journal of Fluid Mechanics 54(3), 565–575 (1972).

[22]
Khan, M. H., Sooraj, P., Sharma, A., and Agrawal, A.
Experimental Thermal and Fluid Science 93, 257–271

[23]
Proudman, I. and Pearson, J.
Journal of Fluid Mechanics 2(3), 237–262 (1957).

[24]
Kaplun, S. and Lagerstrom, P.
Journal of Mathematics and Mechanics 6(5), 585–593

[25]
Thom, A.
Technical report, HM Stationery Office,  (1928).

[26]
Frössling, N.
Fysiografiska Saellskapets Handlingar 51(NACA-TM-1432)

[27]
Dennis, S. and Walker, J.
Journal of Fluid Mechanics 48(4), 771–789 (1971).

[28]
Cheng, W., Pullin, D., Samtaney, R., Zhang, W., and Gao, W.
Journal of Fluid Mechanics 820, 121–158 (2017).

[29]
Sanjeevi, S. K., Dietiker, J. F., and Padding, J. T.
Chemical Engineering Journal 444, 136325 (2022).

[30]
Ouchene, R.
Physics of Fluids 32(7), 073303 (2020).

[31]
Proudman, I.
Journal of Fluid Mechanics 9(4), 593–602 (1960).

[32]
Haider, A. and Levenspiel, O.
Powder technology 58(1), 63–70 (1989).

[33]
Clift, R.
Drops and Particles 117 (1978).

[34]
El Hasadi, Y. M. and Padding, J. T.
AIP Advances 9(11), 115218 (2019).

[35]
Brown, P. P. and Lawler, D. F.
Journal of environmental engineering 129(3), 222–231

[36]
Acrivos, A. and Goddard, J.
Journal of Fluid Mechanics 23(2), 273–291 (1965).

[37]
El Hasadi, Y. M. and Padding, J. T.
Chemical Engineering Science , 118195 (2022).

[38]
Hölzer, A. and Sommerfeld, M.
Powder Technology 184(3), 361–365 (2008).

[39]
Happel, J. and Brenner, H.
Low Reynolds number hydrodynamics: with special applications to
particulate media, volume 1.
Springer Science & Business Media,  (2012).

[40]
Sanjeevi, S. K., Kuipers, J., and Padding, J. T.
International Journal of Multiphase Flow  (2018).

[41]
Ouchene, R., Khalij, M., Arcen, B., and Tanière, A.
Powder Technology 303, 33–43 (2016).

[42]
Fröhlich, K., Meinke, M., and Schröder, W.
Journal of Fluid Mechanics 901 (2020).

[43]
Chen, Y., Jiang, P., Xiong, T., Wei, W., Fang, Z., and Wang, B.
Chemical Engineering Journal 424, 130391 (2021).

[44]
Flett, T. M.
The Mathematical Gazette 42(339), 38–39 (1958).

[45]
Schlichting, H. and Gersten, K.
Boundary-layer theory.
Springer, Berlin Heidelberg, Germany,  (2016).

[46]
Abraham, F. F.
The Physics of Fluids 13(8), 2194–2195 (1970).

[47]
Abadi, M., Agarwal, A., Barham, P., Brevdo, E., Chen, Z., Citro, C., Corrado,
G. S., Davis, A., Dean, J., Devin, M., et al.
arXiv preprint arXiv:1603.04467  (2016).

[48]
Blasius, H.
Boundary layers in liquids with low friction.
printed by BG Teubner,  (1907).

[49]
Melnik, R. and Chow, R.
NASA. Langley Res. Center Aerodynamic Analysis Requiring
Advanced Computers, Pt. 1  (1975).

[50]
Friedlander, D. J., Leake, C., Zhou, A., Wang, X.-Y., and Huynh, H.
Technical report,  (2021).

[51]
Janour, Z.
Technical report,  (1951).

[52]
Lisoski, D.
California Institute of Technology  (1993).

[53]
Bejan, A.
Convection heat transfer.
John wiley & sons,  (2013).

[54]
Najjar, F. M. and Balachandar, S.
Journal of Fluid Mechanics 370, 101–147 (1998).

[55]
Narasimhamurthy, V. D. and Andersson, H. I.
International Journal of Heat and Fluid Flow 30(6),
1037–1043 (2009).

[56]
Fage, A. and Johansen, F.
Proceedings of the Royal Society of London. Series A, Containing
Papers of a Mathematical and Physical Character 116(773), 170–197

[57]
Tian, X., Ong, M. C., Yang, J., and Myrhaug, D.
Journal of Fluids and Structures 49, 149–169 (2014).

[58]
Kravchenko, A. G. and Moin, P.
Physics of fluids 12(2), 403–417 (2000).

[59]
Lee, S. and You, D.
Journal of Fluid Mechanics 879, 217–254 (2019).

[60]
Raissi, M., Wang, Z., Triantafyllou, M. S., and Karniadakis, G. E.
Journal of Fluid Mechanics 861, 119–137 (2019).

[61]
Vlachas, P. R., Arampatzis, G., Uhler, C., and Koumoutsakos, P.
Nature Machine Intelligence 4(4), 359–366 (2022).

[62]
Takami, H. and Keller, H. B.
The Physics of Fluids 12(12), II–51 (1969).

[63]
Imai, I.
Note, no. BN-104  (1957).

[64]
Gerhart, P. M., Gerhart, A. L., and Hochstein, J. I.
Munson, Young and Okiishi's fundamentals of fluid mechanics.
John Wiley & Sons,  (2016).

[65]
Jiang, H. and Cheng, L.
Physics of Fluids 33(3), 034119 (2021).

[66]
Sucker, D. and Brauer, H.
Wärme-und Stoffübertragung 8(3), 149–158 (1975).

[67]
Goldstein, S.
Modern developments in fluid dynamics: an account of theory and
experiment relating to boundary layers, turbulent motion and wakes,
volume 2.
Clarendon Press,  (1938).

[68]
Saha, A. K.
Physics of Fluids 16(5), 1630–1646 (2004).

[69]
Khair, A. S. and Chisholm, N. G.
Journal of Fluid Mechanics 855, 421–444 (2018).

[70]
Jiang, F., Gallardo, J. P., and Andersson, H. I.
Physics of Fluids 26(11), 113602 (2014).

[71]
Feng, Z. and Michaelides, E.
Available at SSRN 4435415 .

[72]
Lain, S., Castang, C., García, D., and Sommerfeld, M.
Available at SSRN 4317628 .