3.1. Force coefficients

Extensive 2-D DNS is conducted to quantify the applicable bound $K$ and $\beta$ values of the S–W solution over $\beta$ from 1 to $10^6$ and $K$ from 0.01 to 20. The $K\text {--}\beta \text {--}\varLambda _K$ map based on the present numerical results and (1.2) is presented in figure 2. The DNS results of $C_d$ agree extremely well with the S–W solution at low $K$ values, forming a significant contrast to the results shown in figure 1. The parameter space covered by $\{\varLambda _K\} \approx 1$ is much larger than its experimental counterpart of $\{\varLambda _K\}_\textrm {S} \approx 1$ shown in figure 1. The iso-lines of $\varLambda _K = 1.05$, 1.10 and 1.50, determined based on the DNS results, are plotted in figure 2 to provide a quantitative measure of applicable $K$ and $\beta$ bound values of the S–W solution. Taking $\varLambda _K = 1.05$ as an example, the applicable upper-bound $K$ values is approximately 0.8 for $\beta > 10^2$ and increases with decreasing $\beta$ for $\beta \leq 10^2$. The surprise observation (ii) from figure 1 is not observed in figure 2. The applicable $K$ and $\beta$ bound values of the S–W solution, e.g. those on the iso-line of $\varLambda _K = 1.05$, are much wider than those suggested by Wang (Reference Wang1968). The variation trend of $\varLambda _K$ with $K$ and $\beta$ for $\varLambda _K > 1.05$ is similar to that of $\varLambda _K = 1.05$ and is more sensitive to $K$ than $\beta$. The iso-lines of $\varLambda _K$ shown in figure 2 can be approximately fitted by

(3.1)\begin{equation} \{\beta\}_{\varLambda_K}=\frac{A_0}{(K-K_\infty)^{2}}, \quad \textrm{where }\ 1.05 \leq \varLambda_K \leq 1.5,\ 1 \leq \beta \leq 10^6. \end{equation}

The fitting parameters, $K_\infty$ and $A_0$ in (3.1), are selected based on the asymptotic $K$ values at $\beta \sim 1$ and $\beta \sim \infty$ as $A_0= 0.6 exp{(4\varLambda _K)}-35.326$ and $K_{\infty }=-571.5 exp{(-6.16\varLambda _K)}+1.7$.

Figure 2.

Distributions of $\varLambda _K$ over the $K\text {--}\beta$ parameter space. Discrete symbols are the present 2-D DNS cases, while the filled colour, with labels shown on the top right corner, in each symbol indicates the ratio of $C_d$ between DNS and the S–W solution, $\varLambda _K$. The dashed lines are the iso-$\varLambda _K$ lines for $\varLambda _K = 1.05, 1.1$ and 1.50 using correlations of $\{\beta \}_{\varLambda _K}= A_0/({K-K_{\infty }})^{2}$, where $A_0= 0.6 exp{(4\varLambda _K)}-35.326$ and $K_{\infty }=-571.5 exp{(-6.16\varLambda _K)}+1.7$. The straight blue and green solid lines indicate the Hall line ($K_h$) and Sarpkaya line ($K_r$) provided by Hall (Reference Hall1984) and Sarpkaya (Reference Sarpkaya2002), respectively. The red line indicates the limiting line of $\beta =1/K^2$. The grey shaded areas below $K_h$ line are $\{\varLambda _K\} \approx 1$ region reproduced based on the results given in Sarpkaya (Reference Sarpkaya2006a) and Sarpkaya (Reference Sarpkaya2010).

A detailed comparison of the present DNS results with the S–W solution is provided in figure 3 at $\beta = 1035$ and 11 240 over a range of $K$ values, where experimental results by Sarpkaya (Reference Sarpkaya1986) are available. Since the present focus is on low $K$ values, our discussions on the results shown in figure 3 are limited to $K < 1.8$ only with the following general observations:

1. (i) The present 2-D DNS results of $C_d$ values agree fairly well with the S–W solution at both $\beta = 1035$ and $11\ 240$ with $\varLambda _K \leq 1.05$ for $K \leq {\sim }0.8$ and the deviation from the S–W solution increases with increasing $K$ for $K > 0.8$.

