J Zhao

J S Leontini

D Lo Jacono

J Sheridan

Chaotic vortex induced vibrations

This study investigates the nature of the dynamic response of an elastically mounted cylinder immersed in a free stream. A novel method is utilized, where the motion of the body during a free vibration experiment is accurately recorded, and then a second experiment is conducted where the cylinder is externally forced to follow this recorded trajectory. Generally, the flow response during both experiments is identical. However, particular regimes exist where the flow response is significantly different. This is taken as evidence of chaos in these regimes.

This letter presents the results of an investigation of the dynamic response of an elastically mounted cylinder immersed in a free stream. In particular, the focus is on the identification of low-dimensional chaos in the upper branch flow regime [START_REF] Khalak | Fluid forces and dynamics of a hydroelastic structure with very low mass and damping[END_REF][START_REF] Govardhan | Modes of vortex formation and frequency response of a freely vibrating cylinder[END_REF] of vortex-induced vibration (VIV). For the dynamically more complicated case of a flexible cylinder, chaos has been identified as a generic feature of the response, [START_REF] Modarres-Sadeghi | Chaotic response is a generic feature of vortex-induced vibrations of flexible risers[END_REF] Hover, Techet, and Triantafyllou [START_REF] Hover | Forces on oscillating uniform and tapered cylinders in crossflow[END_REF] and Hover [START_REF] Hover | Three-dimensionality of mode transition in vortex-induced vibrations of a circular cylinder[END_REF] show experimental evidence of low spanwise correlation and high levels of disorder in the upper branch, and three-dimensional simulations from Lucor, Foo, and Karniadakis [START_REF] Lucor | Vortex mode selection of a rigid cylinder subject to VIV at low mass-damping[END_REF] show a variation in the organization of the vortices along the span.

Morse and Williamson [START_REF] Morse | Fluid forcing, wake modes, and transitions for a cylinder undergoing controlled oscillations[END_REF] used data from experiments on an externally controlled sinusoidally oscillating cylinder to construct a low-order model of the interaction between the body and flow in the upper branch. It was concluded that the upper branch could consist of slow variation between the 2S and 2P o vortex shedding modes or intermittent switching between the 2P o and 2P modes. The mode naming follows the convention outlined in Williamson and Roshko. [START_REF] Williamson | Vortex formation in the wake of an oscillating cylinder[END_REF] The 2S mode consists of two single vortices shed per oscillation cycle, the 2P mode consists of two pairs of vortices shed per oscillation cycle, and the 2P o mode consists of two pairs of vortices per oscillation cycle where the second vortex of each pair is much weaker than the first. A map of the occurrence of these modes produced by Morse and Williamson, [START_REF] Morse | Fluid forcing, wake modes, and transitions for a cylinder undergoing controlled oscillations[END_REF] along with results from fully coupled experiments, is shown in Figure 1. The results of the current study augment and extend the conclusions of Morse and Williamson [START_REF] Morse | Fluid forcing, wake modes, and transitions for a cylinder undergoing controlled oscillations[END_REF] regarding the interaction of these modes in the upper branch.

A defining feature of low-dimensional chaos is the rapid divergence of trajectories in phase space that initially start close together. For an elastically mounted circular cylinder, it is often useful to consider the system as two coupled oscillators; one is the elastically mounted body, the other is the wake. The body motion is one phase space parameter, and a scalar representing the wake (often the lift force) can be used as another. If low-dimensional chaos is present, small perturbations to one parameter (the body motion) should result in large deviations of the other (the lift force).

In this system, the fluid and the structure are coupled by the fact that the fluid velocity must be equal to the velocity of the cylinder at its surface. The cylinder motion is driven by the fluid forcing on the body, which comes from integrating the pressure and shear stresses of the fluid on the cylinder surface.

The dependence of the lift force on the cylinder velocity ẏ, and therefore the cylinder displacement y, means that whether y follows a trajectory governed by the coupled system or is somehow FIG. 1. Plot of the amplitude response of a freely vibrating circular cylinder with m * = 2.66, with an overlay of the map of the wake regimes by Morse and Williamson. 7 The open circles represent the maximum amplitude at each normalized wavelength value of λ * = U * / f * , where f * is the normalized frequency close to the natural frequency of the system. The black solid symbols represent the locations where controlled motion experiments are conducted for comparison. The boundaries between the wake modes of Morse and Williamson 7 are indicated by the dashed lines. externally controlled, is unimportant. A particular trajectory y will result in a particular flow and lift force history. However if the flow is chaotic, even infinitesimally small perturbations to this trajectory y will result in large changes in the lift force history.

Here, this relationship is tested explicitly. Free vibration experiments have been conducted and the trajectory y accurately recorded at a sampling frequency of 1000 Hz. Following these, a second set of experiments has been conducted, where the motion recorded during the free vibration experiments has been fed back to a highly accurate motion control system, forcing the cylinder to follow the same trajectory as during the free vibration experiments. The results show that generally, the forces recorded during the free vibration and controlled motion experiments are very similar. However, for a case where the cylinder response is in the upper branch, a major deviation is recorded. This deviation is indicative of the presence of chaos, as any error in the y trajectory from the free vibration experiment in the controlled motion experiment results in diverging trajectories.

The experiments were conducted in a free-surface recirculating water channel with a test section of 4000 mm in length, 600 mm in width and, 800 mm in depth.

For the free vibration experiments, a test rigid circular cylinder was vertically mounted on a low-friction air bearing system which was placed on top of the water channel, allowing linear oscillatory movement of the cylinder transverse to the free-stream. Full details of this air bearing setup, control of end conditions, and the facility in general can be found in Nemes et al. [START_REF] Nemes | The interaction between flow-induced vibration mechanisms of a square cylinder with varying angles of attack[END_REF] and Zhao et al. [START_REF] Zhao | Fluid structure interaction of a square cylinder at different angles of attack[END_REF] The test cylinder used was a carbon fibre tube having an outside diameter of D = 25 mm and a wall thickness of 1.5 mm. The immersed length of the body was L = 620 mm, giving an aspect ratio of L/D = 24.8 and a displaced water mass of m d = πD 2 L/4 = 306.2 g. The total mass of the oscillating system was m = 813.2 g, resulting in a mass ratio of m * = m/m d = 2.66. Such a low mass ratio is required for the upper branch to occur. [START_REF] Williamson | Vortex-induced vibration[END_REF] Free decay tests were conducted individually in air and quiescent water to measure the natural frequencies of the system, f na = 0.874 Hz in air and f nw = 0.740 Hz in quiescent water. The structural damping ratio with consideration of the added mass effect was found to be

ζ = c/ (2  k(m + m A )) = 2.74 × 10 -3
, where c is the structural damping coefficient, resulting in a low massdamping ratio of m * ζ = 7.28 × 10 -3 , where m A = m d . The independent parameter of the reduced velocity was defined as U * = U/( f nw D), where U is the freestream velocity. The dynamic response of free VIV was investigated over a range of 2.5 < U * < 15 (1450 < Re < 7500).

Shown in Figure 1 is the amplitude response of the free vibration tests plotted in an amplitudewavelength plane overlaid on the map of wake regimes from Morse and Williamson [START_REF] Morse | Fluid forcing, wake modes, and transitions for a cylinder undergoing controlled oscillations[END_REF] for a cylinder undergoing sinusoidal forced oscillations. For comparison, the wavelength of free vibration is given by λ * = U * / f * , where f * is the normalized predominant oscillation frequency. The classic feature delimiting the upper and lower branches is the high and (relatively) low amplitudes of the oscillation in each branch, and this difference is clearly shown in Figure 1.

Four cases were selected for comparison with controlled motion experiments, and they are marked by the black solid symbols on Figure 1. These cases were selected to test the response in each of the major response regimes of VIV. One case is on the initial branch (U * = 4.09), one at the very beginning of the upper branch (U * = 4.93), one in the middle of the upper branch close to where the maximum amplitude is recorded (U * = 5.91), and one at the beginning of the lower branch (U * = 6.10).

In the controlled motion experiments, a real-time feedback closed-loop motion control system was employed to drive the test cylinder. This control system consisted of a hardware-in-theloop (HIL) control board (model: Q4, Quanser, Inc., Canada) embedded in a personal computer (PC) desktop, running on a combined software platform of MATLAB R2010a & Simulink and Real-Time Windows Target (MathWorks), and QUARC 2.1 (Quanser, Inc., Canada). The motion plant consisted of a digital servo drive (model: DPRANIE-015A400, Advanced Motion Controls, USA), a brushless DC rotary motor (model: SM231AL-NMSN, Parker Hannifin, USA) with high resolution incremental encoder feedback device (5000 pulses per revolution or 20 000 counts per revolution post quadrature) attached to a linear ball-screw actuator (12.7 mm in lead and 150 mm in stroke, model: ERS50-B02LA20-FSR150-A, Parker Hannifin, USA). The feedback control algorithm used was a Proportional-Integral-Derivative (PID) method, [START_REF] Callender | Time-lag in a control system[END_REF] which resulted in motion tracking errors within ±5 µm (0.0002D). Errors in the velocity and acceleration could be momentarily large (around 30% of the intended instantaneous value), particularly when the actuator changed direction. However, the error in the velocity was typically less than 5% of the intended instantaneous value, while the error in the acceleration was typically less than 10% of the intended instantaneous value.

Full details of the measurement techniques employed for the cylinder displacement, the transverse lift force, and the flow vorticity fields in the near wake of the body can be found in Zhao et al. [START_REF] Zhao | Fluid structure interaction of a square cylinder at different angles of attack[END_REF] Figure 2 shows time histories of the cylinder trajectory and the lift force for both the free vibration and controlled motion experiments for each of the four test cases. The most pertinent result shown in Figure 2 is the time history of the lift force for U * = 5.91. While all the other cases show a near-perfect match between the force history of the free vibration and controlled motion experiments, this one shows a significant deviation. The magnitude and frequency are similar, however there is a clear difference in the phase of the two signals. This mismatch leads to the conclusion that the flow is chaotic in the middle of the upper branch.

