On the ridge of instability in ferrofluidic Couette flow via alternating magnetic field

There is a huge number of natural and industrial flows, which are subjected to time-dependent boundary conditions. The flow of a magnetic fluid under the influence of temporal modulations is such an example. Here, we perform numerical simulations of ferrofluidic Couette flow subject to time-periodic modulation (with frequency \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Omega _H$$\end{document}ΩH) in a spatially homogeneous magnetic field and report how such a modulation can lead to a significant Reynolds number Re enhancement. Consider a modified Niklas approximation we explain the relation between modulation amplitude, driving frequency and stabilization effect. From this, we describe the system response around the primary instability to be sensitive/critical by an alternating field. We detected that such an alternating field provides an easy and in particular accurate controllable key parameter to trigger the system to change from subcritical to supercritical and vice versa. Our findings provide a framework to study other types of magnetic flows driven by time-dependent forcing.

To date, studies of ferrofluid under alternating magnetic fields are relatively rare and if conducted special attention has been given to their heat behavior 29,30 . Most prominent observation for ferrofluidic flows under alternating fields is the fact that sufficiently high modulation frequency field will force a faster rotation of the particles. This spins up the fluid and thus reduces it's viscosity while in contrast a static field hindrance the free rotation of the magnetic material resulting in an increase of viscosity. In literature occasionally referred to as negative viscosity of the ferrofluid 20,31 . Other works also focussed on the dependence of particle agglomeration under rotating field 32 .
Knowing about the stabilization effect of an applied magnetic field, the question arises how does the system response under an alternating magnetic field for control parameters around the bifurcation threshold of the first instability. To understand the dynamics and system response while "walking" along this ridge of instability between sub-and supercritical states is a main focus of the present work.
Here we numerically study modulated ferrofluidic Couette flow within a wide range modulation frequency and amplitude and observe a significant enhancement in system stability (in Re ≈ 220% ). Further we demonstrate that such an alternating magnetic field provides an easy controllable and quite accurate way to balance the system and walk along the narrow ridge of instability. As such it allows to drive the system to be subcritical or supercritical.

Results
System parameters. In TCS (Fig. 1a) the flow strength is represented in terms of the Reynolds number Re = ω i r i d/ν (the ratio between inertia and viscous forces), which is a very well suited parameter to describe the driving of the system 33 . Here r i is the non-dimensionalized radius and ω i the angular velocity of the inner cylinder. No-slip boundary conditions are used on the cylinder surfaces. The system can be characterized in the cylindrical coordinate system (r, θ, z) by the velocity field u = (u, v, w) and the corresponding vorticity field ∇ × u = (ξ , η, ζ ) . The radius ratio of the cylinders, is kept fixed at 0.5. The time, and length scales are made dimensionless by diffusion time d 2 /ν and gap width d. The pressure in the fluid is normalized by ρν 2 /d 2 .
In the periodically modulated TCS, we give a sinusoidal modulation signal to the external magnetic field (parallel to the system symmetry (z) axis, uniform in space and harmonic in time) as H z = [H S + H M sin (� H t)]e z . As earlier reported such a pure axial magnetic field does not change the system symmetry and only shift the stability thresholds 23,24 . The magnetic field H and the magnetization M are conveniently normalized by the quantity √ ρ/µ 0 ν/d , with free space permeability µ 0 . By using a modified Niklas approach 19,24 the effect of the magnetic field and the magnetic properties of the ferrofluid on the velocity field can be characterized by the (time dependent) Niklas function (see "Methods" section for details) with three control parameters, s z,S being the static contribution of the driving, s z,M the modulation amplitude, and H the modulation frequency. See "Methods" section for more details.
Explored parameter space. We explore the parameter space within s z,S ∈ [0, 1] and s z,M ∈ [0, 1] . The trajectories I and II shown in the parameter space of Fig. 1b represent pure static and pure alternating magnetic fields, respectively. Point A presents the parameters for supercritical flow (Taylor Vortex flow, TVF) at Re = 100 . The trajectories III and IV highlight the parameters at which we provide a more detailed study around the onset of instability at point B for Re = 80 (cf. Fig. 5).