2. (ii) The experimental $C_d$ at $\beta = 1035$ agrees well with the S–W solution and DNS results for $K \leq 0.7$, deviates considerably from the S–W solution and DNS results for $K > 0.7$, reaches a minimum value at $K \approx 1.6$ and increases with increasing $K$ for $K > {\sim }1.6$.

3. (iii) The experimental $C_d$ at $\beta = 11\ 240$ shows significant departures from the S–W solution and DNS results within the range of for $K > 0.8$. The ratio between $\{C_d\}_{Exp}$ and $\{C_d\}_{\text {S--W}}$ increases from ${\sim }1.2$ at $\beta = 1035$ to ${\sim }5$ at $\beta = 11\ 240$, as shown by the grey circles in figure 3(b).

4. (iv) The present 2-D DNS results of $C_m$ values agree fairly well with the S–W solution and experimental results at both $\beta = 1035$ and $11\ 240$.

Figure 3.

Distributions of $C_d$ and $C_m$ with respect to $K$ at (a) $\beta = 1035$ and (b) $\beta = 11\ 240$. The grey lines, black solid symbols and grey hollow symbols represent the force coefficients obtained from the S–W solution, present 2-D DNS and experimental results (Sarpkaya Reference Sarpkaya1986) respectively. The grey and blue dashed lines are the critical $K$ values at corresponding $\beta$ values obtained in terms of the $K_r$ and $K_h$ lines (Hall Reference Hall1984; Sarpkaya Reference Sarpkaya2002).

Corresponding instantaneous vorticity contours to the time instant marked as a filled circle in the inset are depicted in figure 4 for cases at $K = 0.1, 1.2, 1.8$ and 2.0 and $\beta = 1035$, to illustrate the variation of general flow features at different $K$ values. The symmetric features of 2-D flows described here use the nomenclature proposed by Elston, Blackburn & Sheridan (Reference Elston, Blackburn and Sheridan2006). The flow holds an $x$-reflection symmetry condition, i.e. $\omega _z(x, y,t) = -\omega _z(x, -y, t)$ ($\omega _z$ is the vorticity of 2-D flows), about $y/D = 0$ at small $K$ values in figure 4(a,b). As $K$ is increased, e.g. $K = 1.8$ in figure 4(c), the far tips of the shear layers on the cylinder become asymmetric and finally lead to vortex shedding and transition to turbulent flow at $K = 2.0$ in figure 4(d).

Figure 4.

Flow characteristics represented by vorticity contours for cases at (a) $(K, \beta ) = (0.1, 1035)$, (b) (1.2, 1035), (c) (1.8, 1035) and (d) (2.0, 1035) at a selected instant, where the sampling phase is marked as a red circle filled with black in the sinusoidal velocity signal in the inset of (a).

The increasing deviation of the present DNS results from the S–W solution for $K > 1.0$ is primarily due to the influence of flow separation, which will be discussed in detail later on. Sarpkaya (Reference Sarpkaya1986) attributed the large deviations of the experimental $C_d$ from the S–W solution observed in the case of $\beta = 1035$ to the development of 3-D instabilities and the hysteresis effect. We subsequently checked the hysteresis effect through separate DNS tests at $\beta = 1035$ and 1380 by gradually increasing and decreasing $K$ values and failed to identify a noticeable difference in $C_d$ from those DNS with zero initial conditions. The influence of 3-D instabilities on $C_d$ will be discussed in the subsection to follow. Since the minimum $K (= 0.8)$ value tested was relatively large for the case of $\beta = 11\ 240$, the large deviation observed did not attract much attention in Sarpkaya (Reference Sarpkaya1986). We suspect, based on our 2-D and 3-D DNS results, that the large deviation was induced by experimental uncertainties.

3.2. Three-dimensionality