It should also be noted that the case at the beginning of the upper branch (U * = 4.93) shows no evidence of chaos. In fact, this case oscillates periodically, with a 2S wake mode. This periodic oscillation in the 2S mode in the upper branch is in contrast to the results presented in Govardhan and Williamson 2 that show a 2P o mode over the entire branch. Morse and Williamson 7 also concluded from controlled sinusoidal motion experiments that the 2S mode should only occur in competition with the 2P o mode in the upper branch. It is possible that this stable regime of the 2S mode occurs over a very small range of U * , or that the fully coupled oscillation departs (even very slightly) from purely sinusoidal, allowing the 2S mode to become stable.

To further investigate the apparent chaos in the middle of the upper branch, Figure 3 presents portions of the time histories of the total phase φ total and the vortex phase φ vortex for the free vibration and controlled motion experiments for the case in the middle of the upper branch (U * = 5.93). The total phase φ total is defined as the phase difference between the lift force and the cylinder This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 49.127.17.151 On: Tue, 23 Dec 2014 22:48:14 displacement, whereas the vortex phase φ vortex is defined as the phase between the vortex lift force and the cylinder displacement. The vortex force is defined using the method of Lighthill [START_REF] Lighthill | Fundamentals concerning wave loading on offshore structures[END_REF] and can be thought of as the force due to vortex shedding. The instantaneous values of φ total and φ vortex were found by finding the instantaneous value of phase of the force and displacement signals using a Hilbert transform and subtracting the phase of the force from the phase of the displacement.

Govardhan and Williamson, [START_REF] Govardhan | Modes of vortex formation and frequency response of a freely vibrating cylinder[END_REF] Carberry, Sheridan, and Rockwell, [START_REF] Carberry | Controlled oscillations of a cylinder: forces and wake modes[END_REF] and Morse and Williamson 7 have shown that the values of φ total and φ vortex are linked to the wake mode present. The 2S mode has φ total ≃ φ vortex ≃ 0 • , the 2P o mode has φ total ≃ 0 • and φ vortex ≃ 180 • , and the 2P mode has

φ total ≃ φ vortex ≃ 180 • .
Using these definitions, the phase history presented in Figure 3 shows that the free vibration case in the middle of the upper branch switches between the 2S and 2P o wake modes, as suggested by Morse and Williamson. [START_REF] Morse | Fluid forcing, wake modes, and transitions for a cylinder undergoing controlled oscillations[END_REF] However during the controlled motion experiment, the flow switches between the 2S, 2P o , and 2P wake modes. Phase averaged PIV images of these wake modes are also presented in Figure 3. The PIV images were obtained on a plane z/D = 4.5 from the cylinder free end. Previous results from Hover, Techet, and Triantafyllou 4 and Hover [START_REF] Hover | Three-dimensionality of mode transition in vortex-induced vibrations of a circular cylinder[END_REF] have shown that the wake mode can be different at the same instant at different points along the cylinder. It is possible that this occurred during the current experiments. However, sampling at z/D = 4.5 from the cylinder end (where a platform technique was used to mitigate the impact of the channel boundary layer) should provide a clear picture of the modes present and the amount of time the flow spends in each mode. The snapshots used to construct each image were selected based on the phase history, e.g., the images for the 2S mode were selected from instants where φ total ≃ φ vortex ≃ 0. The similarities between the 2S and 2P o modes from both the free vibration and controlled motion experiments are clear. However, the appearance of the 2P mode in the controlled motion experiments only is also This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 49.127.17.151 On: Tue, 23 Dec 2014 22:48:14 clear. These distinct modes and the apparent sudden switching between them illustrated by the time histories of φ total and φ vortex are evidence of chaos via mode competition. The competition between these modes for oscillations of the typical amplitude and frequency experienced in the upper branch is consistent with the data from the wake mode map presented by Morse and Williamson [START_REF] Morse | Fluid forcing, wake modes, and transitions for a cylinder undergoing controlled oscillations[END_REF] and shown in Figure 1. Further evidence of chaos via mode competition is found by looking at the time the flow spends in each mode. The presence of each mode at a given instant can be determined using the relationship between the instantaneous total and vortex phases outlined above. To further simplify, the phases This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 49.127.17.151 On: Tue, 23 Dec 2014 22:48:14 have been binned into two classes, i.e., a phase in the range -90 • < φ ≤ 90 • is treated as φ = 0 • and a phase in the range 90 • < φ ≤ 270 • is treated as φ = 180 • . Each instant can then be classified as belonging to the 2S, 2P o , or 2P mode, as indicated by the shaded regions in the example time series presented in Figure 3.

Figure 4 presents data from 300s-long segments of time histories for free vibration cases over a range of U * in the upper branch. Figure 4(a) presents a stacked histogram of the time spent in each mode as a fraction of the total length of the 300s segment as a function of U * . The plot shows that initially the flow is in the 2S mode practically all the time. Increasing U * sees a regime where the flow switches between the 2S and 2P o modes, with the proportion of time in the 2P o mode increasing with increasing U * . Around the middle of the upper branch, the flow remains in the 2P o mode practically all the time. Further increasing U * sees a regime where the flow switches between the 2P o and 2P modes, with the proportion of time in the 2P mode increasing with increasing U * . With the transition to the lower branch, the flow remains permanently in the 2P mode.

Figure 4(b) shows the mean, minimum, and maximum segment lengths as a proportion of the total time for each mode as a function of U * . As expected from Figure 4(a), initially the average segment length of the 2S mode is the whole measurement time. At high values of U * in the lower branch, the average segment length of the 2P mode is practically the whole measurement time. However, in the middle of the upper branch which is dominated by the 2P o mode, the 2P o segment lengths can still be quite variable, indicated by an average segment length of less than 1 and a significant variation between the minimum and maximum segment lengths. This indicates that the 2P o mode is never seen as a completely stable periodic state.

The overall picture presented by Figure 4 is that at the beginning of the upper branch, the dynamics are governed by a competition between the 2S and 2P o modes. At the end of the upper branch, the dynamics are governed by a competition between the 2P o and 2P modes. This may go some way to explaining why all three modes are found in the controlled motion experiment-small deviations from the free vibration trajectory cause the flow to visit some regions of parameter space where there is 2S-2P o competition and others where there is 2P o -2P competition.

To conclude, the argument for the appearance of chaos in the middle of the upper branch is summarized. Experiments were performed where the recreated trajectory matched the free vibration trajectory to within ±5 µm (0.0002D). Only for the case in the middle of the upper branch were differences in the force recorded, the large change due to a small perturbation indicative of chaos.

Further evidence for chaos due to mode competition is that in the middle of the upper branch the flow switches between at least two distinct wake modes. Mode competition and switching were implied by Morse and Williamson. [START_REF] Morse | Fluid forcing, wake modes, and transitions for a cylinder undergoing controlled oscillations[END_REF] Here, it has been shown there are two competitions in the upper branch (2S-2P o competition at the beginning of the branch and 2P o -2P competition at the end) and the statistics of this competition presented. Consistent with this finding, it has also been shown that a stable, periodic 2S wake mode exists at the very beginning of the upper branch, which has not been predicted by any previous analysis.

FIG. 2 .

 2 FIG. 2. Time histories in terms of periods of oscillation (τ) of cylinder displacement y/D for free vibration (black dashed line) and controlled motion (red solid line) and the lift coefficient C y for free vibration (black dashed line) and controlled motion (red solid line), for (a) U * = 4.09 (initial branch), (b) U * = 4.93 (beginning of upper branch), (c) U * = 5.91 (middle of upper branch) and (d) U * = 6.10 (lower branch). Only the case in the middle of the upper branch (figure (c)) shows a discrepancy between the lift force of the free vibration and controlled motion experiments.

FIG. 3 .

 3 FIG. 3. Upper: total phase φ total (black open circles) and vortex phase φ vortex (blue solid circles) for the free vibration and controlled motion case in the middle of the upper branch (U * = 5.91). Evidence of switching between two wake modes (2S and 2P o modes) is found during free vibration, while evidence of switching between three wake modes (2S, 2P o , and 2P modes) is found in the corresponding controlled motion experiment. The background of segments of each mode are shaded: 2S (dark grey), 2P o (light grey), 2P (white). Lower panes show images of the wake modes constructed via phase averaging. The vorticity contour levels shown are normalized by ω * z = ω z D/U.

FIG. 4 .

 4 FIG. 4. (a) A histogram of the portion of time spent in each of the modes 2S( ), 2P o ( ), and 2P ( ) in the upper branch during free vibration experiments. The flow moves from a 2S mode, to a mix of 2S and 2P o , to a mix of 2P o and 2P, to a 2P mode in the lower branch. (b) Mean (marked by points), minimum, and maximum segment lengths (marked by range bars) of the modes 2S( ), 2P o ( ), and 2P(•). For cases in the middle of the upper branch that spend almost all of their time in the 2P o mode, the flow does not settle to a stable periodic oscillation.

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 49.127.17.151 On: Tue, 23 Dec 2014 22:48:14

The authors would like to acknowledge the financial support of the Australian Research Council (ARC) through Grant No. DP110102141 under the Discovery program, and the Centre National de la Recherche Scientifique (CNRS) through Grant No. PICS161793 under the Projet International de Coopération Scientifique.


Justin S Leontini

David Lo

Mark C Thompson

Stability analysis of the elliptic cylinder wake

Keywords: instability, vortex streets, wakes/jets

