coefficient-of-drag.org

Estimating the Coefficient of Drag for Real-Time Applications

\newcommand{\bm}[1] {{\bf{#1}}}

Keywords: fluid, physics, drag, video-game, casual-papers[fn:0]

[fn:0] This paper is written more informally so that it gets written rather than not.

Introduction

The aim of this paper is to estimate one part of the drag forces of a rigid body in a liquid. The rigid body’s surface is defined by a closed mesh of triangles and assumed to be convex. Let us begin with the drag equation from wikipedia \cite{drag}:

\begin{align} F = \frac{1}{2}~ρ~C_D~A~{v}^2 \end{align}

@@latex:\noindent@@ where \(ρ\) is the density of the fluid, \(C_D\) is the drag coefficient, \(A\) is the cross sectional area, and \(v\) is the relative velocity of the object in fluid.

Estimating Drag Coefficient

The drag coefficient \(C_D\) takes into account many factors and is often determined experimentally. In particular it takes into account the shape and direction of an object as shown in Figure fig:drag, which we’re going to try and estimate.

Shape Descriptors

Sphericity

Sphericity describes how close a particle’s shape is to that of a sphere. A sphere has a sphericity of 1 while anything that isn’t a sphere is less than 1.

\begin{align} ψ = \frac{π^\frac{1}{3} (6 V)^\frac{2}{3}}{A} \end{align}

@@latex:\noindent@@ Sphercity was used by \cite{SALMAN:1988ur} to estimate the coefficient of drag leading to this expression:

\begin{align} C_D &= 5.31 - 4.88 ψ \end{align}

@@latex:\noindent@@ This is a good approximation; however, we can see from figure fig:drag that the shape and its orientation matter. Let’s try to incorporate that difference if somewhat crudely.

Consider a Triangular Mesh

For a given mesh of $m$ triangles each triangle $i$ has three vertices in $\mathbb{R}^3$: $a_i$, $b_i$, $c_i$; a normal $n_i$; and an area $A_i$. The normal and area can be derived from the vertices.

\begin{align} n’ &= (b - a) × (c - a)
n &= \frac{n’}{|n’|} \ A &= \frac{|n’|}{2} \end{align}

Volume of Mesh

First we calculate the volume of a tetrahedron $T$ that’s formed from a triangle and a reference point $p$ that’s outside of the mesh.

\begin{align} T &= \frac{1}{6}~[(a - p) × (b - p)] ⋅ (c - p) \end{align}

@@latex:\noindent@@ Then we add all the tetrahedrons that face away from point $p$ and subtract all the ones that face it.

\begin{align} V &= ∑_{i = 1}^m \text{sign}(d ⋅ n_i) T_i \end{align}

\begin{align} \text{sign}(x) &= \begin{cases} -1 & x < 0 \ 1 & x \geq 0 \end{cases} \end{align}

Cross Sectional Area

We can estimate the cross sectional area $C$ for a given direction $d$ thusly:

\begin{align} cos_i &= d ⋅ n_i
C &= ∑_{i = 1}^m \text{step}(-cos_i)~|cos_i|~A_i \ \text{step}(x) &= \begin{cases} 0 & x < 0 \ 1 & x \geq 0 \end{cases} \end{align}

@@latex:\noindent@@ The area for each triangle is scaled by the proportion $cos_i$ which is the cosine of the angle between the normal $n_i$ and the direction $d$. The step function ensures only the faces directed toward $d$ are counted

Surface Area Facing a Direction $d$

We can also calculate the surface area $S$ facing a particular direction $d$.

\begin{align} S &= ∑_{i = 1}^m \text{step}(-cos_i)~A_i \end{align}

A Directional Sphericity

We want to be able to discriminate from a flat, dragging plane, and a smooth sphere. Let’s consider the case of a plane, its cross sectional area and its surface area will be the same. A sphere will have a larger surface area than cross sectional area. Let’s call our measurement $φ$.

\begin{align} φ = \frac{C}{S} \end{align}

@@latex:\noindent@@ For a plane $φ$ equals 1 and for a sphere $φ$ is greater than 1. For a cup with a concave surface, we may want $φ$ to be negative; however, this work will constrain itself to dealing with convex shapes only.

A Model

Going back to the equation by \cite{SALMAN:1988ur}.

\[ C_D = 5.31 - 4.88 ψ \]

@@latex:\noindent@@ It had the following model.

\[ C_D = a + b ψ \]

@@latex:\noindent@@ We’re going to amend that model with the following parameters that capture something of the shape of an object’s front $φ(d)$ versus its back $φ(-d)$.

\[ C_D = a + b ψ + c φ(d) + e φ(-d) \]

@@latex:\noindent@@ Using the coefficients of drag given in the above figure and the objects found in figure fig:objects, a fit was found using the Mathematica program in figure code:drag.

(* data entry is {volume, surface area, sphericity, front, back, C_D } *)

data = 
 { (* sphere *)         {0.51, 3.11, 1.00, 0.50, 0.50, 0.47},
   (* angled cube *)    {1.00, 6.00, 0.81, 0.35, 0.35, 0.80},
   (* cube *)           {1.00, 6.00, 0.81, 0.20, 0.20, 1.05},
   (* long cylinder *)  {1.55, 7.80, 0.83, 0.11, 0.11, 0.82},
   (* short cylinder *) {0.77, 4.67, 0.87, 0.20, 0.20, 1.15},
   (* cone *)           {0.26, 2.52, 0.78, 0.63, 1.00, 0.50},
   (* half sphere *)    {0.26, 2.32, 0.84, 0.50, 1.00, 0.42} };

model = a  + b sphericity + c front + e back;
fit = FindFit[data, model, {a, b, c, e}, {volume, surfaceArea, sphericity, front, back} ];
modelf = Function[{volume, surfaceArea, sphericity, front, back}, Evaluate[model /. fit]];
variables = data[[All, 1 ;; 5]];
original = data[[All, 6]];
estimated = Map[modelf @@ # &, variables];
absoluteDifference = Abs[estimated - original];

parameter	original	new
a	5.31	2.10
b	-4.88	-1.19
c	0	-0.25
e	0	-0.54

shape	$C_D^1$	$C_D^2$	$| C_D^1 - C_D^2|$
sphere	0.47	0.52	0.05
angled cube	0.80	0.86	0.06
cube	1.05	0.98	0.07
long cylinder	0.82	1.03	0.21
short cylinder	1.15	0.91	0.24
cone	0.50	0.48	0.02
half sphere	0.42	0.44	0.02

Other Considerations

The long and short cylinder are the most poorly estimated coefficients of drag. This is probably because the length of the shape is not really taken into account by any of the measures. A possible extension of this work could take into account a shape’s eccentricity in three dimensions. Or use something like the “Corey shape factor” as referenced by Tran-Cong et al. \cite{TranCong:2004cb} which is defined as the ratio of the shortest particle axis to the square root of the product of the other two axes.

bibliographystyle:plainurl bibliography:buoyancy.bib