To check the influence of 3-D effect on $C_d$, a limited number of 3-D DNS are conducted at $K$ values that are smaller than the corresponding $K$ values on the iso-line of $\varLambda _K = 1.5$ with $\beta = 200\text {--}20\, 950$. The results from 3-D DNS at $\beta$ from 200 to 20 950 are compared with the S–W solution and 2-D DNS results in table 1. The difference between $C_d$ values predicted by 3-D and 2-D DNS does increase with increasing $K$ and $\beta$ values. For instance, the difference increases from 0.50 % to 7.68 % as $K$ is increased from 0.8 to 1.4 in the case of $\beta = 20\, 950$. For a constant $K$ value, e.g. $K=1.2$, the difference increases from 0.30 % at $\beta = 1035$ to 2.69 % at $\beta = 20\, 950$. The present results appear to contradict with previous findings reported in the literature (Nehari, Armenio & Ballio Reference Nehari, Armenio and Ballio2004; Rashid, Vartdal & Grue Reference Rashid, Vartdal and Grue2011) that the 3-D effect does not affect the drag coefficient significantly. This contradiction arises because previous studies were mainly concerned with flows at small $\beta$ values, while the large differences between 2-D and 3-D DNS results are observed at either high $\beta$ or $K$ values in the present study. For the cases with small $\varLambda _K$ values investigated in the present study, e.g. $\varLambda _K \leq 1.1$, the difference between $C_d$ values predicted by 2-D and 3-D DNS is less than 2 %, suggesting 2-D DNS is sufficient for the parameter space bounded by $\varLambda _K \leq 1.1$. For this reason, the following discussions are based on 2-D results. Further research efforts are recommended to investigate the potential cause of the significant difference between 2-D and 3-D DNS at large $K$ and $\beta$ values.

Table 1.

Comparisons of $C_d$ values in the S–W solution, 2-D and quasi-3-D simulations, i.e. $\{C_d\}_{\text {S--W}}$, $\{C_d\}_{\text {2-D}}$ and $\{C_d\}_{\text {3-D}}$ respectively. The values in the brackets behind the $\{C_d\}_{\text {3-D}}$ are the relative differences compared with the corresponding 2-D results.

3.3. Flow separation

We speculate that the large $\varLambda _K$ values are mainly induced by flow separations around the cylinder surface because the S–W solution does not take into account the influence of boundary layer separation on the solution. To examine the influence of flow separation on the $\varLambda _K$ values, the spatio-temporal variations of flow separation on the upper surface of the cylinder are first quantified. The separation point is defined as the location where the vorticity on the cylinder surface changes sign and is measured by the separation angle ($\alpha _s$) relative to the front stagnation point of the cylinder, $(x/D, y/D) = (-0.5, 0)$ (see the inset in figure 5a).

Figure 5.

Variations of the instantaneous separation angle ($\alpha _s$) as a function of phase angle $\theta$ at (a) $\beta = 1035$ and (b) $K = 1$. Corresponding evolutions of free-stream velocity (red solid line) and acceleration (green dashed line) are plotted at the bottom of (b). (c) Evolution of boundary layer separation along the cylinder, represented by streamlines (line with arrows) and vorticity contours, at $(K, \beta ) = (1.2, 1035)$. The yellow and red contours represent positive and negative vorticities around the cylinder surface. The streamlines with blue and green colours represent the positive and negative signs of the shear stress.

The influence of $K$ on the variation of $\alpha _s$ with $\theta$ is examined for a number of $K$ values at $\beta = 1035$ in figure 5(a). Corresponding evolutions of free-stream velocity (red solid line) and acceleration (green dashed line) are plotted at the bottom of figure 5(b). For $K < 1.0$, the separation initially develops at $(x/D, y/D) = (0.5, 0)$ during the deceleration phase of the free-stream velocity ($\theta = 90^\circ$–$180^\circ$, where $\theta = 2{\rm \pi} t/T$ is the phase angle of incoming flow). It then propagates towards the upstream, featured by a decrease of $\alpha _s$ with increasing $\theta$, and terminates by merging with the upstream separation bubble for $\theta > 135^\circ$. The increase of $K$ induces an early onset of flow separation and a reduction of propagation rate, $\partial \alpha _s/\partial \theta$, in the phase space. For instance, the critical phase angle ($\theta _{cr}$) at which separation is initiated equals $122^\circ$ at $K = 0.4$ and decreases to $92^\circ$ at $K = 1.0$. A general trend observed in figure 5(a) is that $\alpha _s$ decreases almost linearly with $\theta$ for $\alpha _s > 90^\circ$ and drops sharply for $\alpha _s< 90^\circ$. Another flow separation bubble develops from the upstream stagnation point at $\theta \approx 135^\circ$, quickly climbs upwards and merges with the separation bubble originated from the downstream. The major flow separation features observed at different $\beta$ values with $K = 1$ in figure 5(b) are similar to those observed in figure 5(a). The average propagation rate of the separation point along the cylinder surface appears to be less affected by $\beta$ than $K$. Although $\theta _{cr}$ may be sensitive to $\beta$, the overall spatio-temporal development patterns of $\alpha _s$ at different $\beta$ values are similar. For example, $\alpha _s$ at $\beta =20\,950$ remains near $180^\circ$ for $70^\circ < \theta < 90^\circ$ before it propagates upstream at a similar rate to those at other $\beta$ values. The early occurrence of $\theta _{cr}$ in this case has little influence on $\varLambda _K$, as shown in figure 2.