This paper presents the results of numerical stability analysis of the wake of an elliptical cylinder. Aspect ratios where the ellipse is longer in the streamwise direction than in the transverse direction are considered. The focus is on the dependence on the aspect ratio of the ellipse of the various bifurcations to three-dimensional flow from the two-dimensional Kármán vortex street. It is shown that the three modes present in the wake of a circular cylinder (modes A, B and QP) are present in the ellipse wake, and that in general they are all stabilized by increasing the aspect ratio of the ellipse. Two new pertinent modes are found: one long-wavelength mode with similarities to mode A, and a second that is only unstable for aspect ratios greater than approximately 1.75, which has similar spatiotemporal symmetries to mode B but has a distinct spatial structure. Results from fully three-dimensional simulations are also presented confirming the existence and growth of these two new modes in the saturated wakes.

Introduction

The wake formed behind a cylinder of elliptical cross-section can be thought of as one example of the wider class of bluff-body wake flows. For many phenomena in this class of flows, the canonical case is that of a circular cylinder. It is well known that the circular cylinder wake progresses from creeping flow to a steady separated flow at around Re = 5 [START_REF] Taneda | Experimental investigation of the wakes behind cylinders and plates at low Reynolds numbers[END_REF]). The flow then undergoes a transition from a steady to periodic state at Re 47 [START_REF] Dušek | A numerical and theoretical study of the first Hopf bifurcation in a cylinder wake[END_REF][START_REF] Le Gal | Hysteresis in the forced Stuart-Landau equation: application to vortex shedding from an oscillating cylinder[END_REF]) apparently governed by a supercritical Hopf bifurcation [START_REF] Provansal | Bénard-von Kármán instability: transient and forced regimes[END_REF][START_REF] Thompson | The Stuart-Landau model applied to wake transition revisited[END_REF] where the Reynolds number Re = UD/ν, and U is the free stream velocity, D is the cylinder diameter and ν is the kinematic viscosity. Beyond this Re, the wake is characterized by the periodic two-dimensional vortex shedding of the Kármán vortex street. Excellent overviews of the regimes of flow past the circular cylinder are provided by [START_REF] Roshko | Perspectives on bluff body aerodynamics[END_REF] and [START_REF] Zdravkovich | Flow Around Circular Cylinders: Volume 1: Fundamentals[END_REF].

Further increases in Re see the Kármán wake subsequently become unstable to three-dimensional perturbations. At approximately Re = 190, the advecting vortex street becomes unstable to mode A [START_REF] Williamson | The existence of two stages in the transition to three-dimensionality of a cylinder wake[END_REF][START_REF] Barkley | Three-dimensional Floquet stability analysis of the wake of a circular cylinder[END_REF], characterized by a waviness of the vortex cores and speculated to be caused by a cooperative elliptic instability of these cores [START_REF] Landman | The three-dimensional instability of strained vortices in a viscous fluid[END_REF][START_REF] Leweke | Three-dimensional instabilities in wake transition[END_REF][START_REF] Thompson | The physical mechanism of transition in bluff body wakes[END_REF][START_REF] Ryan | Three-dimensional transition in the wake of bluff elongated cylinders[END_REF][START_REF] Leontini | Three-dimensional transition in the wake of a transversely oscillating cylinder[END_REF]. (See also [START_REF] Giannetti | Structural sensitivity of the secondary instability in the wake of a circular cylinder[END_REF] for an analysis in terms of structural stability.) Further increases in Re see the onset of mode B, characterized by fine-scale streamwise vortices growing in the braid shear layers between the wake vortices. Stability analysis predicts mode B to begin at Re = 260 [START_REF] Barkley | Three-dimensional Floquet stability analysis of the wake of a circular cylinder[END_REF], however in experiments [START_REF] Williamson | Three-dimensional wake transition[END_REF] and direct numerical simulations [START_REF] Thompson | Three-dimensional instabilities in the wake of a circular cylinder[END_REF] mode B is observed for Re 230. The presence of mode A appears to have a destabilizing influence on mode B [START_REF] Barkley | Bifurcation theory for three-dimensional flow in the wake of a circular cylinder[END_REF][START_REF] Sheard | A coupled Landau model describing the Strouhal-Reynolds number profile of the three-dimensional wake of a circular cylinder[END_REF]. The three-dimensional numerical study of [START_REF] Henderson | Nonlinear dynamics and pattern formation in turbulent wake transition[END_REF] provides a thorough overview of both modes A and B. A third mode is predicted by stability analysis to occur around Re = 377 [START_REF] Blackburn | Symmetry breaking of two-dimensional time-periodic wakes[END_REF], known as mode QP (for quasi-periodic), but clear evidence of it has not been found experimentally, most likely due to the presence and dominance of mode B.

The symmetry group argument presented by [START_REF] Blackburn | Symmetry breaking of two-dimensional time-periodic wakes[END_REF] states that there are only three generic three-dimensional spatiotemporal wake states that can develop from the two-dimensional Kármán wake. The symmetries of modes A, B and QP each represent one of these states. It might therefore be expected that the transition to three-dimensional flow in the wake of a cylinder may be taken as representative of that for any body that produces a Kármán wake. However, a body of research has been built up regarding three-dimensional transition for flows past other bodies and modifications to the basic flow past a cylinder. Square cross-sections [START_REF] Robichaux | Three-dimensional Floquet instability of the wake of a square cylinder[END_REF][START_REF] Sheard | Cylinders with square cross-section: wake instabilities with incidence angle variation[END_REF], elliptical leading-edge plates [START_REF] Ryan | Three-dimensional transition in the wake of bluff elongated cylinders[END_REF], transverse-oscillating cylinders [START_REF] Leontini | Three-dimensional transition in the wake of a transversely oscillating cylinder[END_REF]) and rotational oscillating cylinders [START_REF] Lo | Modification of three-dimensional transition in the wake of a rotationally oscillating cylinder[END_REF] are some examples of such flows where the symmetry of the resulting base flows are the same as the Kármán vortex street flow past a cylinder, that is spatiotemporally symmetric such that u(x, y, t) = u(x, -y, t + T/2), where u is the velocity field and T is the period of the vortex shedding in the wake. Many of these studies show that the order of inception of modes A, B and QP can change. They also show that other modes can be found, albeit with one of the same three spatiotemporal symmetries.

Beyond some critical value of Re, the wake of an ellipse is the Kármán vortex street. The critical value of Re for the onset of vortex shedding is approximately proportional to the aspect ratio [START_REF] Jackson | A finite-element study of the onset of vortex shedding in flow past variously shaped bodies[END_REF]. The frequency of vortex shedding, represented by the Strouhal number St = f D/U, where f is the frequency of vortex shedding and D is the length of the body across the flow, is also approximately proportional to the aspect ratio [START_REF] Kuo | The formation of elliptical wakes[END_REF][START_REF] Modi | Near-wakes of elliptic cylinders in subcritical flow[END_REF]. This is consistent with the theory of a universal Strouhal curve from [START_REF] Roshko | On the wake and drag of bluff bodies[END_REF], that states as the separating shear layers on each side of the body get closer to each other, St should increase. Given this universality, it should be expected that, at least for small deviations from the circular cylinder geometry, modes A, B and QP will still exist, albeit perhaps somewhat modified. The results of this paper confirm that this is the case, but for larger departures from a circular cylinder, new modes also arise.

There appears to be little published information regarding the development of three-dimensional flow in the wakes of elliptical cylinders. [START_REF] Mittal | Effect of three-dimensionality on the lift and drag of nominally two-dimensional cylinders[END_REF] reported on two-and three-dimensional simulations of the flow past cylinders and ellipses with a focus on the Re = 525 case for aspect ratio Γ = 2 (Γ = major/minor axis ratio). The two-dimensional simulations resulted in a peak lift almost twice as large as the three-dimensional simulations, indicating the flow is truly three-dimensional for these conditions. More recently [START_REF] Thompson | Low-Reynolds-number wakes of elliptical cylinders: from the circular cylinder to the normal flat plate[END_REF] examined the three-dimensional transitions for elliptical cylinders with aspect ratios less than one, finding transitions similar to those that occur for a normal flat plate at low aspect ratios, and convergence towards the circular cylinder transition sequence as the aspect ratio became close to unity.

In this paper, the variations of the various three-dimensional modes as a function of aspect ratio Γ are studied. Stability analysis of the two-dimensional periodic wakes of ellipses with 1 Γ 2.4 is undertaken. The case Γ = 1 corresponds to the canonical circular cylinder.

Section 2 presents an overview of the techniques employed both for the base flow simulations and the stability analysis. Section 3.1 first presents results analysing the onset of vortex shedding and the Kármán wake. Then, § § 3.2 and 3.3 present results of Floquet stability analysis of these Kármán wakes. These results show the modification of modes A, B and QP as Γ is changed, as well as the appearance of two new modes, A and B. Sections 3.5 and 3.6 present some details of these new modes, and § 3.7 presents limited results of fully three-dimensional simulations confirming the existence and growth of these new modes. Finally, some concluding remarks are presented in § 4.

Methodology

Governing equations and Floquet stability analysis

The two-dimensional base flows for this paper were found by solving the nondimensionalized incompressible Navier-Stokes equations,

∂u ∂τ = -(u • ∇)u -∇P + 1 Re ∇ 2 u, ∇ • u = 0, (2.1)
where u is the velocity field non-dimensionalized by the free stream velocity U, τ = tU/D is the non-dimensional time where t is time and D is the length of the body across the flow, P = p/(0.5ρU 2 ) is the non-dimensional pressure where p is the pressure and ρ is the fluid density and Reynolds number Re = UD/ν, where ν is the kinematic viscosity.

Apart from the Reynolds number Re, the only other variable controlling the flow state is the aspect ratio of the ellipse Γ = L/D, where L is the length of the body in the streamwise direction.