Stability behavior.
Static magnetic fields ( s z,M = 0). Such static fields have been studied in detail in numerous works [22][23][24][25] with the common result that any applied magnetic field regardless it's orientation stabilizes the CCF basic state. Thus the bifurcation thresholds for primary instability (of TVF) are moved to larger Re with increasing field strength s z,S (Fig. 2a). Without magnetic fields, i.e. s z = 0 the critical value is Re 0 c = 68.8 . Worth mentioning that other axial wavenumber (here k = 3.927 ) will lead to other critical Reynolds numbers Re 0 c .
(1) Modulated magnetic fields ( s z,M = 0). Similar to increasing field strength s z,S for static magnetic fields, also an increase in the modulation amplitude s z,M stabilizes the system. Figure 2b presents the surface of Re c (over (s z,S , s z,M )-plane) which is convex in all points in any direction. Thereby the quantity of stabilization also increases with increasing modulation amplitude (Fig. 2c). Although this behavior remains qualitative the same for alternating magnetic fields with different static contributions s z,S , the relative effect weakens with increasing the static contribution. For parameters in Fig. 2 the maximum stability enhancement in Re is about 220% , comparing the system in absence of any magnetic field with alternating magnetic field at (s z,S = 1 = s z,M ) . The inset in Fig In terms of stability one can summarize, that the system reacts to an alternating modulation of the magnetic field similar as increasing the magnetic field strength in the static case. The stronger stabilization with increasing modulation amplitude originates from the static field behavior in particular it's non-linear grows with power of 2 ( Fig. 2a). Thus during one modulation period the system experience a stronger stabilization effect while the modulation amplitude is above the average field strength in comparison to the de-stabilization in the other half period. As result (for high frequency) the stabilization within an alternating magnetic field corresponds to a static field strength, which lies above the mean value of the alternating field. For the same reason the stabilization also grows with increasing modulation amplitude (Fig. 2c). Figure 3 illustrates the stable forward bifurcating branches of TVF solutions for different modulation amplitudes s z,M as indicated ( s z,S = 0 ) of the magnetic field. The onsets corresponds to the critical curve in absence of any magnetic field s z,S = 0 ( Fig. 2 in main paper). Being supercritical the dominant mode amplitudes |u 0,1 | grow in well known square root manner. For better comparison we also consider use the relative distance µ = Re(s z,S )/Re 0 c − 1 (inset), with Re 0 c (depending on system parameters, e.g. the axial wavenumber k) being the Reynolds number for which the flow becomes supercritical in the absence of a magnetic fields.