The $\theta _{cr}$ value approaches $135^\circ$ at small $K$ values, e.g. $(K, \beta ) = (0.01, 1035)$ in figure 5(a), with a propagation rate close to infinity. This case can be considered equivalent to no flow separation because the separation point moves from the back stagnation point to the front stagnation point instantly at $\theta \sim 135^\circ$, which is identical to the critical phase angle in the Stokes solution of oscillatory flow over a flat plate when the wall shear stress changes sign(Stokes Reference Stokes1851). Furthermore, the leading term of the friction in-line force ($F_s$) in the S–W solution by Wang (Reference Wang1968) (see (3.6) in § 3.4) also shows a reversal of skin friction at $135^\circ$. This behaviour is expected as the flow with small $K$ and large $\beta$ is unlikely to feel the curvature effect induced by the cylinder. Since the flow velocity near the cylinder surface is in phase with the shear stress on the cylinder surface and leads the free-stream velocity by $45^\circ$, the reversal of the boundary layer flow direction at $\theta = 135^\circ$ is responsible for the above flow behaviour. It is worth noting that the point $(K, \beta ) = (0.01, 1035)$ is bounded by the $\beta = 1/K^2$ line in figure 2, where no flow separation is implied by Wang (Reference Wang1968).

The temporal evolution of flow separations at $(K, \beta ) = (1.2, 1035)$ is visualised through streamlines near the cylinder surface, as shown in figure 5(c). The blue and green colours of the streamlines represent the positive and negative signs of the shear stress. The red and orange contours represent positive and negative vorticities ($\omega _z$) around the cylinder at each phase angle. Although the separation point undergoes substantial development along the circumferential direction on the cylinder surface, the size of the separation bubble in the radial direction is rather limited in this case.

The above observations clearly show that flow separation occurs on the cylinder surface even at very low $K$ values, as shown in figure 5(a), which is different from the assumption made by Wang (Reference Wang1968). It is the spatio-temporal extent of the flow separation that is dependent on $K$ values. For example, the flow separation develops near the end of the deceleration phase of local flow near the cylinder (around $\theta = 135^\circ$) and vanishes shortly after the reversal of local flow near the cylinder at $K = 0.2$ and 0.4. The onset of flow separation advances forward in the phase space and vanishes at approximately the same phase after the reversal of local flow as $K$ is increased. To quantify the influence of $K$ on flow separation, a measure for the extent of flow separation over a half-oscillation period is proposed as follows:

(3.2)\begin{equation} \varGamma=\frac{1}{{\rm \pi}^2}\Big \{\int_{0}^{\rm \pi}({\rm \pi}-\alpha_s(\theta))\mathrm{d}\theta+\int_{0}^{\rm \pi}\alpha_s^*(\theta)\mathrm{d}\theta\Big \}, \end{equation}

where $a_s^*$ and $a_s$ represent the front and back separation bubbles, respectively; $\varGamma$ effectively measures a normalised spatio-temporal extent of flow separation (the shaded area shown in figure 5a) during a half-oscillation period. The results shown in figure 6(a,b) suggest that $\varGamma$ increases almost linearly with $K$ over a range of $\beta$ values and is less sensitive to $\beta$. This variation trend of $\beta$ with $K$ appears to correlate surprisingly well with the variation trend of $\varLambda _K$ with $K$ shown in figure 2. Correlation between $\varGamma$ and $\varLambda _K$ over a series of $K$ and $\beta$ values is plotted in figure 6(c) and can be represented by

(3.3)\begin{equation} \varLambda_K \sim 1+101.03\varGamma^{3}. \end{equation}

Figure 6.