The stability of two-dimensional periodic solutions of these equations is then assessed using Floquet stability analysis. Essentially, this consists of solving for the growth of a three-dimensional periodic perturbation of a given spanwise wavelength λ from one period to the next. This growth is represented by the Floquet multiplier, µ f , which can be a complex quantity. The magnitude of µ f dictates the stability of the flow: if |µ f | > 1, the flow is unstable, indicating that perturbations will grow. Critical conditions occur at neutral stability, where

|µ f | = 1.
The growth of perturbations is governed by the linearized Navier-Stokes equations,

∂u ∂τ = -(u • ∇)u -(u • ∇)u -∇P + 1 Re ∇ 2 u , ∇ • u = 0,    (2.2)
where u is the three-dimensional perturbation velocity field, and P is the perturbation pressure field.

In this general form, u is a three-dimensional vector field and P is a threedimensional scalar field. These quantities can be further subdivided by decomposing into Fourier modes in the z-direction. Following [START_REF] Barkley | Three-dimensional Floquet stability analysis of the wake of a circular cylinder[END_REF], if the perturbations are written as u (x, y, z, τ ) = û cos(2π/λ)z, v cos(2π/λ)z, ŵ sin(2π/λ)z , P (x, y, z, τ ) = P cos(2π/λ)z, (2.3) then the perturbation equations shown in (2.2) can be written as

∂ û ∂τ = - û ∂U ∂x + v ∂U ∂y + U ∂ û ∂x + V ∂ û ∂y - ∂ P ∂x + 1 Re ∂ 2 û ∂x 2 + ∂ 2 û ∂y 2 -(2π/λ) 2 û ∂ v ∂τ = - û ∂V ∂x + v ∂V ∂y + U ∂ v ∂x + V ∂ v ∂y - ∂ P ∂y + 1 Re ∂ 2 v ∂x 2 + ∂ 2 v ∂y 2 -(2π/λ) 2 v ∂ ŵ ∂τ = -U ∂ ŵ ∂x + V ∂ v ∂y -(2π/λ) P + 1 Re ∂ 2 ŵ ∂x 2 + ∂ 2 ŵ ∂y 2 -(2π/λ) 2 ŵ ∂ û ∂x + ∂ v ∂y + (2π/λ) ŵ = 0.                              (2.4)
These equations are explicitly a function of the wavelength λ, meaning the twodimensional perturbation field û is a function of the Reynolds number Re and the wavelength of the perturbation in the z-direction λ.

Floquet stability analysis works by finding the eigenvalues µ f of some linear operator L that maps the perturbation field from one period to the next, i.e. û(x, y, λ, Re, t + T) = L( û(x, y, λ, Re, t)).

(2.5)

If L has eigenvalues where |µ f | > 1, then the perturbation will grow over one period of oscillation and the flow is unstable. Here, L need never be explicitly formed: its action is equivalent to integrating equations (2.4) forward in time by one period.

Values of the Floquet multiplier µ f and the corresponding eigenvectors (or Floquet modes, or simply modes) are a function of the flow variables Re and Γ and the wavelength λ. Therefore, the analysis task becomes a matter of calculating the leading modes (those with the largest corresponding µ f ) for each combination of Re, Γ and λ. This process is well described in [START_REF] Barkley | Three-dimensional Floquet stability analysis of the wake of a circular cylinder[END_REF], Tuckerman & Barkley (2000), [START_REF] Ryan | Three-dimensional transition in the wake of bluff elongated cylinders[END_REF] and Rao et al. (2013b).

Numerical method

Both the base flow and perturbation equations (2.1) and (2.2) were solved using a spectral-element code [START_REF] Karniadakis | Spectral/HP Methods for Computational Fluid Dynamics[END_REF]. The spatial domain is split into quadrilateral finite elements that employ high-order Lagrange polynomials as shape and weighting functions. These polynomials are associated with Gauss-Lobatto-Legendre quadrature points which leads to the efficient calculation of integrals required for solving the equations in the weak form. Time integration was executed using a three-step time-splitting scheme [START_REF] Karniadakis | High-order splitting methods of the incompressible Navier-Stokes equations[END_REF][START_REF] Thompson | Hydrodynamics of a particle impact on a wall[END_REF], employing a third-order Adams-Bashforth scheme for the convective substep, a centred scheme for the pressure correction and an implicit (theta-modified) Crank-Nicholson scheme [START_REF] Canuto | Spectral Methods in Fluid Dynamics[END_REF] for the diffusion term.

In the current implementation, the elements making up the mesh can have straight or constant-radius curved surfaces. The elements at the body boundary were therefore curved with a curvature calculated to minimize the difference between the constantradius curve and the true elliptical surface. The error was generally imperceptible (the root mean square (RMS) of the error in the y-position was smaller than 10 -4 D for all meshes), and the favourable comparison with previously published results of [START_REF] Jackson | A finite-element study of the onset of vortex shedding in flow past variously shaped bodies[END_REF] presented below in § 3.1 confirms this method presents an accurate representation of an ellipse.

A domain extending 20D upstream, 46D downstream and 30D to either side of the centre of the body was employed, giving a blockage ratio of less than 2 %. Here D refers to the length of the body across the flow. For the velocity field of the base flow, a free stream Dirichlet condition was set at the upstream and transverse boundaries and a Neumann condition set at the downstream boundary. A no-slip Dirichlet condition was set at the body surface. A Neumann condition for the gradient normal to the boundary calculated from the Navier-Stokes equations was set for the pressure at all boundaries [START_REF] Karniadakis | High-order splitting methods of the incompressible Navier-Stokes equations[END_REF] except the outflow boundary where the pressure was set to zero. The perturbation velocity was set to zero at all boundaries except the outflow boundary where a zero normal gradient condition was imposed on the perturbation velocity components.

The analysis was conducted by first running the simulation of the base flow only, until the base flow had saturated to a periodic state (typically after around 40 vortex shedding cycles). Then, the perturbation velocity field was initialized to random noise of amplitude 10 -3 , and the simulation continued, integrating both the base flow and the perturbation field forward in time. The perturbation field was stored once per period, and these saved fields used to find the leading eigenvectors and corresponding eigenvalues.

This was done using a Krylov subspace method, and an Arnoldi decomposition employed to resolve the full complex leading multipliers. The saved perturbation fields (typically five to eight successive fields) formed the Krylov subspace. The Arnoldi decomposition was then used to project these vectors into an orthonormal space. The resulting orthonormal vectors approximate the eigenvectors (or Floquet modes) of the system, and the coefficient matrix associated with the orthonormalization has eigenvalues that approximate the leading eigenvalues associated with the Floquet modes. The implementation here is given by [START_REF] Mamun | Asymmetry and Hopf-bifurcation in spherical Couette flow[END_REF] and Tuckerman & Barkley (2000), and has been successfully applied in previous studies of bluff-body wakes such as stenotic pipes [START_REF] Griffith | Pulsatile flow in stenotic geometries: flow behaviour and stability[END_REF], rotationally oscillating cylinders [START_REF] Lo | Modification of three-dimensional transition in the wake of a rotationally oscillating cylinder[END_REF], sphere-wall impact [START_REF] Thompson | Sphere-wall collisions: vortex dynamics and stability[END_REF] and constantly rotating cylinders (Rao et al. 2013a,b). [START_REF] Jackson | A finite-element study of the onset of vortex shedding in flow past variously shaped bodies[END_REF].

A resolution test on this mesh was conducted by varying the order of the shape functions p. The maximum lift force on the ellipse of Γ = 2 at Re = 400 was found to vary by less than 0.5 % when varying p from 5 to 8. All of the subsequent simulations were performed using p = 8.

Results

Characteristics of the base flows

The primary focus of this paper is the loss of stability of the Kármán wake to three-dimensional perturbations at a given Re. This assumes that the flow has already undergone transition from a steady flow to the periodic Kármán wake.

Data for the onset of this periodic shedding are provided in table 1. The value of Re, as well as the value of St at the onset of vortex shedding are provided, and compared to the values presented in [START_REF] Jackson | A finite-element study of the onset of vortex shedding in flow past variously shaped bodies[END_REF]. To find the value of Re at the onset of vortex shedding, a series of simulations was run at increments of Re of 1, and it was observed whether small perturbations grew or decayed. Growth was measured by monitoring the lift force coefficient C l = F L /(1/2)ρU 2 D on the body, where F L is the lift force per unit length, and ρ is the fluid density. It was then assumed that the transition to periodic flow could be adequately modelled as a supercritical Hopf bifurcation [START_REF] Provansal | Bénard-von Kármán instability: transient and forced regimes[END_REF][START_REF] Dušek | A numerical and theoretical study of the first Hopf bifurcation in a cylinder wake[END_REF][START_REF] Thompson | The Stuart-Landau model applied to wake transition revisited[END_REF], governed by the normal form equation

dA dτ = σ A + l c |A| 2 A, (3.1)
where A is the amplitude of the oscillation (here taken as the amplitude of the lift force coefficient), τ is the non-dimensional time tU/D where t is time, σ is the initial linear growth rate and l c is the Landau coefficient, which governs the saturation of the oscillation to a stable periodic state. This amplitude equation is often referred to as a Stuart-Landau model [START_REF] Stuart | On the non-linear mechanics of hydrodynamic instability[END_REF]. By dividing by A, (3.1) can be written as

d log(A) dτ = σ + l c |A| 2 . (3.2)
The values of σ and l c for a given simulation were then determined from a leastsquares fit to the amplitude data of the form shown in (3.2). The values of σ as a function of Re were then fitted and the fit extrapolated to find the value of Re where σ = 0.

To find the frequency at the onset of vortex shedding, the frequency of the saturated state of the same set of simulations was measured, giving St as a function of Re. A curve was then fitted to these data, and the curve extrapolated to find the frequency at the value of Re corresponding to σ = 0.

