Dynamics of heavy and buoyant underwater pendulums
Varghese Mathai, Laura Loeffen, Timothy Chan, Sander Wildeman
TThis draft was prepared using the LaTeX style file belonging to the Journal of Fluid Mechanics Dynamics of heavy and buoyant underwaterpendulums
Varghese Mathai † , Laura A. W. M. Loeffen † ,Timothy T. K. Chan , and Sander Wildeman ‡ School of Engineering, Brown University, Providence, RI 02912, USA. Physics of Fluids Group and Max Planck Center for Complex Fluids,Faculty of Science and Technology, University of Twente, P.O. Box 217, 7500 AE Enschede,The Netherlands. Department of Physics, The Chinese University of Hong Kong, Shatin, Hong Kong. Institut Langevin, ESPCI, CNRS, PSL Research University, 1 rue Jussieu, 75005 Paris,France.(Received xx; revised xx; accepted xx)
The humble pendulum is often invoked as the archetype of a simple, gravity driven,oscillator. Under ideal circumstances, the oscillation frequency of the pendulum is inde-pendent of its mass and swing amplitude. However, in most real-world situations, thedynamics of pendulums is not quite so simple, particularly with additional interactionsbetween the pendulum and a surrounding fluid. Here we extend the realm of pendulumstudies to include large amplitude oscillations of heavy and buoyant pendulums in afluid. We performed experiments with massive and hollow cylindrical pendulums inwater, and constructed a simple model that takes the buoyancy, added mass, fluid(nonlinear) drag, and bearing friction into account. To first order, the model predicts theoscillation frequencies, peak decelerations and damping rate well. An interesting effect ofthe nonlinear drag captured well by the model is that for heavy pendulums, the dampingtime shows a non-monotonic dependence on pendulum mass, reaching a minimum whenthe pendulum mass density is nearly twice that of the fluid. Small deviations fromthe model’s predictions are seen, particularly in the second and subsequent maximaof oscillations. Using Time- Resolved Particle Image Velocimetry (TR-PIV), we revealthat these deviations likely arise due to the disturbed flow created by the pendulum atearlier times. The mean wake velocity obtained from PIV is used to model an extra dragterm due to incoming wake flow. The revised model significantly improves the predictionsfor the second and subsequent oscillations.
Key words: nonlinear dynamical systems, particle/fluid flow, wakes
1. Introduction
The first known study on pendulums dates back to Galileo Galilei in 1605. The storygoes that Galileo observed the lamplighter pushing one of the swaying chandeliers in thePisa cathedral. Galileo timed the swings with his pulse and concluded that, althoughthe amplitude decreased, the time of each swing remained constant. Half a century laterthis observation inspired Christiaan Huygens to invent the pendulum clock, which until1930 has set the standard for accurate timekeeping. Based on Newton’s laws of motionHuygens could derive that an ideal, frictionless pendulum has a period T = 2 π (cid:112) L/g for † Equally contributed authors ; ‡ [email protected] a r X i v : . [ phy s i c s . f l u - dyn ] J a n Mathai, Loeffen, Chan, and Wildeman small amplitudes, where L is its length, and g the acceleration due to gravity (Huygens1986). However, as is well known to any clockmaker, real pendulums are often far fromideal. Even in vacuum, and when friction is negligible, one has to take into account thefinite amplitude and mass distribution of the pendulum in order to correctly predictits period. When a pendulum swings in a fluid such as air or water, additional fluidforces have to be taken into account. In fact, it was through careful observations ofdeviations from ideal pendulum motion that the nature and strength of many of thesefluid forces were brought to light. Concepts like buoyancy, added mass and fluid friction(viscosity) emerged concurrently with the study of pendulums swinging in different fluidenvironments (Stokes 1851).In recent years there has been a renewed interest in the use of simple pendulums toprobe fluid-structure interactions. Neill et al. (2007) conducted experiments on steel andbrass pendulums oscillating in different fluids. Their results agreed fairly well with themotion predicted by taking the buoyancy and added mass into account, implying thatviscous corrections played only a minor role in the oscillation frequency. In an extensionto this, Bolster et al. (2010) conducted experiments on a pendulum oscillating in a viscousfluid. In all of the above studies, the researchers focused on small amplitude vibrationsof dense pendulums in viscous fluids, where the dynamics could be reasonably modelledusing the ideal flow approximation, with the inclusion of a linear drag in the equationof pendulum motion. Many studies have focused on the flow disturbances induced by anoscillating body in a fluid. For one-dimensional oscillations at low to moderate Reynoldsnumbers, the effect of a body’s oscillation on the wake around it has been studied in somedetail (Tatsuno & Bearman 1990; Tatsuno & Karasudani 1993). Various types of wakeshave been quantified depending on the body oscillation frequency, velocity, fluid viscosityand the characteristic length scale of the problem. More recent studies have addressedthe influence of flow disturbances on pendula stability. Obligado et al. (2013) foundthat the equilibrium of a pendulum facing an incoming flow displays both bi-stabilityand hysteresis. Further, flow disturbances due to incoming turbulence can modify thestability map by changing the drag acting on the pendulum.Here we extend the scope of pendulum-fluid interaction studies to the regime of largeamplitude motions of heavy and buoyant cylindrical pendulums in a dense fluid. We haveperformed controlled experiments in which we simultaneously track the motion of thependulum and the flow field around it, in a regime where the Reynolds number (basedon the cylinder diameter) Re ∼ O (10 ) . The control parameter in our experiment isthe cylinder to fluid density ratio m ∗ ≡ ρ p /ρ f The ratio of the cylinder diameter tothe pendulum length
D/L , the fluid properties ( ρ f and ν ) and the release angle ( θ =90 ◦ ) are fixed. We vary m ∗ , covering both buoyant ( m ∗ < ) and heavy ( m ∗ > )pendulums. We build on the knowledge of buoyancy, added mass and inertial drag toconstruct a basic model for the pendulum motion, which we compare in detail to theexperiments. We find that a reduced added mass as compared to the ideal flow case leadsto better agreement between experiment and model predictions. We use insights fromparticle image velocimetry (PIV) measurements to model a higher-order fluid-structureinteraction, by which initially induced fluid motions affect the force on the pendulum atlater instants.
2. Experiments
The experiments were conducted using a cylindrical pendulum, with length L =203 mm, diameter D = 32 mm and span W = 300 mm. This resulted in a cylinderspan to diameter ratio W/D ≈ . The cylinder attached to the pendulum arm wasconstructed from polyvinyl chloride (PVC) pipe sealed with two plastic caps. Masses
800 mm mm )b()a((c) L = 203 mm D = 32 mm R = 3 mm L = 203 mm R = 3 mm D = 32 mm retaw ni muludnep tnayouBretaw ni muludnep yvaeH g g (d) Figure 1.
Schematics of the experimental setup for (a) heavy and (b) buoyant pendulums.(c) Cylinder mass was varied by inserting brass disks into the hollow cylinder. (d) Experimentalsetup used for performing PIV measurements. could be inserted inside the pipe (see Fig. 1(c)), which let us study both heavy andbuoyant pendulums in water.The pendulum was placed in a water tank (800 mm ×
400 mm ×
500 mm), and pivotedon a metal rod (radius R = 3 mm), which was fixed to the walls using suction cups. Thewire connecting the cylinder to the pivot axis had a diameter of 1.2 mm. This was thickenough to not bend during oscillations, while also providing low resistance to the flow.The cylinder ends remained at least . D away from the side walls. The lowest positionof the heavy cylinder during its motion was D above the base of the water tank, andsimilarly, the highest position of the buoyant cylinder during its motion was D belowthe water level. Schematics of the basic experimental arrangement for the heavy andbuoyant cases are shown in Fig. 1(a) & (b), respectively.The release mechanism consisted of a semicircular cavity with two curved metal wiresegments. The cylinder was held against the semicircular cavity, and locked in place bythe metal wires. The release was done by pulling the wire, which ensured a gentle releaseof the pendulum without creating disturbances in the water. By inverting these parts,the buoyant pendulums could also be released without disturbing the water. The releaseposition was at a ◦ angle from the vertical. A waiting interval of 15 minutes ensuredthat any residual flow from a previous run had damped out.A high-speed camera (Photron Fastcam 1024PCI) was used to record the pendulummotion. The cylinder movement was detected using a circle detection method (Mathai et al. <
1% of the pendulum diameter. Additionally, we performed Time-Resolved Particle Image Velocimetry (TR-PIV) to measure the flow field surrounding thependulum. The flow was seeded with fluorescent polystyrene tracer particles (diameter ≈ µ m). The seeding particle size was small enough to have both a small Stokesnumber and a small Stokes/Froude ratio (Mathai et al. Mathai, Loeffen, Chan, and Wildeman to have well-lit tracer particles all around the cylinder. A double-frame camera (ImagersCMOS) was used at a frame rate of 25 fps. The particle seeding density for PIV wassuch that an interrogation window contained around 6 particles, which is consideredgood for accuracy. The inter-pulse time ∆t was varied from 1–3 ms. ∆t was optimisedby ensuring that most of the tracer particles remained within the interrogation window,with a maximum movement ≈ / th of the window width. The PIV analysis was donewith two stage processing and 50 percent overlap using LaVision software.
3. Model equation of motion
The equation of motion of a pendulum oscillating in a fluid can be expressed byNewton’s law in the angular form:
I d θdt = (cid:88) τ, (3.1)Here, θ is the angular position of the pendulum with respect to the equilibrium position, I the moment of inertia of the system and τ the net torque acting on the system (seeFig. 2(a)). The pendulum mass and moment of inertia come mostly from the cylinder ata distance L from the rotation axis. Since L/D > , the cylinder moment of inertiaabout its mass centre I cm = mD / is negligible compared to the pendulum moment ofinertia mL . The effect of the pendulum accelerating through the fluid is modelled bymeans of an added mass m a . Therefore, the effective mass of the system can be writtenas m eff = m + m a . With these approximations, we model the system as a point masswith a moment of inertia: I p = m eff L .The torques acting on the system can be written as (cid:80) τ = (cid:80) r × F , where F are theindividual forces, and r the corresponding arm lengths. The forces acting on the systemare gravity F g ≡ ρ p V g , buoyancy F B ≡ ρ f V g , fluid drag F D and a bearing friction F f ≈ µ f F N (see Fig. 2(a)). Here, ρ p and ρ f are the particle and fluid densities, respectively, V is the cylinder volume, µ f is the friction coefficient and F N ≈ ( F B − F g ) cos θ is thenormal force acting at the bearing. The tension due to the centrifugal force is negligible.The Galileo number Ga ≡ (cid:112) gD | m ∗ − | /ν . The velocity scale in Ga is a gravitationalvelocity v g = (cid:112) gD | − m ∗ | . The Reynolds number is defined as Re = v max D/ν , where v max is the measured maximum velocity of the cylinder in the first swing. We note thatGa ∼ O (10 − ) gives a predictive estimate of the Re range in the experiment (seeTable 1). The drag force may be written as F D = ρ f A p C D v p , where C D is the dragcoefficient, A p is the projected area of the cylinder and v p ≡ Ldθ/dt is the instantaneouscylinder velocity (Batchelor 2000). By defining m ∗ ≡ ρ p /ρ f and ˜ t ≡ t (cid:112) g/L , the equationof motion can be written in non-dimensional form as m ∗ eff d θd ˜ t = − k sin θ − c (cid:12)(cid:12)(cid:12)(cid:12) dθd ˜ t (cid:12)(cid:12)(cid:12)(cid:12) dθd ˜ t − h | cos θ | sgn (cid:18) dθd ˜ t (cid:19) , (3.2)where m ∗ eff = ( m ∗ + m ∗ a ); k = | m ∗ − | ; c = C D A p L V ; and h = µ f | m ∗ − | RL .In order to solve eq 3.2, m ∗ a , C D and µ f have to be set. It is well known that addedmass coefficients m ∗ a for oscillating cylinders depend on a variety of parameters, includingthe frequency of oscillation, distance to boundaries, free surfaces, etc. Most studies havefocused on the case of cylinders near a free surface and for small oscillations (Dong1978; Konstantinidis 2013; Koo & Kim 2015; Tatsuno & Bearman 1990; Mathai et al. m ∗ a can deviate significantly (b)(a) t (s) θ (r a d ) -π/80π/8π/43π/8π/2 m* = 4.98 m* = 1.16
42 860 (p1) (p2) L R D Force diagram
Figure 2. (a) Schematic of the cylindrical pendulum at an instant when it swings to the right.Blue arrows indicates the forces on the pendulum: buoyancy F B , weight F g , drag F D , normalforce due to the pendulum arm F N , and a bearing friction F f . The dashed arrow and v p denotethe instantaneous velocity of the cylinder. (b) Angular position ( θ ) vs time for m ∗ = 4 . and m ∗ = 1 . cases. (p1) & (p2) – peak deceleration points. Table 1.
Ga, and Re for a selected number of mass density ratios m ∗ used in the experiment.Heavy Buoyant m ∗ × − × − from the two-dimensional potential flow added mass coefficient. Moreover, the fact thatthe cylinder span is finite will induce a three-dimensional (3-D) flow near the cylinderends, allowing some fluid to move to the sides. Therefore the cylinder motion is expectedto accelerate less fluid in the 3-D case as compared to 2-D potential flow. Nevertheless,to begin with, we choose the 2-D potential flow added mass coefficient of an infinitecylinder, m ∗ a = 1 . .The drag coefficient range in our experiments could be estimated from the range ofRe (see table 1). We begin with a simplified mean drag coefficient C D ≈ . based on theReynolds number range in our experiments (Lienhard 1966). Note that the influence ofReynolds number and vortex shedding on C D are not included at this stage, and will bediscussed later. µ f is expected to lie in the range [ . − . ] for lubricated steel-on-steelcontact. We estimated µ f ≈ . , based on tests performed using very heavy pendulumsin air, since in these cases the drag and added mass forces were at a minimum.Fig. 2(b) shows θ vs t predicted by the model for heavy pendulums with m ∗ = 1 . and m ∗ = 4 . . Clearly, the amplitude decay and oscillation frequency depend on the massdensity ratio. With this in mind, we explore the predictions for a range of m ∗ , coveringheavy and buoyant underwater pendulums.
4. Results
Model predictions
We first present some typical pendulum motions predicted by the model equation ofmotion. In Fig. 3(a), we show a contour plot of the angular position θ vs time t forcontinuous variation of m ∗ . We notice a clear asymmetry about m ∗ = 1 line, whichseparates the heavy cases from the buoyant cases. For instance, the oscillation frequencyis slightly higher for a buoyant cylinder as compared to a heavy cylinder that has identical Mathai, Loeffen, Chan, and Wildeman t ( s ) -π/32-π/640π/64π/32 θ (r a d ) (a) m* (b) π/2 ω p e nd s / d a r( ) -π-2π/3-π/30π/3 m* = 0.02 m* = 1.98 θ (rad) Figure 3. (a) Contour-plot of θ vs time t and mass-density ratio m ∗ . (b) Phase portraits fromnumerical solution of eq. 3.2 showing the evolution of angular velocity ω pend vs angular position θ for m ∗ = 0 . and m ∗ = 1 . cases. Note that these cases were chosen with identical driving: | F B − F g | . Phase portraits obtained for available experimental data (not shown here) showsimilar behaviour. driving | F B − F g | . On the buoyant side ( m ∗ < ), the oscillations damp faster withdecreasing m ∗ , owing to the decreasing inertia. Other aspects of the pendulum motionare altered in going from m ∗ > to m ∗ < . In Fig. 3(b), we show the phase portraits oftwo cases, one buoyant ( m ∗ = 0 . ) and other heavy ( m ∗ = 1 . ), with identical driving | F B − F g | . The low mass-density ratio of m ∗ = 0 . chosen here corresponds to a cylindermade from a very light material such as expanded polystyrene (Mathai et al. (2015)).The buoyant pendulum with ( m ∗ = 0 . ) accelerates quickly, reaching a high angularspeed in a short period of time. At the same time, owing to its low inertia, which in thiscase comes entirely from the fluid (or m ∗ a ), the motion damps out quickly. Therefore, thepeak angular displacement at first minimum θ is lower for the buoyant pendulum ascompared to the heavy pendulum ( m ∗ = 1 . ), while the peak angular velocity reached ω peak is higher for the buoyant cylinder.4.2. Model vs experiment
We now present a one-to-one comparison between the model and experiment. InFig. 4(a), we show the normalised frequency of oscillation. As expected, the frequencyof oscillation has a clear dependence on the mass ratio m ∗ . When the pendulum massdensity is close to the fluid density (or m ∗ → ), the frequencies are low. For heavycylinders, the frequency increases with m ∗ , and asymptotically approaches the frequencyof a large amplitude simple pendulum f ∗ π/ . The predictions of the model (eq. 3.2) aregiven by the blue and red dashed curves. The blue curve corresponds to the predictionwhen the potential flow added mass m ∗ a = 1 . is used. This slightly underpredicts theoscillation frequency. A best fit to the data is obtained for m ∗ a = 0 . . Thus, for C D = 1 . , m ∗ a = 0 . and µ f = 0 . , the frequency predictions of the model are in good agreementwith the experiments for both buoyant and heavy cylinders.It is interesting to compare the predictions of the model (Eq. 3.2) when the drag andthe bearing friction are ignored altogether. The solid green curve in Fig. 4(a) shows thefrequency prediction for large amplitude oscillations when the potential flow added massalone is accounted for. This underpredicts the oscillation frequency. However, when thesame undamped model is used for small amplitude oscillations, the predictions improvesignificantly (see solid blue curve). This is because the nonlinear drag is highly effectivein quickly reducing the oscillation amplitudes to modest values. A further improvement -1 Model ( m a * = 1.0)Model ( m a * = 0.53) Experiment
Model ( m a * = 1.0)Model ( m a * = 0.53) Experiment C D = 1.2 (Lienhard, 1966) )b()a( Undamped ( , m a * = 0.53)Undamped ( , m a * = 1.0)Undamped ( , m a * = 0.53) Figure 4. (a) Normalized oscillation frequency f ∗ from experiment and model (eq. 3.2). Here, f is averaged over the first four oscillations, and normalised by f sp = π (cid:112) gL . For the full equationof motion, using potential flow added mass m ∗ a = 1 . leads to an underprediction of f ∗ . For areduced m ∗ a = 0 . , the predictions of eq. 3.2 match well with the experiments. Remarkably,an undamped model for small amplitude simple pendulum, and with the same m ∗ a = 0 . ,reproduces the curve, providing evidence that the nonlinear drag plays only a weak role in theoscillation frequency f ∗ . The green curve shows that an undamped model with large initialamplitude θ = π/ underpredicts the oscillation frequency. (b) Peak pendulum deceleration, a max vs m ∗ . is obtained upon using a reduced added mass of m ∗ a = 0 . . Neill et al. (2007) had notedthat linear drag had only a small effect on the oscillation frequency. Here we found thatthe same holds for the nonlinear drag term.The frequency of oscillation represents an averaged quantity for the pendulum motionconsidered here. A more sensitive quantity would be the peak deceleration a max ofthe oscillating pendulum, which signifies the instant when the largest force acts onthe pendulum. Following the initial instants of release of the cylinder, a peak in thedeceleration is experienced when the cylinder reaches the end of the first swing, i.e. pointssuch as (p1) & (p2) in Fig. 2(b). In Fig. 4(b), we compare the experimentally measuredpeak deceleration with the predictions of the model (eq. 3.2). The peak decelerationchanges significantly with m ∗ , and the model predictions are in reasonable agreement withthe experimental measurements. Interestingly, the value of the added mass coefficienthas only a minor role in a max . Instead, the fluid drag has a leading role in the peakdeceleration of the pendulum. Thus, the frequency curve (Fig. 4(a)) provides a probethat singles out the added mass effect, while the deceleration curve (Fig. 4(b)) probesmainly the nonlinear drag’s effect.Figure 3(a) also provides insight into an interesting aspect of the oscillation decayfor heavy pendulums. For m ∗ > cases, the oscillations decay quicker with increasing m ∗ . This occurred despite the higher inertia of heavier pendulums, and seems counter-intuitive based on our common knowledge of pendulums oscillating in air. However, thistrend cannot continue to large m ∗ , since in the limit of very large inertia ( m ∗ (cid:29) )the decay rate should reduce with m ∗ . This implies the existence of an optimal heavypendulum ( m ∗ opt > ) for which one observes the quickest damping. We will provide anapproximate model to predict m ∗ opt , corresponding to the quickest damped pendulum.4.3. Amplitude decay
The initial damping is determined mainly by the nonlinear drag; therefore, we ignorethe bearing friction force in the equation of motion (eq. 3.2). This yields: m ∗ eff d θd ˜ t = Mathai, Loeffen, Chan, and Wildeman (d) (a) (c) (b)
Figure 5. (a) Time to decay 60% of the initial amplitude ˜ τ vs mass density ratio m ∗ .(b) Amplitude envelope θ m vs m ∗ after a time τ ref = 3 π (cid:112) L/g . For the heavy pendulums, anoptimal m ∗ and optimal f ∗ are visible in both experiment and simulations. (c) ˜ τ vs normalisedoscillation frequency f ∗ = f/f sp , where f sp = π (cid:112) gL . (d) θ m vs f ∗ at τ ref = 3 π (cid:112) L/g . Notethat (c) & (d) show heavy cases only. The basic model overpredicts the decay time and amplitudeenvelope for all cases. − k sin θ − c (cid:12)(cid:12) dθd ˜ t (cid:12)(cid:12) dθd ˜ t . In the absence of fluid drag ( c = 0 ), the solution is given by θ (˜ t ) = θ cos ω ˜ t , with ω = (cid:113) k/m ∗ eff , and θ the initial angle. If the damping constant c is small such that the pendulum amplitude is approximately constant over one period,we can assume that the solution is still approximately the same, but with a maximumamplitude θ m (cid:0) ˜ t (cid:1) which slowly decays in time. Consequently, the energy in the oscillationsevolves as E (cid:0) ˜ t (cid:1) ≈ k θ m (cid:0) ˜ t (cid:1) . The overall decay rate of the energy can be equated to theaverage work done by the drag force over one period T = 2 π/ω : dEd ˜ t = ω π (cid:90) πω F D (cid:0) ˜ t (cid:1) ˙ θ (cid:0) ˜ t (cid:1) d ˜ t (4.1)Using E (cid:0) ˜ t (cid:1) ≈ k θ m (cid:0) ˜ t (cid:1) and F D = − c ˙ θ , we obtain a differential equation forthe slowly varying envelope θ m ( t ) of the form: dθ m d ˜ t ≈ −B ( m ∗ ) θ m , where B ( m ∗ ) = c π (cid:113) | m ∗ − | ( m ∗ + m ∗ a ) . This yields the time to decay to a certain fraction ( − α ) of the initialamplitude θ as ˜ τ (1 − α ) = α/ (1 − α ) θ B ( m ∗ ) . For m ∗ > , and using the best fit m ∗ a = 0 . , B has a maxima at m ∗ opt ≈ . . Inother words, for a given θ the damping is the quickest when the pendulum is nearlytwice the density of the fluid. Alternately, one could plot the amplitude envelope θ m after a dimensionless time ˜ τ ref ≡ τ ref √ L/g for different m ∗ . This yields an expression θ m = B ( m )˜ τ ref +(1 /θ ) . In Fig. 5(a)–(b) we show comparisons between the model and experimentfor ˜ τ and θ m , respectively, vs mass-density ratio m ∗ . The black symbols denote themodel predictions by solving the full equation of motion (eq. 3.2). The black curve showsthe prediction by the approximate model described above. Note that we first fit anenvelope to the amplitude decay, using the maxima (peaks) of θ vs t curve, i.e a curve thatgrazes through the maxima of the amplitude curve. The optimal damping at m ∗ opt ≈ . is clearly visible in both the simulations and the reduced model. The blue data pointsare from experiment. A similar behaviour with a minima at an optimal m ∗ is visiblein the experiments. At the same time we note that both ˜ τ and θ m are lower in theexperiment. Further, the optimal m ∗ value in experiment is slightly lower as comparedto the model predictions. For very heavy pendulums ( m ∗ (cid:29) ), the damping B ∝ /m ∗ ,which leads to the linear relation ˜ τ ∝ m ∗ . This is also the commonly encounteredsituation for pendulums oscillating in air.The above discussed non-monotonic damping behaviour may be understood as theinterplay between the higher speeds achieved by the heavier pendulums, and the nonlineargrowth of the drag in proportion to the square of the speed. The phenomenon observedhere may hold some analogies to the added damping for oscillating cylinders (Dong 1978),which is usually expressed in terms of the oscillation frequency. Therefore, in Fig. 5(c)–(d)we plot ˜ τ and θ m , respectively, as a function of the normalised oscillation frequency f ∗ . The θ m plot shows a clear minima at f ∗ ≈ . . Similar to the m ∗ plots (Fig. 5(a)–(b)),the model overpredicts both ˜ τ and θ m , which occurs despite using the optimal addedmass coefficient m ∗ a = 0 . determined earlier. Clearly, the origin of these must lie in theinadequacy of the simplified drag model used ( C D = 1 . ). Nevertheless, the constant C D case can be viewed as the most basic model, which captures many essential features ofthe damping. Using a higher value of C D seems the logical choice. In the following wewill show that accounting for a few phenomenological effects can improve the predictions.We will introduce these modifications to the model and equation of motion.4.4. Drag corrections
Firstly, the influence of the instantaneous Re on C D can be taken into account.This leads to a faster decay at small amplitudes (or small v p ), since C D is higher atlower Reynolds numbers (Hoerner 1965). Secondly, the forces due to vortex sheddingbehind the cylinder could be modelled into the drag term. The drag coefficient in thepresence of vortex shedding may be approximated as C Dvs = C D (1 + k sin ω vs t ) , where k ∼ . (Govardhan & Williamson 2000), and ω vs = 2 π St v p /D is the vortex sheddingfrequency. Here St is the Strouhal number (Govardhan & Williamson 2005), v p thecylinder velocity and D the cylinder diameter. This leads to a marginal reduction inthe predicted peaks of the oscillations.Fig. 6(a)–(b) compares the experiment and model predictions for the angular position θ vs time t for a heavy pendulum with m ∗ = 4.98. The oscillation frequency matches wellbetween the experiment and the model. However, deviations are seen in the second peakof oscillation, marked as (r) in Fig. 6(b). The deviation persists for all of the subsequentoscillations. This suggests that the major factor causing the deviations between the modeland experiment is still missing in the equation of motion.4.5. PIV and wake corrections
To understand the origin of the deviation in the second and subsequent peaks ofthe oscillations, we employed high-speed PIV to quantify the flow around the cylinder0
Mathai, Loeffen, Chan, and Wildeman t (s) (a) (r) U x / v ω / ( v / D ) θ (r a d ) θ (r a d )
32 mm 32 mm32 mm
32 mm32 mm 32 mm t (s)
ExperimentOriginal model (q) (r)(q) (b)(e)(d)(c) U x / v U x / v ω / ( v / D ) ω / ( v / D ) (p) Figure 6. (a)Angular position θ vs time t for m ∗ = 4 . . The red curve shows the modelpredictions. (b) shows zoom-in of the second swing. The model overpredicts the maximumamplitude (point marked as (r) in (a) & (b)). (c) – (e) Normalised horizontal velocity U x /v ( left )and normalised vorticity ω/ ( v /D ) ( right ) at the three instants marked in subfigure (a) as (p) , (q) , and (r) , respectively. Here, v is the measured maximum speed of the pendulum. during its motion. Fig. 6(c) shows the normalised velocity and vorticity fields at a timeinstant when the cylinder is moving towards the right during its first swing (point (p) in Fig. 6(a)). The wake behind the cylinder is clearly visible, while the flow ahead ofthe cylinder has negligible vorticity i.e. nearly irrotational (Fig. 6(c)( right )). The wakeis unsteady and is expected to induce unsteady forces on the cylinder during its swing.One can expect that the mean drag on the cylinder at this instant is fairly predicted1 -1 -2 -1 slope ~ -0.8 -0.60.40.200.2-0.4 (a) (b) D D D Figure 7. (a) Schematic showing the estimation of mean wake velocity U f from a samplePIV flow field. A sample video of the PIV flow field is provided as supplemental material. U f iscomputed inside a D × D window located at a selected angular position θ sel . (b) Decay of wakevelocity U f in time at various angular positions θ sel . t pass is the time when the cylinder passed θ sel . The initial oscillation of the wake shows clear indication of vortex shedding. However, afterthis phase the wake decay is smooth. Inset shows the same plot on log-log scale. While at shorttimes U f shows oscillations, the long-time decay shows a power law decay. by the nonlinear drag term F D in the equation of motion (eq. 3.2). Consequently, thefirst swing of the pendulum is captured well by our model. Next we focus at a laterinstant in time (point (q) in Fig. 6(a) & (b)), when the cylinder has just completedits first swing and is returning to the left. Fig. 6(d) shows the velocity and vorticityfields at this instant. In this case, however, the cylinder is moving towards a disturbedbackground flow. The horizontal velocity plot (see left-side plot of Fig. 6(d)) indicatesthat the cylinder faces an incoming flow to the right. The effect of this disturbed flow isnot taken into account in our model. Therefore, as the cylinder travels further throughthe disturbed flow, it experiences a greater resistance than what is predicted by thedrag term in our model. Fig. 6(e) shows the flow field at the instant when the cylindercompleted its reverse swing (point (r) in Fig. 6(a) & (b)). By this time, the original wakein front of the cylinder appears to have dissipated.One of the approximations of our model is that the cylinder velocity v p nearly equalsthe relative velocity v rel between the cylinder and the flow. However, the PIV snapshotsreveal a strong mean flow after the first swing of the pendulum. This disturbed flow in thewake could result in v rel deviating significantly from v p , which is the likely reason for theoverprediction of the maxima of the second and subsequent oscillations. We look at themean fluid velocity U f at different angular positions in the cylinder wake. Fig. 7(a) showsdetails of the D × D window at a selected angular position θ sel used in the estimationof U f . The evolution of U f in time t − t pass is shown in Fig. 7(b). Here, t pass is the timewhen the cylinder passed θ sel . The wake decay shows some interesting characteristics.Initially, U f oscillates, which is a clear indication of vortex shedding in the cylinder wakeduring the initial moments after the cylinder has passed θ sel . We also see variations inthe peak for different angular positions θ sel . This part of the curve is highly unsteadyand difficult to model. However, later in time U f decays in a gradual and monotonic way.The inset to the figure shows the same plot on log-log scale. The long-time decay appearsto follow a power law decay.Several studies have addressed the decay of wakes for rising and falling particles (Wu &Faeth 1993; Bagchi & Balachandar 2004). These revealed many aspects of the temporaland spatial decays of wakes behind spheres in quiescent fluids, and in turbulent flows.2 Mathai, Loeffen, Chan, and Wildeman
Alméras et al. (2017) disentangled the wake decay behind rising bubbles into near wakeand far fields, with their characteristic decay rates. The wake decay observed in Fig. 7also represents the combined effect of the local decay of the wake and the advection ofthe wake with a velocity U f , i.e. dU f dt = ∂U f ∂t + U f L ∂U f ∂θ . These terms can be very difficultto disentangle for the wide range of m ∗ (or Re) in our experiments. Alternatively, onecan look at the wake velocity at locations just upstream of the cylinder during its returnswing. Fig. 8(a) shows the normalised wake velocity ˜ U f as a function of normalisedangular position θ ∗ for four different mass ratios. Here U f is normalised by the pendulumspeed at equilibrium position of each m ∗ case, thus representing a characteristic speedfor that m ∗ . For the heavier pendulums, v p is comparable to the gravitational velocityscale v g , but it deviates as m ∗ is reduced. Similarly, the maximum angle reached by thependulum at the end of its first swing θ max (point ( q ) in Fig. 6(a)(b)) is used in thenormalisation for θ . For clarity, we also show the arrow of increasing time, indicatingthat the return swing starts at θ ∗ → . With these normalisations the data show areasonable collapse. It is remarkable that the data collapse despite the large variation in m ∗ : [1.15, 4.98], which corresponds to a Reynolds number variation of almost a decade(Re ∈ [4000, 34000]). At the initial instant, i.e. when θ ∗ → , the cylinder sees a largeincoming flow velocity. As the pendulum swings back to θ ∗ ∼ , the wake velocity hasalmost completely decayed. Beyond this we observe a surprising increase in the incomingflow velocity. This increase arises from the cylinder wake during the first swing, wherethe pendulum had its highest swing velocity. Since the cylinder velocity at this positionwas large, the wake has not decayed completely.We obtain an approximate fit for the variation of ˜ U f vs θ ∗ (see Fig. 8(a)). With thewake flow modelled using ˜ U f ( θ ∗ ) shown above, we can include the effect of the incomingflow through a modified drag term: F Df ≈ ρ f A p C D ( v p − U f ) , where U f = ˜ U f v p .This can be implemented as a modified torque that replaces the original drag term ineq. 3.2. The modified equation of motion reads: m ∗ eff d θd ˜ t = − k sin θ − c (cid:12)(cid:12)(cid:12)(cid:12) dθ rel d ˜ t (cid:12)(cid:12)(cid:12)(cid:12) dθ rel d ˜ t − h | cos θ | sgn (cid:18) dθd ˜ t (cid:19) , (4.2)where the relative angular velocity dθ rel /d ˜ t = ( v p − U f ) / √ Lg . Note that Eq. 4.2 requiresno new input parameters, and works for the full range of Reynolds number in our ex-periments, i.e. Re ∈ [4000 , . The improvements in the predictions of the amplitudedecay are presented in Fig. 8(b), which shows the error in the predicted amplitude ofthe second peak of oscillation for all mass ratios ( m ∗ > ) studied. The predictions ofthe revised model are shown along with those of the original model. Accounting for thewake flow results in significant improvement in the second peak for all m ∗ cases. Themean error reduces from ∼
5. Conclusions and Discussion
In this work, we have studied the large amplitude oscillations of heavy and buoyantcylindrical pendulums in water. The oscillation frequency and peak deceleration of thependulum depend on the mass-density ratio m ∗ . Accounting for buoyancy and addedmass alone, as done in previous investigations (Neill et al. fit E rr o r ( % ) Original ModelRevised Model
532 4 (a) (b) increasing time m* = 4.98 m* = 3.74 m* = 2.20 m* =1.16 decreasing Figure 8. (a) Decay of the normalised wake velocity ˜ U f vs normalised angular position θ ∗ . Thewake velocity is normalised by the pendulum velocity at equilibrium position during the firstswing v p , and the angular position is normalised by the peak angular position at the end ofthe first swing θ max . Both v p and θ max can be found by solving eq. 3.2, and do not require anyinput from PIV. The normalisation leads to a reasonable collapse of the data despite the widevariation in Reynolds number from Re ∼ m ∗ = 1 . ) to Re ∼ m ∗ = 4 . ).Note that the arrow of time is pointing opposite to the arrow of θ , as indicated in the (a).(b) Error (%) in the peak amplitude after including the history force due to wake flow. Thecoloured bands show the error ranges for the original and revised model. The dashed lines showthe mean errors for the two models, which decreases from 14.8% (for the original model) to3.8% (for the revised model). reasonable agreement with the experimental observations for a wide range of mass ratios.The added mass coefficient in experiment is found to be m ∗ a = 0 . , i.e. lower that thetwo-dimensional (2-D) potential flow value for a cylinder. A reduced m ∗ a as comparedto the 2-D potential flow added mass can be expected, given that the flow is viscousand three-dimensional, and also due to the finite span of the pendulum. This result isconsistent with existing literature on oscillating bodies in viscous fluids (Dong 1978;Konstantinidis 2013; Koo & Kim 2015).The theoretical model presented here provides fair predictions of the oscillation fre-quency and peak decelerations for a wide range of m ∗ . However, the model overpredictsthe peak oscillation amplitudes for the second and subsequent oscillations. Introducinga Reynolds number dependent drag, along with vortex shedding forces, can providemarginal improvements to the model predictions. While modelling C D as a functionof Re improves the later oscillations, the vortex shedding term mainly influences the firstoscillation.Particle image velocimetry (PIV) measurements have revealed that the major factorcausing the deviations between experiment and the simplified model is the disturbed flowsurrounding the cylinder, which was not accounted for in the original model. We haveused insights from PIV measurements to obtain a simple model for the incoming wakeflow for a wide range of mass ratios. With this history force modelled, the revised modelpredictions are significantly improved; the mean error reducing from 14.8 to 3.8%. Thepredictions are also improved for other mass ratios for which PIV experiments were notperformed.Even with the wake history force included, the current model is still quite basic. Inreality, the dynamics is highly nonlinear, with changes in direction, curvilinear trajec-tories and wide variations in instantaneous Re. Under such conditions, exact analyt-ical expressions for the drag, added mass and history forces are not available. Fullyresolved direct numerical simulations (Immersed boundary (Mittal & Iaccarino 2005), or4 Mathai, Loeffen, Chan, and Wildeman
Physalis (Naso & Prosperetti 2010)) can provide better insights into the flow-inducedforces. On a different note, the present experiments have shown the importance of fluiddrag for a bluff body oscillating in a fluid. An interesting extension to the study wouldbe to use a more streamlined body, for which the nonlinear drag would be less dominant.In the absence of major flow separations, we can expect the frequency and dynamics toconform slightly better with a potential flow based added mass.We thank D. Lohse, C. Sun, A. Pandey and S. Maheshwari for fruitful discussions. Wethank G.-W. Bruggert, M. Bos and D. van Gils for technical support. V.M. acknowledgesfunding from the Netherlands Organisation for Scientific Research NWO, and EuropeanCooperation in Science and Technology, Eu-COST action MP1305.
REFERENCESAlméras, Elise, Mathai, Varghese, Lohse, Detlef & Sun, Chao
J. Fluid Mech. , 1091–1112.
Bagchi, Prosenjit & Balachandar, S
J. Fluid Mech. , 95–123.
Batchelor, G. K.
An introduction to fluid dynamics . Cambridge university press.
Bolster, D., Hershberger, R. E. & Donnelly, R. J.
Phys. Rev. E (4), 046317. Dong, RG
Tech. Rep. . CaliforniaUniv., Livermore (USA). Lawrence Livermore Lab.
Govardhan, R & Williamson, C. H. K.
J. Fluid Mech. , 85–130.
Govardhan, R. N. & Williamson, C. H. K.
J.Fluid Mech. , 11–47.
Hoerner, Sighard F
Fluid-dynamic drag: practical information on aerodynamic dragand hydrodynamic resistance . Hoerner Fluid Dynamics.
Huygens, C.
Christiaan Huygens’ the pendulum clock, or, Geometrical demonstrationsconcerning the motion of pendula as applied to clocks . Iowa State Pr.
Konstantinidis, Efstathios
Proc. R. Soc. A (2156), 20130135.
Koo, Weoncheol & Kim, Jun-Dong
International Journal of Naval Architecture and Ocean Engineering (1), 115–127. Lienhard, J. H.
Synopsis of lift, drag, and vortex frequency data for rigid circularcylinders , , vol. 300. Technical Extension Service, Washington State University.
Mathai, V., Calzavarini, E., Brons, J., Sun, C. & Lohse, D.
Phys. Rev. Lett. (2),024501.
Mathai, V., Prakash, V. N., Brons, J., Sun, C. & Lohse, D.
Phys. Rev. Lett. (12), 124501.
Mathai, Varghese, Zhu, Xiaojue, Sun, Chao & Lohse, Detlef
Phys. Rev. Lett. (5), 054501.
Mathai, Varghese, Zhu, Xiaojue, Sun, Chao & Lohse, Detlef
Nat. Commun. (1),1792. Mittal, R. & Iaccarino, G.
Annu. Rev. Fluid Mech. ,239–261. Naso, A. & Prosperetti, A.
New J. Phys. (3), 033040. Neill, D., Livelybrooks, D. & Donnelly, R. J.
Am. J. Phys. (3), 226–229. Obligado, Martin, Puy, Martin & Bourgoin, Mickaël
J. Fluid Mech. . Stokes, G. G.
On the effect of the internal friction of fluids on the motion of pendulums .Pitt Press Cambridge.
Tatsuno, M & Bearman, PW
J. Fluid Mech. ,157–182.
Tatsuno, M & Karasudani, T
Fluid Dyna. Res. (6), 313. Wu, J-S & Faeth, GM
AIAA journal31