Bifurcation behavior.
As a result one can say, increasing the modulation amplitude s z,M moves the onset of instability to larger control parameters (Re, to the right in Fig. 3) and therefore stabilizes the CCF basic state. In addition it also slightly effects/modifies the bifurcation characteristics itself. Rescaling the bifurcation scenario by the corresponding onsets (insets in Fig. 3) one sees that with increasing modulation amplitudes s z,M also the mode amplitudes grow faster, the corresponding slopes become steeper. Similar observation numerical and experimental has been already found for increasing field strength in static magnetic fields [22][23][24] .
A further bifurcation scenario of the dominant mode amplitudes |u 0,1 | for static non-zero, and two different modulated driven magnetic fields (also with finite static contribution) is shown at the inset in Fig. 2a ( H = 100, Re = 80 ). The chosen parameters (cf. horizontal arrows in Figs. 1b, 2a) correspond to those for which in the following the dynamic system response will be investigated. For larger modulation amplitude s z,M = 0.2 , TVF dissapear at smaller corresponding static control parameter s z,S and the system returns to the CCF basic state. This is another confirmation for stabilization effect with increasing modulation amplitude s z,M . In the high-frequency limit, solely the time average of s z (t) affects the stability behavior. Thus, in this limit the stability boundary coincides with a static stability boundary using an equivalent static magnetic Niklas parameter. Note, that this larger than the mean value s z (t) T H (cf. Fig. 4a). For the given field the order parameters for equivalent static driving is s z = 0.245 , which, for the sake of reference is also included (red dashed lines) in Fig. 4 (Note here �s z (t)� T H = 0.2 ). For the modulation with the high frequency H 100 , the flow dynamics is nearly averaged. Variations in the dominant mode amplitude |u 0,1 | are small compared to its mean value. For H = 100 the modulation amplitude |u 0,1 | is barely 0.29% of its time mean (Fig. 4b). A phase shift between the maximum and minimum of field function s z (t) versus the minimum and maximum of the mode amplitudes |u 0,1 | occurs: the latter ones are temporally delayed to the former because of the inertia of the fluid resisting the fast changing accelerating Kelvin force leading to this time lag. Consistently the phase shift decreases with decreasing frequency (best visible for |u 0,1 | in Fig. 4b). Thereby the oscillation amplitudes are increasing with smaller H . The lower the modulation frequency, the closer the oscillation profiles get to the curve of a static magnetic field (red squares). Deviations just persist in the vicinity of the bifurcation threshold, because the dynamics become infinitely slow there.
The inharmonic behavior in the mode amplitudes |u 0,1 | (Fig. 4b) for very low frequency (and the static limit) reflects the increasing effect onto the flow dynamics with increasing field strength s z (t) (Fig. 2). For �s z (t)� = 0.  Figure 4b reflects this by either steeper/larger variation |u 0,1 | for positive modulation amplitude s z,M > 0 as well as a much flatter profile |u 0,1 | for negative modulation amplitude s z,M < 0 . The latter is a direct consequence of smaller variation in Re c with small field strengths (Fig. 2a). Worth mentioning, although not further studied such inharmonic response behavior has been earlier reported in Rayleigh-Bénard system exposed to a time-periodic magnetic field 34 .
Interesting observation is the fact that for low frequencies H , approaching the static state the mode amplitudes |u 0,1 | within one period slightly overshoot the maximum and minimum values of their static counterparts. For high frequencies H 30 the mode amplitudes |u 0,1 | move around the average well within their maximum and minimum limits. It is the inertia of the fluid itself which causes this overshooting.  Small modulation amplitude ( s z,M = 0.1). The system becomes only slightly subcritical over one period. In the high frequencies limit the time averaged magnetic field s z (t) for modulated driving (dashed red line in Fig. 51c) corresponds to a static magnetic field with s z,S ≈ 0.423 (cf. inset in Fig. 2). With decreasing frequency H first the amplitude in the oscillating mode |u 0,1 | continuously increase before at H 0.5 it eventually becomes tem-    www.nature.com/scientificreports/ porally zero indicating that the system is now subcritical. The smaller the frequency H , the longer the system remains subcritical (Fig. 51c). Over one period, for such low frequencies, a fast growth of the mode amplitude |u 0,1 | followed by a relaxing just similarly to values close to the stationary can be observed. With further decreasing H the oscillation profile in the mode amplitudes |u 0,1 | approaches the static scenario. As for full supercritical flow state (Fig. 4), a temporally delay between the extrema of s z (t) and the corresponding extrema (min and max) in the mode amplitudes |u 0,1 | appears.
Large modulation amplitude ( s z,M = 0.2). The system goes deeper into the subcritical regime within one period of driving ( Fig. 5(2)) and as a result remains subcritical in the high frequency limit (see inset in Fig. 2), which is just opposite to the scenario for small modulation amplitude (1). The decay of the mode amplitudes |u 0,1 | in Fig. 5b for larger frequencies H clearly illustrates the subcriticality for here given alternating field. Corresponding equivalent static magnetic field is s z,S ≈ 0.407 . However, with decreasing H the system becomes temporally supercritical. Analogous fast growth of the mode amplitude |u 0,1 | followed by a relaxing close to the stationary state appears. Again, the lower the modulation frequency H , the larger the oscillation amplitudes |u 0,1 | , thereby approaching the limiting value given by the respective stationary solution curve (red squares in Fig. 5(2c)).

Discussion
We have shown that ferrofluidic Couette flow under alternating magnetic field becomes stabilized. The primary instability of TVF is moved towards larger Re, whereby the quantity of stabilization increases with larger modulation amplitude. This is similar to the modification due to static magnetic fields with increasing field strength.
With increasing oscillation frequency the temporal evolution/response in the system decreases. The stability boundary for alternating magnetic field in high-frequency limit corresponds to a static stability boundary which is above the mean of the alternating magnetic field. This results from the fact that during one modulation period the system experience a stronger stabilization effect while the modulation amplitude is above the average field strength in comparison to the de-stabilization in the other half period. For very low modulation frequencies the oscillation profiles approach the stationary curves.
In addition, we found that the system response is selective to driving parameters around the primary instability. As such an alternating magnetic field can force/drive the system to be subcritical or supercritical.
The schematic in Fig. 6 summarizes the non-linear system response based on small and large modulation amplitudes with respect to variation in the driving frequency H . In any case the high frequency limit selects a single solution, the system is either sub-or supercritical. Which of the solution is selected (subcritical or supercritical) depends on the modulation amplitude. For studied parameters, s z,S = 0.4 and modulation amplitudes s z,M = 0.1 and s z,M = 0.2 , respectively (Fig. 5) the selection appears at H ≈ 0.27 . Main characteristics while "surfing the edge of instability" can be described as follows: For small magnetic modulation amplitudes the system is supercritical for high frequencies H . Decreasing the frequency modifies this scenario and the system becomes temporally sub-and supercritical. On the other hand for large modulation amplitudes the system is subcritical in the high frequency limit. However, decreasing the frequency also modifies the system response to be temporally subcritical and supercritical.
The present work highlights the importance of complexer fluids under external driving. As such the variation in frequency of the alternating field provides a very simple and in particular accurate controllable way to trigger the system response to be either subcritical or supercritical. This offers various ways for industrial applications, e.g. focussing on the significant difference in torque between the subcritical CCF basic state and the primary instability of supercritical TVF.