The comparison between the current results for the value of Re at the onset of shedding and those of [START_REF] Jackson | A finite-element study of the onset of vortex shedding in flow past variously shaped bodies[END_REF] is good, differing by a maximum of only 2.3 %. Whilst there is a larger discrepancy in the predicted frequencies, with differences of up to 15 %, it should be noted that the simulations of [START_REF] Jackson | A finite-element study of the onset of vortex shedding in flow past variously shaped bodies[END_REF] were performed on a small domain, and so are affected by substantial blockage, artificially raising St. Support for this comes from more recent studies (e.g. [START_REF] Thompson | The Stuart-Landau model applied to wake transition revisited[END_REF][START_REF] Sipp | Global stability of base and mean flows: a general approach and its applications to cylinder and open cavity flows[END_REF] based on global stability analyses of the steady circular cylinder wake (Γ = 1). These studies report values of Re c = 46.6 and 46.4 at the onset of vortex shedding, and global mode frequencies of St = 0.11777 and 0.1182, respectively. These values are very close to those found here, again with the small differences consistent with different domain blockage ratios, so providing some confidence in the results obtained for this study.

3.2. Variation of the growth rate of modes as a function of Γ , Re = 300 In this section, results for a fixed Re = 300 are presented as a function of Γ . The value Re = 300 is chosen as both of modes A and B in the wake of a circular cylinder are known to be unstable at this value. Figure 1 presents |µ f | as a function of spanwise wavelength λ for six aspect ratios ranging from Γ = 1 (a circular cylinder) to Γ = 2.

For wavelengths λ 5, the Γ = 1 case shows the familiar picture for a circular cylinder; the flow is unstable to mode A over a band of wavelengths focused at λ 3.5D and unstable to mode B for wavelengths of λ 0.8D. Mode QP is found at λ 1.8D, but it is not unstable at this Re, achieving a maximum |µ f | 0.7.

For longer wavelengths, 5 < λ 10, two other modes can be identified. The first of these modes, here named A, appears to be a continuation of the mode A branch. It has the same spatiotemporal symmetries as mode A, and a purely real Floquet multiplier, yet has a slightly different spatial structure, as will be discussed in more detail later. The second of these modes is a quasi-periodic mode with a complex Floquet multiplier, and is here named mode QPL (for long-wavelength quasi-periodic). It has a maximum growth rate for λ 9D.

With moderate increases in Γ , this picture remains more or less the same. These same five modes are the only relevant modes up to Γ = 1.4. Increases in Γ over this range see |µ f | for mode A and mode B decrease. However, |µ f | for mode QP increases, indicating an increase in growth rate. Despite this, |µ f | for mode QP remains below |µ f | = 1, indicating it does not become unstable at this Re. The maximum growth rate for mode A changes little with Γ , however the range of λ over which it is unstable increases, extending to higher λ with increasing Γ . The growth of mode QPL remains largely untouched, with a maximum |µ f | 0.8 occurring at λ 9D.

Further increases in Γ up to Γ = 2 show similar behaviour for modes A and B. Increases in Γ see a reduction in |µ f |. Mode B is stabilized (µ f < 1) at Γ 1.4; mode A is stabilized at Γ 1.8. The |µ f | for mode QP also reduces, and this mode is not resolved for Γ > 1.8. A clear distinction between modes A and A is seen for Γ 1.8; mode A shows a clear local maximum in growth rate at λ 8 for Γ = 1.8, further delineating it from mode A. Mode QPL is again practically unchanged. An important feature is that at Γ = 1.6, another mode with wavelengths around λ = 2.8, a purely real µ f and spatiotemporal symmetry the same as for mode B is found. This mode, here named mode B, sees |µ f | increase with increasing Γ . Mode B is shown to be unstable for Γ 1.8.

Figure 2 shows images of each of the modes, all at Re = 300, in order of increasing characteristic wavelength λ. All of the modes except mode B are shown at Γ = 1.8. Mode B is shown at Γ = 1.4 as µ f becomes so low at Γ = 1.8 that the mode is no longer identifiable. The images show the spanwise perturbation vorticity as filled colour contours, overlaid with black contour lines showing the positions of the base flow vortex cores.

The structures of modes A, B and QP are the same as the structures found in the wake of a circular cylinder [START_REF] Barkley | Three-dimensional Floquet stability analysis of the wake of a circular cylinder[END_REF][START_REF] Blackburn | Symmetry breaking of two-dimensional time-periodic wakes[END_REF]. It should be noted that mode QP has a complex µ f which introduces a second frequency, and as such the image shown in figure 2 is at an arbitrary phase.

Mode B is shown to have a distinct structure. While the symmetry of this mode is the same as mode B, the mode shares many features with mode A, particularly the fact that it seems mainly focused in the wake vortex cores (although there is Mode A has a similar structure to mode A in the two vortices closest to the body, but is distinct from the third vortex (circled with a dashed line) and further downstream. Note that modes QP and QPL are shown at an arbitrary phase as they have a complex µ f . also significant amplitude in the braid shear layers that connect the wake vortices). However, mode B has a finer structure on these vortex cores than mode A. An example of this is found by focusing on the negative vortex formed from the top of the body that is about to be shed into the wake immediately behind the body. Traversing across the centre of this approximately elliptical vortex along the short axis sees the sign of the perturbation vorticity change twice in mode B (black-whiteblack); for mode A, it changes sign once (white-black). This is, of course, consistent with the shorter preferred wavelength of mode B relative to mode A.

A mode found in the wake of elliptic leading-edge plates [START_REF] Ryan | Three-dimensional transition in the wake of bluff elongated cylinders[END_REF]), named mode B , shares many of the characteristics of mode B. It too has the same spatiotemporal symmetries as mode B, with a critical wavelength in the vicinity of λ = 2. The spatial structure of the two modes near the body is also quite similar in topology. These facts make it likely that mode B found here and mode B found in the wake of elliptic leading edge plates are due to the same physical mechanisms.

Exactly what these mechanisms are is not completely clear, and given the complicated structure of this mode it is possible that there are a number of cooperative mechanisms that contribute to its growth.

Mode A has a structure very similar to mode A near the body. In the images shown in figure 2, the structure of modes A and A are almost the same in the two vortices closest to the body. However, there is a distinct difference from the third vortex and further downstream. For instance, figure 2 shows mode A has one white lobe and one black lobe in the lower portion of the third vortex, and a white lobe in the upper portion of this vortex. Mode A instead has a white lobe and a black lobe, both of which extend over the entire vortex (the vortices in question are circled with a dashed line in figure 2). Further downstream in the wake (but not shown here), mode A displays a spanwise perturbation vorticity distribution strongly characteristic of elliptic instability. i.e. positive and negative lobes aligned at approximately 45 • to the main axis of the base flow vortices. This is not as apparent in the downstream wake for the shorter wavelength mode A. The structure of mode QPL is not dissimilar to mode QP. However, comparison of this mode's structure is more difficult as it is dependent on the phase at which images are obtained. Over the parameter space tested, mode QPL was never seen to become unstable, and therefore further investigation of its structure does not seem warranted.

Variation of critical Re and wavelength for all modes as a function of Γ

Expanding on the results at a single Re = 300 presented in § 3.2, this section presents the critical values of Reynolds number Re c and wavelength λ c for all four modes as a function of Γ . Re c is the value of the Reynolds number at which a given mode is predicted to first become unstable (where |µ f | first exceeds |µ f | = 1 at some value of wavelength); λ c is the wavelength at which this occurs. For a given Γ , calculations at a series of λ in steps of 0.2 for a fixed Re were run to find the maximum |µ f | for a given mode, similar to the data presented in figure 1. The value of Re was then incremented in steps of 20, finding a maximum |µ f | and associated value of λ for each Re. This process was repeated until the maximum |µ f | exceeded |µ f | = 1. The values of Re c and λ c were then found by quadratic interpolation of the surrounding values. Figures 3 and4 present the results of these calculations. This same process was then performed for each aspect ratio. Figure 3 shows a number of features that warrant discussion. First, Re c is almost a linear function of Γ for both modes A and B over the entire range tested. It also shows that, regardless of Γ , the onset of mode A always precedes mode B.

This same figure shows that the relationship is more complicated for mode QP, as Re c first decreases with Γ (for Γ < 1.5) and then increases with Γ . Regardless, Re c for mode QP is always well beyond Re c for mode A, and so it is difficult to predict the influence this mode will have in a physical flow.

For Γ > 1.2, mode A is found to be unstable. Below this value of Γ , mode A coalesces with mode A. Figure 3 shows that Re c is almost equal for modes A and A at Γ = 1.2, and figure 4 shows that the critical wavelengths of the two modes are also the same at this Γ , both providing evidence of the mode coalescence at this Γ .

For Γ > 1.4, mode A is clearly unstable at lower Re than mode A. In fact, with respect to increasing Re, mode A is the first mode to become unstable until Γ = 1.9, where mode B becomes the first to become unstable.

For Γ > 1.75, mode B is found to be unstable. Unlike the other modes, an upper limit is found, such that mode B is only unstable over a finite range of Re. Of course, this is based on a linear analysis that fails to hold once the wake has undergone three-dimensional transition, changing the base flow state. In any case, mode B is the leading instability mode, i.e. it is the first to become unstable with increasing Re, for Γ > 1.9.

The prevalence of these two new modes, A and B, means that for moderate aspect ratios Γ > 1.5, the three-dimensional structure of the flow should be controlled by these modes, at least at Re just above the first onset of three-dimensionality. This is further investigated in § 3.7.

The values of λ c presented in figure 4 show that the change in Γ has the greatest impact on the wavelengths of modes A and B. Figure 4(a) shows that the λ c for mode A first increases mildly, then decreases with increasing Γ , falling to approximately λ c = 3.5 at Γ = 2, from λ c = 4 for the circular cylinder (Γ = 1). Once it is resolved, λ c for mode B decreases with increasing Γ . The value of λ c for mode B is essentially unchanged, and that for mode QP only varies slightly.