Variations of the extent of flow separation $\varGamma$ with respect to (a) $K$ and (b) $\beta$ values. (c) Correlations between $\varGamma$ and $\varLambda _K$ for cases shown in (a,b).

The results shown in figure 6(c) suggest that $\varLambda _K$ correlates extremely well with $\varGamma$. The larger the $\varGamma$, the more influence it has on $\varLambda _K$. The influences of flow separation on pressure and shear stresses on the cylinder surface and hence the drag force are further explored below.

3.4. Time-dependent in-line forces on the cylinder surface

The instantaneous in-line force acting on the cylinder under the oscillatory flow condition ($U=U_m \sin \theta$) is comprised of three components (Stokes Reference Stokes1851; Wang Reference Wang1968; Bearman et al. Reference Bearman, Downie, Graham and Obasaju1985). The first component is the inviscid inertia force, $F_0$, induced by the acceleration of the outer flow, which can be estimated through the potential flow theory as,

(3.4)\begin{equation} F_{0}=\frac{{\rm \pi}^2}{K}\rho U_m^2 D \cos \theta.\end{equation}

The second component considers the contribution of viscous interactions of boundary layers with the cylinder. The viscous boundary layer on the cylinder surface affects the in-line force in two ways (Stokes Reference Stokes1851; Wang Reference Wang1968). The boundary layer profiles determine the distribution of skin friction $\tau _w$ along the cylinder surface. The shear force, $F_s$, can be determined by integrating $\tau _w$ along the cylinder surface as

(3.5)\begin{equation} F_s = \int_{0}^{2{\rm \pi}} \tau_w \sin \alpha\, \textrm{d}\alpha,\end{equation}

where $\alpha$ measures the angle from the cylinder surface relative to $(x/D, y/D) = (-0.5, 0)$. As the growth of the boundary layer is not uniform over the surface of the cylinder due to the curvature of the body, it induces a perturbation to the pressure distribution around the cylinder surface which subsequently alters the in-line force as well (Stokes Reference Stokes1851; Wang Reference Wang1968). This force is denoted as the viscous force, $F_v$, in the present study. In the analytical solution given by Wang (Reference Wang1968), $F_s$ and $F_v$ are identical and can be calculated as

(3.6)\begin{align} F_{v}=F_{s}&=\frac{{\rm \pi}^2}{K}\rho U_m^2 D \left[({\rm \pi}\beta)^{{-}1/2}(\sin \theta+\cos \theta) +({\rm \pi}\beta)^{{-}1}\sin \theta \vphantom{\frac{{\rm \pi}^2}{K}}\right.\nonumber\\ &\quad \left.+\frac{1}{4}({\rm \pi}\beta)^{{-}3/2}(\cos \theta-\sin \theta)\right]. \end{align}

The third component of the in-line force is induced by extensive separations of viscous boundary layer flow and the generation of vortices (Bearman et al. Reference Bearman, Downie, Graham and Obasaju1985), which modify the distributions of pressure and shear stress around the cylinder. This force is similar to the form drag in steady flow past a stationary cylinder and is denoted as $F_p$ in the present study. The value of $F_p$ becomes significant as $K$ exceeds a critical value, depending on $\beta$.