Methods
Direct numerical simulation. DNS for ferrohydrodynamical flow using the Niklas approximation are employed 19,24 . In the present study, we consider in axial direction periodic boundary conditions corresponding to a fixed axial wavenumber k = 3.927 , motivated by experimental findings for the appearance of primary TVF instability in Taylor-Couette flow with outer cylinder at rest 3,33 . No-slip boundary conditions are used on the cylinder surfaces and the radius ratio of inner and outer cylinders, is kept fixed at r i /r o = 0.5 . The DNS are conducted combining a standard, second-order finite-difference scheme in (r, z) with a Fourier spectral decomposition in θ and explicit time splitting. The explored parameter range spans 68 Re 152 , 0 (s z,S , s z,M ) 1 , . For these parameters, the choice of 16 azimuthal modes provides adequate accuracy. We use a uniform grid with spacing δr = δz = 0.02 and time steps δt < 1/3800.

Ferrohydrodynamical equation of motion.
The non-dimensionalized hydrodynamical equations 25,27,35 are given by: On the cylindrical surfaces, the velocity fields are given by u(r i , θ, z) = (0, Re, 0) and u(r o , θ, z) = (0, 0, 0) , with the the Reynolds numbers Re = ω i r i d/ν , where r i = R i /(R o − R i ) is the non-dimensionalized inner cylinder radius.
Equation (2) is to be solved together with an equation that describes the magnetization of the ferrofluid. Using the equilibrium magnetization of an unperturbed state where homogeneously magnetized ferrofluid is at rest and the mean magnetic moment is orientated in the direction of the magnetic field, we have M eq = χH . The magnetic susceptibility χ of the ferrofluid can be approximated with the Langevin's formula 36 , where we set the initial value of χ to be 0.9 and use a linear magnetization law. The ferrofluid studied corresponds to APG933 37 . We consider the near equilibrium approximations of Niklas 19,38 with small ||M − M eq || and small magnetic relaxation time τ : |∇ × u|τ ≪ 1 . Using these approximations, one can obtain 27 the following magnetization equation: where is the Niklas coefficient 19 , µ is the dynamic viscosity, is the volume fraction of the magnetic material, S is the symmetric component of the velocity gradient tensor 27,35 , and 2 is the material-dependent transport coefficient 35 , which we choose to be 2 = 4/5 24,35,39 based on experimental observation. Using Eq. (3), we can eliminate the magnetization from Eq. (2) to obtain the following ferro-hydrodynamical equations of motion 25,27,35 : where F = (∇ × u/2) × H , p M is the dynamic pressure incorporating all magnetic terms that can be expressed as gradients, and s z is the Niklas parameter (Eq. (7)). Note, while in earlier studies considering static magnetic field this is a real parameter, in the present work devoted to alternating magnetic fields this is basically a timedependent function, which we will refer to as Niklas function. To the leading order, the internal magnetic field in the ferrofluid can be approximated as the externally imposed field 25 , which is reasonable for obtaining dynamical solutions of the magnetically driven fluid motion. Equation (5) (∂ t + u · ∇)u − ∇ 2 u + ∇p M = s 2 z ∇ 2 u −   www.nature.com/scientificreports/ is solved with the FTCS (Forward Time, Centered Space) algorithm 40 . Further pressure and velocity fields are iteratively adjusted to each other with the method of "artificial compressibility" 41 .
The pressure correction dp (n) in the nth iteration step being proportional to the divergence of u (n) is used to adapt the velocity field u (n+1) . The iteration loop (Eq. (10)) is executed for each azimuthal Fourier mode separately. It is iterated until ∇ · u has become sufficiently small for each m mode considered-the magnitude of the total divergence never exceeded 0.02 and typically it was much smaller. Time steps were always well below the von Neumann stability criterion and by more than a factor of 3 below the Courant-Friederichs-Lewy criterion.
Hereafter the next FTCS time step is executed. For diagnostic purposes, we also evaluate the complex mode amplitudes f m,n (r, t) obtained from a Fourier decomposition in the axial direction: where k = 2πd/ is the axial wavenumber.