These changes naturally reflect the changes in the base flow. With increasing Γ , the body becomes more streamlined, allowing the flow to remain attached for longer, reducing the wake width and the lateral spacing between vortices. The size and shape of these wake vortices will also be changed, becoming smaller. This is also caused by the reduced acceleration of the boundary-layer fluid as it advects over the ellipse prior to separation, leading to a lower vorticity flux into the wake. As outlined in § 1 (see also [START_REF] Leweke | Three-dimensional instabilities in wake transition[END_REF][START_REF] Thompson | The physical mechanism of transition in bluff body wakes[END_REF], the mode A spanwise perturbation vorticity distribution shows strong characteristics associated with elliptic instability of shed near-wake vortices, and so it seems likely that its growth and preferred wavelength will be influenced by changes in the base flow (vorticity) configuration. The images presented in figure 2 show that mode B also retains a high perturbation amplitude in the advecting vortex cores, and so it seems logical to assume that it too will be dependent on their configuration. On the other hand, mode B scales on the braid shear layers between vortices [START_REF] Leweke | Three-dimensional instabilities in wake transition[END_REF], and this scale apparently is not as susceptible to changes in the body aspect ratio as the wake vortex size and strength.

Figure 4(b) shows that λ c for mode A increases almost linearly with increasing Γ . It also shows that the λ c for mode A can become very long, approaching λ c = 28D at Γ = 2.4. However, the wavelength of highest growth for mode A is a strong function of Re. Therefore, as Re is increased beyond Re c , the wavelength of mode A typically becomes much shorter than λ c . This fact is further investigated and highlighted in § 3.6.

Behaviour of mode A in terms of elliptic instability theory

Leweke & Williamson (1998) employed the elliptic instability theory developed by [START_REF] Landman | The three-dimensional instability of strained vortices in a viscous fluid[END_REF] to the vortices in the wake of a circular cylinder. Even though the original theory was based on the instability of a single vortex with constant vorticity and a constant strain rate causing elliptical streamlines, [START_REF] Leweke | Three-dimensional instabilities in wake transition[END_REF] showed that the behaviour of the instability leading to mode A was at least consistent with this theory. Here, we take the same theory, and show that it predicts both the decrease in λ c with increasing Γ , and the increase in Re c .

The theory developed by [START_REF] Landman | The three-dimensional instability of strained vortices in a viscous fluid[END_REF] predicts the growth rate of perturbations to an elliptic vortex, given the eccentricity and size of the vortex. The eccentricity β = /(2|ω|), where is the strain rate and ω is the vorticity. For a twodimensional flow, the local strain rate in Cartesian coordinates is given by

= ∂u ∂x ∂v ∂y - 1 4 ∂u ∂y + ∂v ∂x 2 , (3.3)
which is the square root of the determinant of the symmetric part of the velocity gradient tensor. For a purely elliptic vortex, both ω and are constant, however both can vary in space in a more complicated flow such as a bluff-body wake, and hence the local eccentricity β can also vary. Figure 5 shows images of the wake for a series of Γ , at Re c for each case, where contours of vorticity are marked with solid lines, and values of the eccentricity are marked with greyscale contours. All of the images are shown at the instant where the lift force on the body is a maximum.

Figure 5 shows that regardless of Γ (and the variation of Re for each case), the eccentricity β is basically constant with a value of around β = 0.6 at the centre of the shed vortices. The eccentricity is not constant over the entire vortex but there is a region at the centre of each vortex where β is not strongly varying. FIGURE 5. The wake at the instant of maximum lift force as a function of aspect ratio, at Re c for each case. Contours of non-dimensional vorticity at levels ±1 are marked with solid lines. Greyscale colour contours represent the local eccentricity β at levels from 0 (black) to 1 (white). The white and black arrows on the example at Γ = 1 show the measurement of the lengths of the minor and major axes l minor and l major of the wake vortices, respectively. The point of intersection of these lines indicates the point where β of each vortex was measured. The vortices are shown to decrease in size with increasing Γ , while the eccentricity at the centre of the vortices remains around a value of β 0.6.

1.7 1.8 1.9 2.0 2.1 1.0 1.2 1.4 1.6 1.8 2.0 2.2 2.4
FIGURE 6. Variation of the size of the wake vortices with Γ . Here, L = l 2 major + l 2 minor where l major and l minor are the lengths along the major and minor axes of the vortex indicated in figure 5. Points are measurements, the dashed line is a linear fit.

The figure also shows that as Γ increases, the size of the vortices in the wake decreases. To quantify this decrease, the vortex furthest from the body in each of the images in Figure 5 was measured along its major and minor axes (where the boundary of the vortex was defined, reasonably arbitrarily, as the contour of |ω| = 1), and a length scale L = l 2 major + l 2 minor defined, where l major and l minor are the lengths along the major and minor axes, respectively. This length scale L is plotted as a function of Γ in figure 6, and it is shown the L decreases approximately linearly with Γ . [START_REF] Leweke | Three-dimensional instabilities in wake transition[END_REF] showed that the elliptic instability theory of [START_REF] Landman | The three-dimensional instability of strained vortices in a viscous fluid[END_REF] predicts wavelengths along the vortex core of

λ = L 2 1 -β tan θ , (3.4)
where θ is the initial angle of incidence of the wavevector which results in the largest growth. For β 0.6, figure 3 of [START_REF] Landman | The three-dimensional instability of strained vortices in a viscous fluid[END_REF] shows that θ is basically a constant, and for β = 0.6 (as measured at the centre of the vortices here) θ 7π/24 53 • . If both β and θ are constant, (3.4) shows that λ is proportional to L. Figure 6 shows that L decreases linearly with Γ , so it is expected that λ should also decrease linearly. Figure 4 shows that this is indeed the case, at least over the range 1.2 < Γ < 2.4.

Concerning the critical Reynolds number, [START_REF] Landman | The three-dimensional instability of strained vortices in a viscous fluid[END_REF] showed that the growth rate of the elliptic instability is made up of two components; an inviscid growth rate σ i and a viscous decay rate (effectively a negative growth rate) σ v . When time is non-dimensionalized by D/U, this viscous decay rate is given by

σ v = 1 Re 2π λ 1 -β cos 2 θ (1 -β) cos 2 θ . (3.5)
At neutral stability, the overall growth rate is zero, such that σ i = σ v . It therefore also seems reasonable to assume that at neutral stability, the viscous decay rate is constant. If σ v is constant, and since both β and θ have been shown to be constant, then from

(3.5) Re c ∝ λ -2 ∝ (1 -ζ Γ ) -2
, where ζ is a constant. It is therefore expected that Re c should increase with Γ , which seems consistent with the trend of Re c observed and plotted in figure 3. The fact that the elliptic instability theory presented by [START_REF] Landman | The three-dimensional instability of strained vortices in a viscous fluid[END_REF] seems to at least be consistent with the trends in λ c and Re c observed during the Floquet stability analysis is further evidence that an elliptic instability plays a role in the development of mode A.

3.5. Calculation of the upper limit of instability for mode B

The upper limit of instability for mode B is difficult to compute. figure 3 shows that this upper limit occurs at much higher Re than the Re c of mode A. The Floquet multiplier of mode A can therefore be quite large, while by definition the Floquet multiplier for mode B at the upper limit is |µ f | = 1. This means that it can become very difficult to isolate mode B when searching for leading eigenvalues and eigenvectors using a Krylov subspace-based method such as that employed here. Of course, if the flow is already three-dimensional through amplification and saturation of other three-dimensional modes, as is the case here, the upper limit of positive growth for mode B will not have strong physical relevance, since the base flow for the stability analysis is no longer two-dimensional.

To estimate the location of this upper limit, the maximum Floquet multipliers have been calculated for mode B at lower values of Re, and the trend extrapolated to find where |µ f | = 1. The data for aspect ratios 1.8 Γ 2.4 are presented in figure 7. The figure plots the maximum value of |µ f | for all wavelengths for mode B as a function of Re. For all values of Γ , it is clear that there is a turning point, and that beyond this turning point the maximum |µ f | decreases with further increases in Re. For the lower values of Γ , e.g. Γ = 1.8, it is possible to resolve mode B beyond the upper limit, and the value of Re at this upper limit can be found from interpolation. However, for higher values of Γ this upper limit needs to be estimated by extrapolating the data.

The relationship between modes A and A

As highlighted in § 3.3, modes A and A appear to coalesce for Γ < 1.2. Figure 2 shows that they are quite similar in spatial structure, especially in the vortices that are close to the body. Only in the vortices further downstream do they differ enough to be identified as distinct modes. As such, care needs to be taken when performing the stability analysis to ensure the two modes are successfully distinguished.

Figure 8 shows the magnitude of the Floquet multiplier |µ f | as a function of wavelength λ for modes A and A at Γ = 1.4 and Γ = 2.0, at a series of values of Re. The change in character of the curves with increasing Re for the Γ = 1.4 case is clear.

At Re = 220, modes A and A have clearly defined 'branches'. Mode A shows a peak in |µ f | at λ 3.5, and |µ f | for mode A increases with increasing λ beyond this.

At Re = 240, modes A and A have clearly defined branches, and there is a peak for |µ f | for mode A around λ = 7.

By Re = 260, the distinction between the branches of mode A and A disappears, and there seems to be a single curve of |µ f | as a function of λ. Modes A and A in this case can only be (weakly) distinguished by inspecting the perturbation fields similar to those presented in figure 2. This also means that, like the case at Re = 220, there is no value of λ at which a peak in |µ f | occurs.

The data for Γ = 2 presented in figure 8(b) show similar issues. While the branches of mode A and A are clearly defined for all of the values of Re shown, lower values of Re result in curves that do not have a clearly defined peak for mode A. Implicitly, figure 8 also shows how the characteristic wavelength of mode A varies with Re. It is shown for both values of Γ that the maximum value of growth changes from occurring at very large wavelengths at lower Re, to relatively small wavelengths at higher Re. In particular, this means that the wavelength of the mode which is actually observed in a flow at a value of Re beyond Re c may be markedly shorter that λ c . This point is shown explicitly in the three-dimensional simulations presented in § 3.7.

The usual method of finding the critical values Re c and λ c involves finding |µ f | max at specified values of λ as a function of Re, and then interpolating for the value of Re where |µ f | max = 1. As |µ f | max is not defined with a distinct peak, applying this method is more difficult for mode A. Instead, for a given value of λ, the value of Re where |µ f | = 1 is found via interpolation. These values are then plotted as a function of λ. The minimum value of Re where |µ f | = 1, and the corresponding value of λ, are identified, and these are Re c and λ c , respectively. An example of the data from this process, for Γ = 1.4, is shown in figure 9.

The exact physical mechanism that leads to instability of mode A is not clear. Due to its structural similarity to mode A near the body, it seems fair to assume that the mechanisms that cause mode A (which is at least partially due to an elliptic instability) also play a role in mode A. However, the wavelength selection and modification in the wake away from the body is possibly due to a number of cooperative mechanisms, and at this stage any description of them is only speculative. 3.7. Three-dimensional simulations confirming the existence of modes A and B The data from the stability analysis, particularly those shown in figure 3, indicate that for higher values of Γ , the two new modes B and A should be the first modes observed with increasing Re.

Full three-dimensional simulations at Γ = 2 have been conducted to investigate this at two values of the Reynolds number: Re = 300 and Re = 350. These simulations are based on extending the two-dimensional spectral-element code and meshes into the spanwise direction using a Fourier series representation as described by [START_REF] Karniadakis | Three-dimensional dynamics and transition to turbulence in the wake of bluff objects[END_REF] and [START_REF] Thompson | Three-dimensional instabilities in the wake of a circular cylinder[END_REF]. The simulations were initialized Visualizations of the flow using isosurfaces of the λ 2 criterion defined by [START_REF] Jeong | On the identification of a vortex[END_REF] for vortex identification are shown in figure 10.

The image in figure 10(a) corresponds to Re = 300, which figure 3 indicates should result in both modes A and B being unstable. However, figure 4 indicates that λ c for mode A for Re = 300 occurs at λ c 20. Therefore, to isolate the modes, the domain for this simulation was restricted to 9.6D so that only mode B could exist within the domain and thus grow. The spanwise domain size was chosen to allow four mode B wavelength to fit within the domain. The image shows that indeed a mode with a wavelength λ = 2.4D grows, which is very close to the wavelength predicted for mode B by the stability analysis as shown in figure 4.

At this Reynolds number, the saturated state of mode B shown in figure 10(a) shows that it manifests as a relatively weak modulation of the forming vortices. This modulation attenuates as the vortices progress downstream.

The simulation at Re = 350 was conducted on a longer domain of 22D. The larger domain was used in this case to allow modes with a wider range of wavelengths to develop, and to observe which modes would naturally dominate the saturated state. Figure 3 shows that modes A, B and A are all linearly unstable at this Reynolds number. As the flow evolves towards saturation, figure 10(b) shows the emergence of a mode with a wavelength λ 7.5D. This is mode A. The wavelength of the mode in the three-dimensional simulation (λ = 7.5D) is much shorter than the critical wavelength for mode A (λ c = 20D). However, mode A first becomes critical at Re = 290, and the three-dimensional simulation was conducted at Re = 350. Stability calculations beyond the critical the critical Re show that the fastest-growing wavelength of mode A rapidly reduces with increasing Re. In fact, stability calculations at Re = 350 predict a wavelength of around 7.5D, similar to that measured in the three-dimensional simulation. The image shows that mode A results in an initial sinuous modulation of the wake vortices, and the rapid development of a complicated series of streamwise vortices at a distance of only a few diameters. Further downstream again, there appears to be complicated vortex pairing and nonlinear interactions occurring.

Interestingly, even though mode B is unstable first with respect to increasing Re, it appears that the wake is dominated by a wavelength attributable to mode A at Re = 350. This is not entirely inconsistent with the predictions from the stability analysis. Figure 7 shows that the Floquet multiplier of mode B never reaches large values, and actually begins to decrease with further increases in Re. Therefore, beyond the initial transition the Floquet multiplier of mode A can exceed that of mode B and dominate the wake.

Concluding remarks

Floquet stability analysis of the flow past an elliptical cylinder has shown that for aspect ratios Γ < 1.4, the three-dimensional transition sequence (mode A, B, QP) is the same as for the circular cylinder, albeit occurring at higher Re with increasing Γ . At Γ > 1.2, the mode A stability branch effectively splits into two segments with two distinct peaks in the Floquet multiplier curves. In fact, the longer-wavelength branch gives mode A, which is the most unstable mode for 1.4D λ 1.9D. For Γ > 1.75, another new mode, named mode B, is found. This mode has the same spatiotemporal symmetry as mode B, but maintains high amplitude in the advecting wake vortex cores, similar to that found for mode A. In terms of structure, mode B has an apparent higher radial wavenumber on the wake vortices than mode A, consistent with its shorter preferred wavelength. A similar mode called mode B was previously found for elongated cylinders with elliptical noses and a flat base [START_REF] Ryan | Three-dimensional transition in the wake of bluff elongated cylinders[END_REF].

Here, this mode is the first occurring for Γ 1.9D.

Based on the stability analysis at a given aspect ratio, at higher Reynolds numbers the mode A and A Floquet branches merge, with the mode A Floquet multiplier regaining dominance. In addition, the maximum value of the mode B multiplier is limited and indeed decreases as the Reynolds number is further increased beyond some point. Thus, beyond the initial transition, it is not clear how strong the influence of the two new modes will be on the fully developed wake state. But in any case, while the circular cylinder transition sequence is maintained for aspect ratios close to unity, further away, this scenario is altered considerably.

FIGURE 1 .

 1 FIGURE 1. Floquet multiplier, µ f , as a function of wavelength λ, at Re = 300 for: (a) Γ = 1.0 (circle); (b) Γ = 1.2; (c) Γ = 1.4; (d) Γ = 1.6; (e) Γ = 1.8; ( f ) Γ = 2.0. Each symbol designates a mode: mode A ( u ); mode B ( p ); mode QP ( E ); mode B (×); mode A ( f ); mode QPL (♦).

FIGURE 2 .

 2 FIGURE 2. Images of the different modes at Γ = 1.8 (except for mode B shown at Γ = 1.4) shown in order of increasing characteristic wavelength. The white and black contours represent spanwise perturbation vorticity over normalized levels of ±0.1. The solid lines show contours of the base flow vorticity at levels ±1. The structure of modes A, B and QP is very similar to their counterparts found in the wake of a circular cylinder. Mode B has a distinct structure, but shares the same spatiotemporal symmetries as mode B.Mode A has a similar structure to mode A in the two vortices closest to the body, but is distinct from the third vortex (circled with a dashed line) and further downstream. Note that modes QP and QPL are shown at an arbitrary phase as they have a complex µ f .

FIGURE 4 .

 4 FIGURE 4. The value of wavelength at marginal stability λ c for (a) mode A ( u ), mode B ( p ), mode QP ( E ), mode B (×) and (b) mode A ( u ), mode A ( f ) as a function of Γ . Lines are best fits to the measured values.

FIGURE 7 .

 7 FIGURE 7. The magnitude of the maximum Floquet multiplier for mode B as a function of Re, for different Γ . Dashed lines show extrapolation to |µ f | = 1 for Re where mode B cannot be resolved due to the much higher Floquet multiplier of mode A.

FIGURE 8 .

 8 FIGURE 8. Magnitude of the Floquet multiplier |µ f | as a function of wavelength λ for (a) Γ = 1.4 and (b) Γ = 2.0. Symbols mark mode A ( u ) and mode A ( f ). The Re for each curve is indicated on the figure. For Γ = 1.4 shown in (a), a peak in the curves for mode A is difficult to discern and the mode A and A curves appear to merge at higher Re. For Γ = 2.0, the two modes are more clearly delineated.

FIGURE 9 .

 9 FIGURE 9. Plot of Re where |µ f | = 1 for mode A, Γ = 1.4. Points represent values interpolated from the data (such as that presented in figure 8), the solid line represents a quadratic fit to these values. The minimum of this curve gives the critical values of Re and λ.

TABLE 1 .

 1 Values of Re and St at the onset of periodic vortex shedding from the current study, compared with those of

	Γ	Re	St	Re (Jackson 1987) St (Jackson 1987)
	1.0 46.559 0.11885 1.2 50.970 0.11974 1.4 56.247 0.12116 1.6 62.174 0.12305 1.8 68.877 0.12582 2.0 75.793 0.12768	45.503 50.586 56.478 --76.794	0.13626 0.13766 0.13943 --0.14644

† Email address for correspondence: justin.leontini@gmail.com

Acknowledgements

The authors would like to acknowledge the financial support of the Australian Research Council (ARC) through grant number DP110102141 under the Discovery program, and the Centre National de la Recherche Scientifique (CNRS) through grant number CNRS/PICS 49586. J.S.L. would like to acknowledge the financial support of the ARC through an Australian Postdoctoral Fellowship.


Xuan-Vu Phan

Laurent Ferro-Famil

Michel Gay

Yves Durand

Marie Dumont

COMPARISON BETWEEN DMRT SIMULATIONS FOR MULTILAYER SNOWPACK AND DATA FROM NOSREX REPORT

INTRODUCTION

This paper presents a multilayer snowpack Electromagnetic Backscattering Model (EBM), based on Dense Media Radiative Transfer (DMRT). This model is capable of simulating the interaction of electromagnetic waves (EMW) at X-band and Ku-band frequencies with multilayer snowpack. The airsnow interface and snow-ground backscattering components are calculated using the Integral Equation Model (IEM), Fung et al. [START_REF] Fung | An update on the iem surface backscattering model[END_REF], whereas the volume backscattering component is calculated by the solution of Vector Radiative Transfer (VRT) equation at order 1. We have applied these models using measurement data from NoSREx report [START_REF] Lemmetyinen | Technical assistance for the deployment of an x-to ku-band scatterometer during the nosrex ii experiment[END_REF], which includes SnowScat data in X-band and Ku-band, TerraSAR-X acquisitions and snowpack stratigraphic profiles. The results of model simulations show consistency with the radar observations, and therefore allow the EBM to be used in various applications, such as data assimilation [START_REF] Phan | 3dvar multilayer assimilation of x-band sar data into a detailed snowpack model[END_REF],

DENSE MEDIA RADIATIVE TRANSFER FOR MULTILAYER SNOWPACK

The DMRT model simulates the multilayer snowpack backscattering coefficient in three components:

σ 0 = σ as + σ v + σ g (1)
where:

σ as and σ sg represent respectively the air-snow and the snow-ground interfaces backscattering, calculated using the Integral Equation Model (IEM) Fung et al. [START_REF] Fung | An update on the iem surface backscattering model[END_REF].

σ v represents the snowpack volume backscattering, calculated by resolving the VRT equation, Tsang et al. [START_REF] Tsang | Theory of microwave remote sensing[END_REF].

The physical parameters of a multilayer snowpack such as thickness, snow density and optical diameter of each layer are considered in the calculation of each layer's effective permittivity ε ef f , which is based on the Strong Fluctuation Theory (SFT). Detailed equations can be found in Huining et al. [START_REF] Huining | Effective permittivity of dry snow in the 18 to 90 ghz range[END_REF]. Next, the extinction coefficient κ e can be derived: κ e = 2k 0 Im √ ε ef f . This factor is then used in the calculation of various scattering mechanisms occurring during the propagation of EMW through a multilayer snowpack, which is categorized into 4 types: (1) transmission between two layers, (2) attenuation by the snow particles, (3) scattering and (4) coherent recombination.

Considering a snowpack made of n distinct layers, where θ k is the incidence angle and d k is the thickness of layer k, the total contribution of the volume backscattering mechanism σ v can be written as follows:

σ v = 4π cos θ 0 n k=1 Att up (k -1)T (k-1)k 1 -exp - 2κ k e d k cos θ k 2κ k e P k T k(k-1) Att down (k -1) (2)
where Att represents the attenuation matrix, T represents the Fresnel trasmission matrix and P represents the phase matrix. The expressions of these factors can be found in Phan et al. [START_REF] Phan | 3dvar multilayer assimilation of x-band sar data into a detailed snowpack model[END_REF].

VALIDATION OF MODEL USING DATA FROM NOSREX REPORT

In this study, data from NoSREx report [START_REF] Lemmetyinen | Technical assistance for the deployment of an x-to ku-band scatterometer during the nosrex ii experiment[END_REF] are used to evaluate the EBM (Electromagnetic Backscatter Model). In the following comparisons, it is worth mentioning that the snowpack stratigraphic profiles are not taken directly at the SnowScat radar footprint and the spatial variability of the snow stratigraphy is high. Consequently, the observations can only be used as a global reference of the evolution of snowpack and not as a precise reference of the snowpack conditions, which are measured by the radar.

In a first part, for comparing EBM model, we consider a single layer snowpack. The snow density is equal to 200 kgm -3 , mean density for the whole period in study and the effective snow grain radius is equal to 0.5 mm. We use the snow water equivalent measured by the gamma wave instrument (GWI) [START_REF] Lemmetyinen | Technical assistance for the deployment of an x-to ku-band scatterometer during the nosrex ii experiment[END_REF]. The ground roughness is fixed σ sg = 1 cm and l sg = 8 cm. Figure 1 shows the time series of simulated backscatter at Ku-and X-band using the EBM in comparison with the SnowScat measurements from Sodankyla in winter 2010/2011. For both cases, Ku-and X-band, the EBM calculations are consistent with SnowScat measurements. However, it is important to note that the model is very sensitive to the grain radius, therefore if we take another value, the results will be changed. Fig. 2. Example of snowpack stratigraphic profiles measurements [START_REF] Lemmetyinen | Technical assistance for the deployment of an x-to ku-band scatterometer during the nosrex ii experiment[END_REF].

In a second part, we take account the snowpack stratigraphic profiles (Figure 2). The snow grain size is measured using visual analysis, which can be subject to human errors. In order to obtain the optical diameter to be used in the EBM, we have divided the measured snow grain size by a factor equal to 2.2 (experimental value, in future studies, we will prefer to use the optical diameter). The ground roughness is fixed σ sg = 0.7 cm and l sg = 8 cm. At this point, we assume that measurements of snow stratigraphy are wrong due to human error. We attempt to improve these measures using a data assimilation process Phan et al. [START_REF] Phan | 3dvar multilayer assimilation of x-band sar data into a detailed snowpack model[END_REF] (not described in this paper).

Figure 5 shows the modifications of grain sizes produced by the assimilation process. The initial profiles come from in-situ snowpack stratigraphic measurements [START_REF] Lemmetyinen | Technical assistance for the deployment of an x-to ku-band scatterometer during the nosrex ii experiment[END_REF] (before vertical white line) and the final profiles result from assimilation (after vertical white line). One can see an overall decrease in grain size, which seems consistent with the finding that visual analysis tends to overestimate this value. Figure 6 shows changes in density produced by the assimilation process. The initial profiles come from in-situ measurements [START_REF] Lemmetyinen | Technical assistance for the deployment of an x-to ku-band scatterometer during the nosrex ii experiment[END_REF] (before vertical white line), thus the final profiles result from assimilation (after vertical white line) . On observation, there is very little change in the density. We have done the assumption that the density measurement errors are low, therefore we have taken a low value in the error covariance matrix. Finally, we used the snowpack modified profiles by assimilation at 30 0 to compare the σ 0 measured by SnowScat at X-band and simulated by EBM at incidence angle of 40 o .

Figure 7 shows that the EBM simulations with assimilation process (red circle) are very consistent with SnowScat measurements (red line).

CONCLUSION

In conclusion, the EBM simulations of backscattering coefficients and measurements from ground-based radar SnowScat as well as TerraSAR-X satellite show the same tendencies at certain sequences of measurement. The validation of this EBM model establishes a relationship between SAR observations and snow stratigraphic measurements (or physical model of snowpack evolution), which allows the development of different approaches, such as data assimilation of SAR data into a detailed snowpack model Phan et al. [START_REF] Phan | 3dvar multilayer assimilation of x-band sar data into a detailed snowpack model[END_REF]. Changes are then made to the initial profiles through a process of assimilation to take into account human error measurement. These promising results also suggest the necessity of an experimental platform that allows automatic and continuous measurements of ground-based radar, coupled with regular observation of the snowpack evolution and SAR satellite acquisitions, in order to validate the electromagnetic model as well as physical models of snowpack and the assimilation process.

Fig. 1 .

 1 Fig. 1. Time series of measured by SnowScat [2] and simulated backscatter using EBM at Ku-(blue) and X-band(red), VV polarized, 40 o incidence angle.

Fig. 3 .

 3 Fig.3. Comparison of σ 0 observed by SnowScat at X-band frequency, TerraSAR-X and simulation of the EBM using insitu snowpack stratigraphic measurements[START_REF] Lemmetyinen | Technical assistance for the deployment of an x-to ku-band scatterometer during the nosrex ii experiment[END_REF].

Figure 3

 3 Figure 3 represents the comparison of the backscattering coefficients σ 0 observed by SnowScat at X-band frequency, TerraSAR-X and simulation of the EBM using in-situ snowpack stratigraphic measurements [2]. One can see that the simulated backscattering coefficients (red cross) varies around the backscattering coefficients measured by SnowScat. The high variations in the value of simulated backscattering coefficients are due to the high sensitivity of the EBM model in the grain size. The tendencies of EBM simulations and SnowScat measurements seem consistent in the period of early March. However, we cannot reach a conclusion with respect to SnowScat and TerraSAR-X data, due to the large spatial variation of snowpack conditions.

Fig. 4 .

 4 Fig. 4. Comparison of σ 0 measured by SnowScat at X-band and Ku-band frequencies and simulation of the EBM using in-situ snowpack stratigraphic measurements [2], at incidence angles of 30 o (top) and 40 o (bottom).

Figure 4

 4 Figure 4 shows the comparison of the backscattering coefficients σ 0 observed by SnowScat at X-band and Ku-band fre-

Fig. 5 .

 5 Fig. 5. Modification of grain sizes profiles (color blue -> red, corresponding to 0 mm -> 2 mm) produced by the assimilation process to converge to SnowScat measurement at X-band (incidence angle of 30 o ).

Fig. 6 .

 6 Fig. 6. Modification of density profiles produced by the assimilation process to converge to SnowScat measurement at X-band (incidence angle of 30 o ).

Fig. 7 .

 7 Fig. 7. Comparison of σ 0 measured by SnowScat at X-band frequencies (red line) and simulation with EBM at incidence angle of 40 o (red circle) using snowpack assimilated profiles at incidence angle of 30 o .

This work has been funded by GlaRiskAlp, the European project (2010-2013) on glacial hazards in the Western Alps, MeteoFrance, the French national meteorological service and Energy de France (EDF).