Since $F_v$, $F_p$ and $F_i$ are determined by the pressure distribution around the cylinder surface, they cannot be separated easily in analysing DNS results. In order to quantify the influence of $F_p$ and to compare DNS results of $F_v$ with the S–W solution, the following approximations are made to quantify $F_v$ and $F_p$ in DNS. We assume $F_i$ in DNS at low $K$ value is identical to the inviscid solution, $F_0$ ((3.4)), and the remaining part, $F'= \int _{0}^{2{\rm \pi} } p_w \cos \alpha \, \textrm {d}\alpha - F_0$, constitutes an approximation to $F_v + F_p$. Here, $F'$ is expected to reduce to $F_v$ and be identical to $F_s$ when the spatio-temporal extent of flow separation is insignificant at small $K$ values. The contribution of $F_p$ to $F'$ will be large when the spatio-temporal extent of flow separation is significant at large $K$ values. To validate the above understanding, corresponding variations of $\{F_s\}_{DNS}$ and $\{ F'\}_{DNS}$, based on DNS results at ($K, \beta ) = (0.1, 1035$), (1.4, 1035) and (2.0, 1035), are depicted in figure 7(a–c), together with $\{F_s\}_{\text {S--W}}$ based on the S–W solution ((3.6)). The red solid line represents the scaled free-stream velocity. As expected, evolutions of $\{F_s\}_{DNS}$ and $\{F'\}_{DNS}$ are identical and agree very well with $\{F_s\}_{\text {S--W}}$ for the case at ($K, \beta ) = (0.1, 1035$) in figure 7(a). The good agreement observed in this case is because the contribution of flow separation is negligible ($F_p \sim 0$) at $K = 0.1$. The value of $F'$ reduces to $F_v$ and is identical to the S–W solution ((3.6)). As $K$ is increased to 1.4 and 2.0, $\{F_s\}_{DNS}$ is slightly larger than $\{F_s\}_{\text {S--W}}$. The value of $\{F'\}_{DNS}$, on the other hand, deviates significantly from $\{F_s\}_{\text {S--W}}$ both in amplitude and phase. The large deviations observed in those two cases are suspected to be primarily due to the $F_p$ induced by flow separation. Since $F_v = F_s$ based on the S–W solution, the approximation of $F'\approx F_v + F_p \approx F_s + F_p$ would allow us to infer the influence of flow separation on $F_p$ by quantifying the influence of flow separation on $\{F_s\}_{DNS}$ over a half of flow oscillation period, through a measure proposed in the present study:

(3.7)\begin{equation} \varPsi=\frac{\displaystyle \int_{0}^{T/2}\{F_s\}_{DNS}\, \textrm{d}t}{\displaystyle \int_{0}^{T/2}\{F_s\}_{\text{S--W}}\, \textrm{d}t}. \end{equation}

The variations of $\varPsi$ values as a function of $K$ at $\beta = 1035$ and 20 950 shown in figure 7(d) suggest the spatio-temporal extent of flow separation has little effect on $F_{s}$ along the cylinder surface. For instance, flow separation induces only a $\varPsi \sim 1.06$ difference on $F_s$ at $(K, \beta ) = (2.0, 1035)$, whereas the measured $C_d$ is around $\varLambda _K = 1.43$ times of that by the S–W solution. The large $\varLambda _K$ value in this case is clearly caused by the form drag component $F_p$.

Figure 7.

Time evolutions of $F_{s}=\int _{0}^{2{\rm \pi} } \tau _w \sin \alpha \, \textrm {d}\alpha$ (orange short-dashed line) and $F'= \int _{0}^{2{\rm \pi} } p_w \cos \alpha \, \textrm {d}\alpha - F_0$ (blue dashed line) obtained in terms of DNS results at (a) $(K, \beta ) = (0.1, 1035)$, (b) (1.4, 1035) and (c) (2.0, 1035) over one flow oscillation period. The grey solid lines represent the corresponding shear force calculated based on (3.6) in the S–W solution, i.e. $\{F_s\}_{\text {S--W}}$. The red solid line represents the scaled free-stream velocity $U$. (d) Variations of $\varPsi$ values for cases at $\beta = 1035$ and 20 950.

3.5. General form of the Morison equation

The validity of the conventional Morison equation (1.1a–c) is questionable under the condition of low $K$ values where the drag component of $F_v+F_s$ dominates the total drag. The first term in (1.1a–c), proportional to $U^2$, represents the form drag component $F_p$ (Sumer Reference Sumer and Fredsøe1997). The S–W solution shows that $F_s+F_v$ is inversely proportional to $U$, rather than to $U^2$. Based on the results shown in § 3.4, $F_s+F_v$ generally dominates the total drag prior to the onset of extensive flow separations on the cylinder surface and $F_p$ becomes important after the extent of flow separation increases beyond a certain level. It is plausible to infer that there will be a parameter space where $F_s+F_v$ is comparable to $F_p$. Since $F_s$ is in phase with $\partial u_\alpha /\partial r |_{r=0.5}$, where $\alpha$ and $r$ are the polar coordinates, it contributes to the drag only. The solution of the laminar oscillatory flow above a flat plate shows that $\tau _w$, $u_\alpha$ and $\partial u_\alpha /\partial r |_{r=0.5}$ on the plate lead the free-stream velocity $U$ by ${\rm \pi} /4$. It is not difficult to deduce from (3.6) and figure 7(a) that this relationship also applies to oscillatory flow around a circular cylinder, where both $F_s$ and $F_v$ are linearly proportional to $U$ with a phase lead of ${\rm \pi} /4$.

Based on above reasoning, the use of (1.1a–c) under small $K$ conditions is fundamentally incorrect because it cannot account for the linearity and phase lead of $F_s+F_v$ to $U$. Applications of (1.1a–c) will not only lead to an incorrect $C_d$ but also poor predictions of time-dependent drag. To rectify the problem, a general form of the Morison equation is proposed in the present study as

(3.8)\begin{equation} F_x=\frac{1}{2}\rho D C_{d1}U_m^2\sin(\theta+{\rm \pi}/4)+\frac{1}{2}\rho D C_{d2}U_m^2\sin(\theta)|\sin(\theta)|+\frac{\rm \pi}{4}\rho D^2C_m\frac{\textrm{d}U}{\textrm{d}t}, \end{equation}

where $U=U_m\sin {\theta }$, $C_{d1}$ and $C_{d2}$ are the viscous and form drag coefficients, respectively. The first, second and third terms of (3.8) represent $F_s + F_v$, $F_p$ and $F_i$ respectively. The above equation recovers to the conventional Morison equation (1.1a–c) when $F_p$ significantly outweighs $F_s + F_v$ at large $K$ values. Assuming the form drag component is negligible under the applicable condition of the S–W solution, we get the analytical form of $C_{d1}$ based on the S–W solution as $\{C_{d1}/C_d\}_{\text {S--W}} = 16/(3\sqrt {2}{\rm \pi} ) \approx 1.2$ and the asymptotic form of $C_{d1}$ for $\beta \gg 1$ reads $\{KC_{d1}\sqrt {\beta }\}_{\text {S--W}} \approx 31.49$. This outcome suggests that (1.2) underestimates the peak value of the hydrodynamic damping force by 1.2 times. Since three terms in (3.8) have ${\rm \pi} /4$, 0 and ${\rm \pi} /2$ phase differences relative to $U$, respectively, a specific procedure is needed to determine the hydrodynamic coefficients in (3.8), based on the total in-line force, $F_x$, measured in DNS or physical experiments.

The following three options are available to determine the force coefficients in (3.8) based on the numerical results from DNS:

Method I, calculate $F_s+F_v$ and hence $C_{d1}$ based on the S–W solution and determine $C_{d2}$ and $C_m$ based on ($\{F_x\}_{DNS}-\{F_s+F_v\}_{\text {S--W}}$).

Method II, let $F_s + F_v=2F_s$, estimate $F_s$ based on DNS results and determine $C_{d2}$ and $C_m$ using ($\{F_x-2F_s\}_{DNS}$).

Method III, calculate $C_{d1}$ and $C_m$ based on the S–W solution and determine $C_{d2}$ based on ($\{F_x\}_{DNS}-\{F_s+F_v+F_i\}_{\text {S--W}}$).

The above three methods are expected to be identical at small $K$ values and yield different results at intermediate $K$ values. To evaluate the goodness of fit by using those methods, cases at $\beta = 1035$ over a range of $K$ values from 0.1 to 8.0 are used to determine a sensible choice of methods I, II and III. Corresponding distributions of $C_{d2}$ and $C_m$ based on methods I–III, together with those using the conventional Morison equation in (1.1a–c), denoted as method 0 hereafter, are shown in figure 8(a). Figure 8(b) shows the variation of a goodness-of-fit parameter $\zeta$ (proposed by Justesen Reference Justesen1989) with $K$. The goodness-of-fit parameter $\zeta$ is determined based on the measured and predicted phase-averaged drag force over one period $T$ as

(3.9)\begin{equation} \zeta = \frac{\displaystyle \int_{T}(\{F_d\}_m-\{F_d\}_p)^2\,\textrm{d}t}{\displaystyle\int_{T}\{F_d\}_m^2\,\textrm{d}t}, \end{equation}

where $\{F_d\}_m$ and $\{F_d\}_p$ are the time-dependent measured and predicted drag components. The lower the $\zeta$, the better the prediction is. The results suggest that method II yields the best fit and the method 0 leads to the worst fit for $K < {\sim }2.0$ at $\beta = 1035$. The results by method I are only slightly worse than those by method II at $K < {\sim }2.0$. The above results are not totally surprising because method II allows the influence of flow separation to be reflected in all three coefficients, while method I ignores the influence of flow separation on $C_{d1}$. Although method I leads to a slightly worse prediction than method II, method I is recommended to quantify $C_{d1}$, $C_{d2}$ and $C_{ m}$ in (3.8) due to the convenience of using the analytical solution of $\{C_{d1}\}_{\text {S--W}}$ in the equation.

Figure 8.

Distributions of (a) $C_{d2}$ (solid circles) and $C_m$ (solid squares), and (b) goodness-of-fit parameter $\zeta$, as a function of $K$ at $\beta = 1035$ based on DNS results. The red dashed line, blue dotted line, grey solid line and green dash dotted line in (a,b) represent the variations using methods 0–III, respectively.

Reconstructions of temporal variations of the drag force based on (1.1a–c) (method 0) and (3.8) (method I) over one $T$ at $\beta = 1035$ are illustrated in figure 9. Four cases at $K = 0.1\text {--}4.0$ are selected to depict the performance of the proposed general form of the Morison equation. By comparing the absolute difference of the measured and predicted drag components, i.e. $\{F_d\}_m - \{F_d\}_p$, (3.8) significantly outperforms (1.1a–c) for $K < {\sim }2.0$. At small $K$, e.g. $K = 0.1$ in figure 9(a), the reconstructed drag agrees fairly well with both measured drag values in DNS and the S–W solution ($\{F_s+F_v\}_{\text {S--W}}$). Due to the development of the form drag component, the reconstructed drag through method I deviates from the S–W solution as $K$ is increased (figure 9 b,c). The general form Morison equation slightly outperforms the conventional one for $K > {\sim }2.0$ (figure 9 d), suggesting (3.8) gradually reduces to (1.1a–c) as $K$ is increased.

Figure 9.

Time evolution of measured drag, $\{F_d\}_m$, and reconstructed drag, $\{F_d\}_p$, at (a) $(K, \beta ) = (0.1, 1035)$, (b) (1.4, 1035), (c) (2.0, 1035) and (d) (4.0, 1035) over one period. Corresponding $\{F_d\}_m$, $\{F_d\}_p$ and $\{F_d\}_m - \{F_d\}_p$ obtained through method 0 (red) in (1.1a–c) and method I (blue) in (3.8) are shown as thick solid, dashed and dotted lines.

Variations of $C_{d1}$ and $C_{d2}$ with $K$ and $\beta$, obtained by using (3.8) through method I, are shown in figure 10(a). As expected, opposite trends of $C_{d1}$ and $C_{d2}$ are clearly observed as $K$ is increased. It is seen in figure 10(a) that $C_{d1}$ is approximately 10 times larger than $C_{d2}$ for the cases of $K < {\sim }1.0$ over a wide range of $\beta$, which means the second term on the right-hand side of (3.8) can be ignored. Both $C_{d1}$ and $C_{d2}$ are equally important at intermediate $K$, e.g. $0.1< C_{d1}/C_{d2} < 10$ at $K= 1.0\text {--}8.0$, where the newly proposed general form of the Morison equation should be utilised. Distributions of $\zeta$ with respect to $K$ over a range of $\beta$ values for method 0 and method I are shown in figure 10(b). The general form of the Morison equation in (3.8) shows an excellent improvement over (1.1a–c) for $K < 2.0$. It is also slightly better than (1.1a–c) for $2.0 \leq K \leq 8.0$. For large $K$ values, e.g. $K > 8.0$, above which asymmetry and vortex shedding occur and the in-line force is largely drag dominated (Bearman et al. Reference Bearman, Downie, Graham and Obasaju1985; Sarpkaya Reference Sarpkaya2010), (3.8) recovers to (1.1a–c) where $C_{d2}$ is significantly larger than $C_{d1}$.

Figure 10.

Variations of (a) $C_{d1}$, $C_{d2}$ and (b) $\zeta$ as a function of $K$ for $\beta = 200 - 10^5$. The solid lines with square symbols in (a) are the viscous drag coefficient, $C_{d1}$, obtained through the S–W solution, while the triangles indicate the form drag coefficient, $C_{d2}$, determined by using method I.