Sustained turbulence and magnetic energy in non-rotating shear flows

0.84 MB·

by  Farrukh Nauman

Sustained turbulence and magnetic energy in non-rotating shear flows
Farrukh Nauman1,2,∗ and Eric G. Blackman1,†
1Department of Physics and Astronomy, University of Rochester, Rochester, NY 14627, USA.
2Niels Bohr International Academy, The Niels Bohr Institute,
Blegdamsvej 17, DK-2100, Copenhagen Ø, Denmark.
(Dated: March 10, 2017)
From numerical simulations, we show that non-rotating magnetohydrodynamic shear flows are
7 unstabletofiniteamplitudevelocityperturbationsandbecometurbulent,leadingtothegrowthand
1 sustenanceofmagneticenergy,includinglargescalefields. Thissupportstheconceptthatsustained
0 magneticenergyfromturbulenceisindependentofthedrivingmechanismforlargeenoughmagnetic
2 Reynoldsnumbers.

r
a
M INTRODUCTION magnetized Keplerian shear is linearly unstable to the
magnetorotationalinstability(MRI)([10–13])whichdoes

8 sustain growth. However, Ref. [9] employed low resolu-
Shear flows are common in nature, both rotating and
tion ideal MHD simulations, and this problem of linear
] non-rotating. Rotationisessentialwhenangularmomen-
h shearwasnotstudiedforconvergence. Theirconclusions
tum support causes the shear, and most studies of field
p also lead to a cognitive dissonance: if turbulence from
- growth have focused on the former. But shear flows in
m which rotation is inessential are also ubiquitous. Exam- linear shear flows were distinctly unable to sustain mag-
netic energy, it would contradict a lesson from stochas-
s plesfromastrophysicsoccurneartheinterfaceofoutflows
a tically forced turbulence where saturated magnetic en-
propagatinginto ambientmedia [1], nearthe interfaceof
l ergy achieves near equipartition with turbulent kinetic
p outward and inward convective plumes in disks or stars
s. [2]andbetweenturbulenteddiesingalaxiesorclustergas energyfor largeenoughmagnetic Reynolds number, Rm
([14, 15]). We are thus motivated to revisit this non-
c [3]. For some azimuthal shear flows, the role of rotation
i rotating magnetized shear problem with more compre-
s may also be minimal [4]. Magnetic fields are common
y hensive simulations.
in all of these contexts. Whether linear shear flows can
h There is also long standing interest in understanding
generate turbulence, amplify magnetic energy, or even
p theroleofshearinthegenerationoflargescalemagnetic
[ produce large scale fields are all questions that interface
into the long standing questions of magnetic field am- fields (e.g., [16–22]). Ref. [16] was the first to show
2 numericallyusingashearingbox,thatthecombinationof
plification in astrophysics and identifying the minimum
v non-helicalstochasticforcingpluslinearshearcanleadto
conditions needed for field amplification in magnetohy-
1
largescale dynamo in a shearingbox. The forcing in the
3 drodynamics(MHD) (fora reviewofdynamotheory,see
simulationsofRef. [16]is suchthat the stochasticpower
5 [5]).
3 inputwasmuchstrongerthantheshearforcing,andscale
Whileamplificationofmagneticfieldsinstochastically
0 separationwasachievedthroughthe use of largevertical
. forced flows with and without shear has been demon- domains. Buttheaforementionedstudiesofnon-rotating
1
strated, questions of nonlinear stability and magnetic
0 linear shear and large scale magnetic field growth have
fieldsustenanceinnon-rotatingmagnetizedshearflowin
7 employedtheadditionalstochasticforcingastheprimary
1 three dimensions have received little attention (for two source of turbulence. This contrasts our present work.
: dimensions, see [6]). Purely hydrodynamic linearly sta-
v In this paper, we study the nonlinear stability of non-
bleshearflowsdoindeedtransitionintoaturbulentstate
i rotating magnetized shear flow using a suite of numeri-
X for a variety of hydrodynamic systems at high Reynolds
calsimulationsinashearingbox,withoutanyadditional
r number, Re(forexample,planeCouetteflow(PCF):[7],
a stochastic forcing. We explore the sustenance of turbu-
pipe flow: [8]) and without a net magnetic flux, the per-
lent state as well as the creation of large scale magnetic
turbed velocities are affected by magnetic field fluctua-
fields. The behavior of velocity and magnetic fields is
tions only to second order. The linear stability problem
studied as a function of box size and the dissipation co-
thus reduces to that of hydrodynamic PCF but can the
efficients. Most significantly we find that linear shear
resulting turbulence sustain magnetic energy?
flows unstable to turbulence do indeed sustain magnetic
Previously,Ref. [9]foundthatwhilemagnetizedlinear energy for large enough magnetic Reynolds numbers.
shear flows do indeed exhibit flow turbulence, magnetic
energy was found to decay. This was interpreted to sug-
gest that linear shear flows may be intrinsically unable METHODS AND RESULTS
togrowfields intheabsenceofrotation. While the Cori-
olis force stabilizes hydrodynamic Keplerian shear flow We perform direct numerical simulations (DNS) em-
(Rayleigh criterion, e,g. [10]), Ref. [9] emphasized that ploying a shearing box setup (with no other forcing) to
2
1100--22 110000
y 1100--44 y 1100--22
g g
r r
ne 1100--66 ne 1100--44
E E
1100--88 RRmm==11020000 1100--66 1x2x1
Kinetic Rm=1500 Kinetic 1x2x4
1100--1100 Magnetic Rm=2000 1100--88 Magnetic 4x8x4
00 220000 440000 660000 880000 11000000 00 220000 440000 660000 880000 11000000
Time (units of S-1 ) Time (units of S-1 )
FIG. 1: Time history plot for different runs. The kinetic FIG. 2: Same as 1 for 1×2×1 (red, middle two lines),
energy is represented by solid lines that are all at the 1×2×4 (green, top two lines), 4×8×4 (blue, bottom
top of the plot. The dotted lines represent magnetic two lines) at Rm=10,000.
energy. The bottom most line (red, dotted) represents
Re=Rm=1,000 that decays. The Re=Rm=1,200
(green, dotted) line is second from bottom, 1,500 (blue, energy decays immediately after reaching the saturation
dotted) thirdfrombottom and2,000(khaki,dotted) for state, (3) Re > Re and Rm > Rm : Both kinetic
crit crit
a domain size of 1×2×1 with 32×64×32 zones. Both and magnetic energy sustain growth. Fig. 1 shows that
kinetic and magnetic energy are normalized by SL2x. Rmcrit ∼1,200. We estimated Recrit ∼750 for both hy-
drodynamicandMHDruns,whichisconsistentwiththe
value found in the hydrodynamic simulations of PCF [7]
study non-rotating magnetized linear shear flow using (notethatthedefinitiontypicallyusedinPCFliterature
the pseudospectral code snoopy [23] ([24]). We define is a factor of 2 smaller than our definition). These crit-
Re = L2xS/ν, Rm = L2xS/η where Lx = 1 is the size ical values have also been verified at a higher resolution
of domain in the ‘x’ direction, η is the magnetic dif- of 64×128×64. Note that the finite lifetime of turbu-
fusivity, ν is the kinematic viscosity and the shear pa- lence as seen in the kinetic energy of Re =Rm = 1,000
rameter, S = 1. We set the magnetic Prandtl number (red) run in fig. 1 is consistent with what has recently
Pm = Rm/Re = 1 for all runs. We initialize our simu- beenfoundinhydrodynamicshearflowexperiments[26].
lations with zero net magnetic flux Bini = B0sinkxxez Ref. [7] suggests that turbulence in linear shear flows
(whereB0 =0.035)andapplyfiniteamplitudeperturba- in small domains (L ,L ∼ L ) exhibit transient chaos,
y z x
tions (LS) to large scales in the velocity [25]. The shear while large aspect ratio domains L ,L ≫ L would in-
y z x
profile Vsh =−Sxey is subtractedout of the totalveloc- stead abruptly transition into steady turbulence. We do
ityandthevelocitythecodesolvesforisV =Vtotal−Vsh: notexploreextendeddomainshereinsoourresultswould
∂V ∂V represent a lower limit on the robustness of turbulence.
+V +∇·(VV +T)=−SV e +ν∇2V,
sh x y
∂t ∂y
(1)
∂B Domain Size effects
=∇×(V ×B)+η∇2B, (2)
∂t
∇·V =0, ∇·B =0, We plotthe ratioofkinetic energyto magnetic energy
in fig. 2 for three different domain sizes (see table I for
where T = (p+B2/2)I −BB. The shear time unit is description)1×2×1(643), 1×2×4(64×64×256)and
1/S. 4×8×4 (2563). The table also has data for a higher
resolution run for the domain 1×2×1, which suggests
convergence for this aspect ratio. The magnetic energy
Critical Re and Rm is nearly an order of magnitude smaller than kinetic en-
ergyin 1×2×1 and1×2×4 whereasin the the largest
In figure 1, we plot the time history of the kinetic and domain,4×8×4thetwoarenearlyequal. Forthe range
magnetic energies for runs with resolution 32×64×32. of domains studied so far, this energy ratio therefore de-
We identify three distinct regimes: (1) Re < Recrit: the pends not only on the aspect ratio Lz/Lx but increases
flowremainslargelylaminarandtheinitialperturbations with box size for a fixed aspect ratio.
die off, (2) Re > Re but Rm < Rm : kinetic en- The velocity profiles hV i, hV i and the magnetic field
crit crit x y
ergygrowsandsustainsforsometimewhilethemagnetic profiles hB i, hB i averaged over xy are plotted for the
x y
3
FIG. 3: xy averagedhV i, hV i, hB i and hB i for the run 1×2×1 for the first 1000 shear times. The x-axis is the
x y x y
time (1/S units) and the y-axis is the vertical domain size in L units.
Domain Size Resolution [B2] ∂hExi/∂zz+ Rm 4×8×4rundisplaysconsiderablelargescaleorganization
1×2×1 64×64×64 0[.V220]4 hS0B.0x1i7z+ 10,000 compared to the corresponding plot for 1×2×4. The
1×2×1 128×128×128 0.211 0.007 10,000 magnetic fields are strongly correlated with the velocity
1×2×4 64×64×256 0.062 0.002 10,000 fields (see Appendix) for 4×8×4 since all of them have
4×8×4 256×256×256 2.413 0.011 40,000 sinusoidalstructureontheboxscale. Moreinterestingly,
the magnetic to kinetic energy ratio is nearly unity (see
TABLE I: Description of the three runs analyzed table I). This is in contrast to the 1×2×4 run where
herein, with one higher resolution run done for the magnetic energy is more than anorder ofmagnitude
convergence test. The third column is the ratio of smallerthanthekineticenergyandthehBxiprofileseems
magnetic to kinetic energies [B2]/[V2] while the fourth to be very noisy. The cross helicity of the 1×2×4 run
column shows that the shear dominates the EMF is close to unity for a significant duration, while that of
derivative term ∂hE i/∂z in the induction equation for 4×8×4 is smaller in comparison and fluctuates about
x
the y-component of the magnetic field ([Q]= volume zero suggesting that the 1×2×4 run is dominated by
average of Q, Q= time average of Q from 100S−1 to kz =1 mode (see Appendix for plots).
200S−1, Qz+ represents the vertical average from z =0 The spatiotemporal profile of hByi in fig. 4 suggests
to z =+L /2, and xy averageof Q is represented by the existence of a cycle period and thus a large scale
z
hQi.). The last column lists the magnetic Reynolds dynamo. In the xy averagedinduction equation:
number for the three different runs.
∂hB i ∂ ∂2
x =− hV′B′ −V′B′i+η hB i
∂t ∂z z x x z ∂z2 x
∂hB i ∂ ∂2
run with 1×2×1 domain in fig. 3. Unlike the magnetic y =ShB i+ hV′B′ −V′B′i+η hB i (3)
∂t x ∂z y z z y ∂z2 y
fields, we see that while the velocities are dominated by
a sinusoidal profile in ‘z’ that varies in time. Further- the only terms that can contribute to the right side
more, the hBxi profile is more noisy than hByi. The of the hByi equation are the EMF term ∂zhExi (where
simulation began with a shear profile Vsh =−Sxey, and hEi = hV′ × B′i, V′ and B′ are fluctuations result-
eventually reached a steady state with additional shear ing from xy averaging) and the ‘Omega’ term ShB i
x
in the z−direction for the xy averaged velocity fields, as seen in eq. 3 (the mean field term contribution,
hViver. sh ∼sinkzz(ex+ey). This structure is a generic ∂(hVi×hBi)i/∂z ∼ 0, where i = x,y). We estimated
feature of hydrodynamic shear flows at and just above theformertoberoughlytwoordersofmagnitudesmaller
Recrit. Recent work on the transition to turbulence sug- than the shear term (table I)
geststhatasthedomainsizeandReareincreased,these We plot the power spectra of the velocity and mag-
structures disappear into ‘featureless’ turbulence [7, 27]. netic field averaged over xy and Fourier transformed in
‘z’ in fig. 5. The predominance of the k = 1 (in units
z
To explore possible generation of an organized mag- of (2πn /L )) mode is consistent with the nearly sinu-
z z
netic field, we plot the magnetic field profilesfor the two soidal profile of hV i seen in fig. 3 and hB i in 4. The
y y
larger domains with L = 4: 1×2×4 and 4×8×4 velocity spectrum for 1×2×4 seems to follow the 1D
z
in fig. 4. The ‘y’ component of the magnetic field, hB i Kolmogorov scaling k−11/3 for intermediate wavemodes,
y
seems to be very coherent in both runs and hB i in the which would imply that these modes represent the iner-
x
4
FIG. 4: xy averagedhB i and hB i for the runs 1×2×4 (left column) and 4×8×4 (right column) for the first
x y
1000 shear times. The x-axis is in the units of 1/S, while the y-axis represents the vertical domain in L units. Note
a strong anti-correlationbetween the hB i and hB i for both runs.
x y
using periodic boundaries and thus the large power ob-
served in box scale structures is an indication that the
102
boundary conditions are strongly influencing the flow
V 1x2x1
y
100 By 1x2x1 dynamics. It remains to be explored whether magnetic
m k−11/3 Vy 1x2x4 fields in such high Re turbulent flows with a featureless
B 1x2x4
Spectru 1100−−42 VByyy 44xx88xx44 vdgeroolowncotithtyeinpthrsoahfitelaererwicneognutbldoaxnsehasloywsthicolawtrhsgeetohsracyatlfreoorotarlgtaiarognneiziassctniaoolnet.finWeelcde-
er essary for dynamo action when a source of velocity fluc-
ow 10−6 tuations is present in a shear flow [29]. Our simulations
P
satisfy their minimum sufficient conditions, although we
10−8 focus on xy averaging rather than the yz averaging of
their case.
10−10
0.1 1.0 10.0 100.0
k
z
FIG. 5: Power spectra of xy averagedazimuthal
CONCLUSIONS
velocity hV i and magnetic field hB i averagedfrom 501
y y
to 600 shear times as a function of k (in units of
z
(2πn /L )). The top (green) two lines represent Using high resolution 3-D simulations of a shearing
z z
1×2×4, while the middle (red) two represent1×2×1 box with a pseudospectral code, we have demonstrated
and 4×8×4 (blue) lines are at the bottom. The power numerically for the first time that not only does shear
spectra are normalized such that P |hQi(k )|2 =[Q2]. driven turbulence sustain for high enough Re, but this
kz z turbulence amplifies and sustains magnetic energy when
Rm is large enough. This contrasts the work of Ref. [9]
whodidnotidentifysustainedgrowthinmagneticenergy
tial range. However, the velocity power spectrum peaks
because their Rm was below the critical value we have
atthe box scale,whichcouldbe due to strong2Dvortex
found. The turbulence emergent in our simulations is
structures [28]. The velocity structures in the smallest
self-sustained by the linear shear and thus distinct from
domain 1 × 2 × 1 do not seem to follow the 1D Kol-
adifferentclassofworkthatemployedstochasticforcing
mogorov scaling, while the largest domain 4×8×4 has
in additionto the non-rotatingshear [16–22]. Structures
a steeper power law behavior for the most part but ap-
in both velocity and magnetic fields at the largest scales
pears to have a flatter power law spectrum closer to the
are seen in our largest domains and we have identified
dissipation scale.
the EMF terms that sustain the latter. Whether the
Since the velocity fields are dominated by box scale velocity structures break into featureless turbulence at
structures, it becomes a subtle matter to define and dis- evenhigherReynoldsnumbersanddomainsizesremains
tinguish small vs. large scale dynamos [29] or system to be explored, but the minimum ingredients derived for
scale dynamo [30]. An additional caveat is that we are large scale field growth [29] are met.
5
1.0
1x2x1 ∗ [email protected]
1x2x4
4x8x4 † [email protected]
0.5 [1] E. Liang, M. Boettcher, and I. Smith,
ApJL 766, L19 (2013), arXiv:1111.3326 [astro-ph.HE].
+ [2] M. S. Miesch and J. Toomre,
z>
b 0.0 AnnualReview of Fluid Mechanics 41, 317 (2009).
⋅<v [3] M. Bru¨ggen and F. Vazza, in
Magnetic Fields in Diffuse Media, Astrophysics and
Space Science Library, Vol. 407, edited by A. Lazarian,
−0.5
E. M. de Gouveia Dal Pino, and C. Melioli (2015) p.
599.
[4] P. S. Cally, Sol. Phys. 194, 189 (2000).
−1.0
0 200 400 600 800 1000 [5] A. Brandenburg, D. Sokoloff, and K. Sub-
Time (units of S−1) ramanian, Space Sci. Rev.169, 123 (2012),
arXiv:1203.6195 [astro-ph.SR].
[6] G.R.Mamatsashvili,D.Z.Gogichaishvili,G.D.Chagel-
FIG. 6: Time evolution of horizontally averagedcross
ishvili, andW.Horton,Phys.Rev.E 89, 043101 (2014),
helicity hv·biz+ vertically averagedfrom z =0 to
arXiv:1409.8543 [physics.plasm-ph].
z =+Lz/2 plotted for 1×2×1 (red, solid), 1×2×4 [7] P.Manneville,European Journal of Mechanics - B/Fluids 49, Part B
(green, dashed), 4×8×4 (blue, dotted) for the first trendsin HydrodynamicInstability inhonourof Patrick
1000S−1. Huerre’s 65th birthday.
[8] T.Mullin,AnnualReview of Fluid Mechanics 43, 1 (2011),
http://dx.doi.org/10.1146/annurev-fluid-122109-160652.
[9] J. F. Hawley, C. F. Gammie, and S. A. Balbus,
Astrophys. J. 464, 690 (1996).
Acknolwedgements:- WethankF.Ebrahimi,P.Bhat,
[10] N. Shakura and K. Postnov, MNRAS 448, 3697 (2015),
and S. Tobias for related discussions. The simulations
arXiv:1412.1223 [astro-ph.HE].
reported on in this paper were done on the Blue Streak [11] E. P.Velikhov,JETP 36, 995 (1959).
cluster hosted by the Center for Integrated Research [12] S.Chandrasekhar,Proceedings of theNational Academy of Science 46
Computing at the University of Rochester. FN acknowl- [13] S. A. Balbus and J. F. Haw-
edgesfundingfromtheEuropeanResearchCouncilunder ley, Astrophys.J. 376, 214 (1991);
E. G. Blackman and F. Nauman,
the European Unions Seventh Framework Programme
Journal of Plasma Physics 81, 395810505 (2015),
(FP/2007-2013)underERCgrantagreement306614. EB
arXiv:1501.00291 [astro-ph.HE].
acknowledges support from grants HST-AR-13916.002
[14] A. A. Schekochihin, S. C. Cowley, S. F. Tay-
and NSF-AST1515648. lor, J. L. Maron, and J. C. McWilliams,
Astrophys. J. 612, 276 (2004), astro-ph/0312046.
APPENDIX
1.0
1x2x1
1x2x4
We plot the cross helicity averaged from z = 0 to 4x8x4
z =+L /2representedbyhv·biz+ (hv·biz−: z =−L /2 0.5
z z
toz =0;[v·b]: volumeaverage)forthethreerunsinfigs.
−
6,7,8. Itappearsthatthe1×2×4runhasacrosshelic- z>
b 0.0
ity dominated by the largest mode and thus has nearly ⋅v
<
maximal cross helicity with the same sign for different
verticalsectionsoftheboxforaconsiderabledurationof
−0.5
time. This is further supported by the velocity profiles
of this run that are also attached in fig. 9. The velocity
profiles for this run are all seemingly locked into kz = 1 −1.0
0 200 400 600 800 1000
state,similartothemagneticfieldprofilesinfig. 4ofthe
Time (units of S−1)
main text. Of the 3 runs, the 1×2×4 is the only one
with a dominant vertical extent (L /L >1). The large
z x FIG. 7: Same as fig. 6 but for hv·biz− vertically
vertical extent seems to be required for the appearance
averagedfrom z =−L /2 to z =0.
of this large scale dominant mode. z
6
[15] N. E. Haugen, A. Brandenburg, and W. Dobler, [27] D. Barkley and L. S. Tuckerman,
Phys.Rev.E 70, 016308 (2004), astro-ph/0307059. Physical Review Letters 94, 014502 (2005),
[16] T. A. Yousef, T. Heinemann, A. A. Schekochi- physics/0403142.
hin, N. Kleeorin, I. Rogachevskii, A. B. [28] H. Xia, D. Byrne, G. Falkovich, and M. Shats,
Iskakov, S. C. Cowley, and J. C. McWilliams, Nature Physics 7, 321 (2011).
Physical Review Letters 100, 184501 (2008), [29] F. Ebrahimi and E. G. Black-
arXiv:0710.3359. man, MNRAS 459, 1422 (2016),
[17] S. Sridhar and K. Subramanian, arXiv:1509.04572 [astro-ph.HE].
Phys.Rev.E 79, 045305 (2009), arXiv:0812.3269. [30] S. M. Tobias, F. Cattaneo, and N. H. Brummell,
[18] D. S. Shapovalov and E. T. Vishniac, The Astrophysical Journal 728, 153 (2011).
Astrophys.J. 738, 66 (2011).
[19] S. M. Tobias and F. Cattaneo,
Nature(London) 497, 463 (2013).
[20] S. Sridhar and N. K. Singh, MNRAS445, 3770 (2014),
1.0
arXiv:1306.2495.
1x2x1
[21] J. Squire and A. Bhattacharjee,
1x2x4
Physical Review Letters 115, 175003 (2015), 4x8x4
arXiv:1506.04109 [astro-ph.SR]. 0.5
[22] J. Squire and A. Bhattachar-
jee, Astrophys.J. 813, 52 (2015),
arXiv:1507.03154 [astro-ph.HE]. ⋅vb] 0.0
[23] http://ipag.osug.fr/~lesurg/snoopy.html. [
[24] G. Lesur and P.-Y. Longaretti, A&A528, A17 (2011),
arXiv:1012.2690 [astro-ph.EP]. −0.5
[25] We do not conduct a detailed study on the critical am-
plitude required to trigger and sustain turbulence and
simplyusevelocityperturbationsofamplitudeLS forall
−1.0
of our simulations. 0 200 400 600 800 1000
[26] B. Hof, J. Westerweel, T. M. Schneider, and B. Eck- Time (units of S−1)
hardt,Nature (London) 443, 59 (2006).
FIG. 8: Same as fig. 6 but for [v·b] volume averaged.
7
FIG. 9: xy averagedhV i and hV i for the runs 1×2×4 (left column) and 4×8×4 (right column) for the first 1000
x y
shear times. Compare fig. 4 from the main paper.