Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets
An interface-tracking numerical algorithm for the simulation of magnetohydrodynamic multiphase / free surface flows in the low-magnetic-Reynolds-number approximation of (Samulyak R., Du J., Glimm J., Xu Z., J. Comp. Phys., 2007, 226, 1532) is described. The algorithm has been implemented in multi-ph...
Збережено в:
| Дата: | 2010 |
|---|---|
| Автори: | , , , , |
| Формат: | Стаття |
| Мова: | English |
| Опубліковано: |
Інститут фізики конденсованих систем НАН України
2010
|
| Назва видання: | Condensed Matter Physics |
| Онлайн доступ: | http://dspace.nbuv.gov.ua/handle/123456789/32126 |
| Теги: |
Додати тег
Немає тегів, Будьте першим, хто поставить тег для цього запису!
|
| Назва журналу: | Digital Library of Periodicals of National Academy of Sciences of Ukraine |
| Цитувати: | Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets / R.V. Samulyak, W. Bo, X. Li, H. Kirk, K. McDonald // Condensed Matter Physics. — 2010. — Т. 13, № 4. — С. 43402:1-12. — Бібліогр.: 32 назв. — англ. |
Репозитарії
Digital Library of Periodicals of National Academy of Sciences of Ukraine| id |
irk-123456789-32126 |
|---|---|
| record_format |
dspace |
| spelling |
irk-123456789-321262012-04-10T12:16:26Z Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets Samulyak, R.V. Bo, W. Li, X. Kirk, H. McDonald, K. An interface-tracking numerical algorithm for the simulation of magnetohydrodynamic multiphase / free surface flows in the low-magnetic-Reynolds-number approximation of (Samulyak R., Du J., Glimm J., Xu Z., J. Comp. Phys., 2007, 226, 1532) is described. The algorithm has been implemented in multi-physics code FronTier and used for the simulation of MHD processes in liquids and weakly ionized plasmas. In this paper, numerical simulations of a liquid mercury jet entering strong and nonuniform magnetic field and interacting with a powerful proton pulse have been performed and compared with experiments. Such a mercury jet is a prototype of the proposed Muon Collider / Neutrino Factory, a future particle accelerator. Simulations demonstrate the elliptic distortion of the mercury jet as it enters the magnetic solenoid at a small angle to the magnetic axis, jet-surface instabilities (filamentation) induced by the interaction with proton pulses, and the stabilizing effect of the magnetic field. Описано інтерфейс-простежувальний числовий алгоритм для моделювання магнітогідродинамічних (MHD) потоків мультифаза/вільна поверхня у наближенні низькомагнітних чисел Рейнольдса (Samulyak R., Du J., Glimm J., Xu Z., J. Comp. Phys., 2007, 226, 1532). Алгоритм застосовується у мультифізичному коді FronTier і використовується для моделювання MHD процесів у рідині та слабко іонізованій плазмі. Числово змодельовано рідкий ртутний струмень, що входить у сильне та неоднорідне магнітне поле і взаємодіє з потужним протонним імпульсом, і порівняно з експериментом. Такий ртутний струмінь є прототипом запропонованого майбутнього прискорювача частинок Мюонний коллайдер/Нейтринна фабрика. У процесі моделювання продемонстровано еліптичну дисторцію ртутного струменя, коли він входить у магнітне поле за малого кута до магнітної осі, нестійкості струмінь-поверхня (філаментацію), індуковані взаємодією з протонними імпульсами, та стабілізувальний ефект магнітного поля. 2010 Article Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets / R.V. Samulyak, W. Bo, X. Li, H. Kirk, K. McDonald // Condensed Matter Physics. — 2010. — Т. 13, № 4. — С. 43402:1-12. — Бібліогр.: 32 назв. — англ. 1607-324X PACS: 47.11, 47.35, 47.55 http://dspace.nbuv.gov.ua/handle/123456789/32126 en Condensed Matter Physics Інститут фізики конденсованих систем НАН України |
| institution |
Digital Library of Periodicals of National Academy of Sciences of Ukraine |
| collection |
DSpace DC |
| language |
English |
| description |
An interface-tracking numerical algorithm for the simulation of magnetohydrodynamic multiphase / free surface flows in the low-magnetic-Reynolds-number approximation of (Samulyak R., Du J., Glimm J., Xu Z., J. Comp. Phys., 2007, 226, 1532) is described. The algorithm has been implemented in multi-physics code FronTier and used for the simulation of MHD processes in liquids and weakly ionized plasmas. In this paper, numerical simulations of a liquid mercury jet entering strong and nonuniform magnetic field and interacting with a powerful proton pulse have been performed and compared with experiments. Such a mercury jet is a prototype of the proposed Muon Collider / Neutrino Factory, a future particle accelerator. Simulations demonstrate the elliptic distortion of the mercury jet as it enters the magnetic solenoid at a small angle to the magnetic axis, jet-surface instabilities (filamentation) induced by the interaction with proton pulses, and the stabilizing effect of the magnetic field. |
| format |
Article |
| author |
Samulyak, R.V. Bo, W. Li, X. Kirk, H. McDonald, K. |
| spellingShingle |
Samulyak, R.V. Bo, W. Li, X. Kirk, H. McDonald, K. Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets Condensed Matter Physics |
| author_facet |
Samulyak, R.V. Bo, W. Li, X. Kirk, H. McDonald, K. |
| author_sort |
Samulyak, R.V. |
| title |
Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets |
| title_short |
Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets |
| title_full |
Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets |
| title_fullStr |
Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets |
| title_full_unstemmed |
Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets |
| title_sort |
computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets |
| publisher |
Інститут фізики конденсованих систем НАН України |
| publishDate |
2010 |
| url |
http://dspace.nbuv.gov.ua/handle/123456789/32126 |
| citation_txt |
Computational algorithms for multiphase magnetohydrodynamics and applications to accelerator targets / R.V. Samulyak, W. Bo, X. Li, H. Kirk, K. McDonald // Condensed Matter Physics. — 2010. — Т. 13, № 4. — С. 43402:1-12. — Бібліогр.: 32 назв. — англ. |
| series |
Condensed Matter Physics |
| work_keys_str_mv |
AT samulyakrv computationalalgorithmsformultiphasemagnetohydrodynamicsandapplicationstoacceleratortargets AT bow computationalalgorithmsformultiphasemagnetohydrodynamicsandapplicationstoacceleratortargets AT lix computationalalgorithmsformultiphasemagnetohydrodynamicsandapplicationstoacceleratortargets AT kirkh computationalalgorithmsformultiphasemagnetohydrodynamicsandapplicationstoacceleratortargets AT mcdonaldk computationalalgorithmsformultiphasemagnetohydrodynamicsandapplicationstoacceleratortargets |
| first_indexed |
2025-07-03T12:39:37Z |
| last_indexed |
2025-07-03T12:39:37Z |
| _version_ |
1836629494512746496 |
| fulltext |
Condensed Matter Physics 2010, Vol. 13, No 4, 43402: 1–12
http://www.icmp.lviv.ua/journal
Computational algorithms for multiphase
magnetohydrodynamics and applications to accelerator
targets
R.V. Samulyak1,2, W. Bo2, X. Li2, H. Kirk3, K. McDonald4
1 Department of Applied Mathematics and Statistics, Stony Brook University, Stony Brook, NY 11794, USA
2 Computational Science Center, Brookhaven National Laboratory, Upton, NY 11973, USA
3 Physics Deparment, Brookhaven National Laboratory, Upton, NY 11973, USA
4 Joseph Henry Laboratories, Princeton University, Princeton, NJ 08544, USA
Received June 19, 2010
An interface-tracking numerical algorithm for the simulation of magnetohydrodynamic multiphase / free surface
flows in the low-magnetic-Reynolds-number approximation of (Samulyak R., Du J., Glimm J., Xu Z., J. Comp.
Phys., 2007, 226, 1532) is described. The algorithm has been implemented in multi-physics code FronTier
and used for the simulation of MHD processes in liquids and weakly ionized plasmas. In this paper, numerical
simulations of a liquid mercury jet entering strong and nonuniform magnetic field and interacting with a pow-
erful proton pulse have been performed and compared with experiments. Such a mercury jet is a prototype
of the proposed Muon Collider / Neutrino Factory, a future particle accelerator. Simulations demonstrate the
elliptic distortion of the mercury jet as it enters the magnetic solenoid at a small angle to the magnetic axis,
jet-surface instabilities (filamentation) induced by the interaction with proton pulses, and the stabilizing effect
of the magnetic field.
Key words: MHD algorithm, multiphase flow, front tracking, mercury target
PACS: 47.11, 47.35, 47.55
1. Introduction
The main driver for computational magnetohydrodynamics is the research in magnetically con-
fined nuclear fusion which has recently been boosted by the International Thermonuclear Exper-
imental Reactor project (ITER, [1]). Numerical algorithms for nuclear-fusion-simulation research
are optimized for highly conductive, fully ionized plasmas. In simplified studies, even the infinite-
conductivity approximation is widely used (ideal MHD approximation [2]). But even in thermonu-
clear fusion devices such as tokamaks, low-conductivity, weakly ionized plasma may still be present
under special conditions despite very high temperature (the nuclear fusion ignition temperature is
of the order of 10 keV or 108 K). In tokamaks, weakly ionized plasma is found in the ablated clouds
of cryogenic fuel pellets. The injection of such frozen deuterium-tritium pellets is considered the
most efficient technique for the fueling of tokamaks [3]. Numerical algorithms designed for highly
conductive plasma have very limited applicability to weakly ionized gases. Another large area of
application of low-conductity MHD span liquid metals and liquid salts. Liquid metals, such as
lithium, can be used in tokamaks as a plasma-facing component protecting tokamak walls. In this
paper, we are exploring another application area of algorithms for liquid-metal MHD: the study of
liquid-mercury-jet targets for future particle accelerators.
In many applications, liquid metals or weakly ionized plasmas have either free surfaces or mov-
ing, geometrically complex interfaces with either non-conductive or highly conductive media. The
presence of interfaces imposes a major challenge on the development of high-quality numerical al-
gorithms. The majority of numerical studies of free-surface MHD flows are based on semi-analytical
treatment of simplified flow regimes [4, 5]. Simplified models have successfully been used for the
description of self-organized filaments in dielectric barrier-glow discharges [6, 7] and the numerical
c© R.V. Samulyak, W. Bo, X. Li, H. Kirk, K. McDonald 43402-1
http://www.icmp.lviv.ua/journal
R. Samulyak et al.
modeling of micro-plasma instabilities [8]. But analytical models have a limited applicability for
complex systems involving strongly coupled multi-physics phenomena.
The numerical resolution of interfaces can be performed using various methods developed for
multiphase fluid dynamics. In our recent work [10], a free-surface MHD algorithm was developed
based on explicit tracking of material interfaces. To the best of our knowledge, the only other
fully numerical treatment of general, free-surface, incompressible liquid flows is implemented in the
HIMAG code [11] using the level-set algorithm for fluid interfaces, the electric-potential formulation
for electromagnetic forces, and the incompressible-fluid-flow approximation. The key feature of our
algorithm is the use of the method of front tracking [12] for the propagation of vapor interfaces. The
FronTier code is capable of tracking and resolving topological changes of large numbers of interfaces
in two- and three-dimensional spaces [13]. In the method of front tracking, the interface is a
Lagrangian mesh moving through a volume-filling rectangular mesh according to the solution of the
corresponding Riemann problem. High-resolution, shock-capturing, Godunov-type solvers are used
to update hyperbolic states in the interior away from interfaces. The explicit treatment of interfaces
typical of the method of front tracking greatly reduces the numerical diffusion and is especially
advantageous for multi-physics problems involving phase transitions. It allows not only to solve
accurately the Riemann problem for the phase boundary, but also to apply different mathematical
approximations in the regions separated by interfaces to account for different material properties
and, if necessary, eliminate fast time scales in numerical simulations. The three-dimensional, front-
tracking algorithm has been recently upgraded with a so-called locally-grid-based tracking [14]
that enabled robust resolution of very complex topological changes in the FronTier code. The
comparison of front tracking and other numerical methods resolving interfaces, in particular the
level-set method, can be found in [13]. The main features of front tracking such as the conservative
property and the ability to resolve complex interfaces are also present in particle methods such
as the smoothed particle hydrodynamics [9]. However the particle-based codes generally require
longer computational time compared to grid-based codes.
Since most of our applications involve strong hydrodynamic waves, we solve the system of
compressible Navier-Stokes equations coupled to electromagnetic forces. A numerical algorithm
for the reduced MHD system by using the incompressible hydrodynamic approximations will be
described in a forthcoming paper.
The front-tracking MHD algorithm has been used for the simulation of tokamak fueling by
the ablation of cryogenic pellets [15, 16] and the dynamics of laser-ablation plumes [17]. In this
paper, we apply the numerical MHD algorithms for the simulation of liquid mercury jet targets for
the Neutrino Factory/Muon Collider. In particle accelerators, liquid metal targets convert intense
proton beams into neutrons, pions or other particles used in fundamental and applied studies. The
aim of the targetry group of the international, multi-institutional Neutrino Factory/Muon Collider
Collaboration [18] is to explore the feasibility of high power targets for future particle accelerators
and, in particular, for the proposed Neutrino Factory and Muon Collider. The target will contain
a series of mercury jet pulses of about 1 cm in diameter and 30 cm in length. Each pulse will
be shot at a velocity of 20 m/s into a 15-20-T magnetic field at a small angle to the axis of the
magnetic field. When the jet reaches the center of the magnet, it interacts with a series of proton
pulses depositing about 100 J/g of energy in the mercury. The aim of numerical simulations is to
describe the hydrodynamic response of the target interacting with proton pulses in magnetic fields,
the understanding of which is of major importance for reliable target design.
We have already performed mathematical modeling, software development and simulations of
liquid mercury jet targets interacting with high-power proton beams in magnetic fields [19, 20]
and made predictions for the targetry experiment at CERN, so-called MERIT. The simulation
predicted strong instabilities and cavitation of the mercury jet interacting with proton pulses at
zero magnetic field and their reduction / stabilizing by a longitudinal magnetic field. Simulation
predictions were qualitatively confirmed by the MERIT experiment conducted in CERN in the
fall of 2007 [21]. In this paper, we use exact input parameters of the MERIT experiment and
quantitatively compare the simulations with experimental data.
The paper is organized as follows. The main governing equations are presented in section 2.
43402-2
Algorithms for multiphase magnetohydrodynamics and applications
In section 3, we describe the front-tracking numerical algorithms for free-surface MHD and the
algorithm for fluid cavitation. Numerical simulations of the mercury-jet target for the Neutrino
Factory / Muon Collider are described in section 4. We conclude this paper with the summary of
our results and perspectives for the future work.
2. Governing equations
We are interested in the description of multiphase or multi-material systems involving con-
ducting fluids or weakly ionized gases interacting with neutral fluids or gases in the presence of
magnetic fields. Interfaces of the phase or material separation are assumed to be sharp (the thick-
ness of the interface is negligible) and, in general, geometrically complex. The numerical simulation
of liquid metals and liquid salts as well as weakly ionized plasmas, which are relatively weak elec-
trical conductors, is difficult using the standard full systems of MHD equations [2]. Fast diffusion of
the magnetic field, caused by low value of electrical conductivity, introduces unwanted small time
scales into the problem. If the time scale of the diffusion of the magnetic field is small compared
to hydrodynamic time scale, the magnetic Reynolds number [4]
ReM =
4πuσL
c2
,
where L is the typical length scale, u is the fluid velocity, and σ is the electric conductivity, is
small. If in addition the eddy-current-induced magnetic field δB is small compared to the external
field B, the full system of MHD equations can be simplified by neglecting the time evolution
of the magnetic field. In this case, the generalized Ohm’s law is used for the calculation of the
current-density distribution instead of the Maxwell equation J = c
4π∇ ×H, where the magnetic
field H and the magnetic induction B are related by the magnetic permeability coefficient µ:
B = µH. The resulting system of MHD equations for compressible conductive fluids or gases in
the low-magnetic-Reynolds-number approximation is
∂ρ
∂t
= −∇ · (ρu), (1)
ρ
(
∂
∂t
+ u · ∇
)
u = −∇P +
1
c
(J×B), (2)
ρ
(
∂
∂t
+ u · ∇
)
e = −P∇ · u+
1
σ
J
2, (3)
P = P (ρ, e), (4)
J = σ
(
−∇ϕ+
1
c
u×B
)
, (5)
∇ · (σ∇ϕ) =
1
c
∇ · σ(u×B), (6)
where ρ, u, P , and e are the density, velocity, pressure, and specific internal energy of the fluid.
Equation (5) is the consequence of Ohm’s law, and the Poisson equation (6) follows from the local
neutrality of fluids ∇ · J = 0 and the Ohm’s law (5).
Equation (4) is the equation of state (EOS) that describes material properties. Three material
substances are present in our simulations: liquid mercury, ambient gas, and mercury vapor inside
cavitation bubbles. The polytropic ideal gas equation of state is used for the ambient gas and
mercury vapor, and the stiffened polytropic EOS [23]
P = (γl − 1)ρ(E + E∞)− γlP∞
with the adiabatic exponent γl = 3.2 and the stiffening constant P∞ = 8 · 1010g/(cm s
2
) is used
for liquid mercury in this study. The stiffened polytropic EOS extends the range of pressure to
negative values in order to account for the transient effect of tension in liquids. With the stiffened
43402-3
R. Samulyak et al.
polytropic EOS, a thermodynamically consistent state can be achieved for P > −P∞ as the sound
speed becomes imaginary below this point: c2 = γ(P −P∞)/ρ. The cavitation algorithm described
in the next section prevents the pressure in liquids from falling below some specified critical pressure
at which the liquid breaks in the form of cavitation bubbles that expand and relieve the tension.
The following boundary conditions must be satisfied at the interface Γ of a conducting fluid
with a dielectric medium:
i) the normal component of the velocity field is continuous across the interface;
ii) the pressure jump at the interface is defined by the surface tension S and main radii of curvature:
∆P |Γ = S
(
1
r1
+
1
r2
)
; (7)
iii) the normal component of the current density vanishes at the interface giving rise to the Neumann
boundary condition for the electric potential
∂ϕ
∂n
∣
∣
∣
∣
Γ
=
1
c
(u×B) · n, (8)
where n is a normal vector at the fluid free surface Γ.
3. Numerical algorithms and implementation
3.1. Front-tracking algorithm for MHD equations
The existence of sharp, geometrically complex material interfaces presents a great challenge to
numerical algorithms. We solve the coupled hyperbolic-elliptic system (1–8) in an operator-splitting
fashion using the method of front tracking. The front-tracking technique, its implementation in
the multiphysics code FronTier, and applications to computational fluid dynamics problems are
described in detail in [12–14] and references therein. The front-tracking algorithm for multiphase
MHD equations in the low magnetic Reynolds number approximation was developed in [10]. In this
paper, we summarize the main ideas of the algorithm and refer the reader to the above references
for details.
Front tracking is an adaptive computational method in which a lower dimensional moving
grid or interface (a triangulated surface in 3D space) represents material interfaces. The interface
contains left and right physical states corresponding to materials on both sides and keeps the
discontinuity sharp. The key feature of the method is greatly reduced numerical diffusion due to
the absence of numerical finite differentiation across the interface. Front tracking is implemented
in a multi-physics code FronTier [13]. FronTier is capable of tracking and resolving the topological
changes of a large number of interfaces in two- and three-dimensional spaces.
Using the operator splitting, we solve the hyperbolic equations first, followed by solving the
elliptic boundary value problem (6) and (8). The time step starts with the advance of the interface.
The propagation of the interface is accomplished by solving the MHD equations in the normal and
tangential directions to each interface point. In the low-magnetic-Reynolds-number approximation,
the MHD equations reduce to Euler equations with external forces that come from the elliptic step.
The rectangular grid, interface, and states for the normal propagation of the interface are shown
schematically in figure 1. After the propagation of the interface points, the new interface is checked
for consistency of intersections. The untangling of the interface at this stage consists in removing
the unphysical intersections, and rebuilding a topologically correct interface. The final phase of
the hyperbolic time step update consists of computing new states on the rectangular spatial grid
using one of the second order shock capturing solvers.
After the hyperbolic step, the elliptic problem for the electric potential (6, 8) is solved using the
embedded-boundary method [24] (see [25] for the latest development of the embedded-boundary
method for 3D moving elliptic-interface problems with front tracking). High-performance, parallel-
software libraries of preconditioners and iterative solvers based on Krylov subspace methods such
43402-4
Algorithms for multiphase magnetohydrodynamics and applications
Figure 1. Rectangular grid, interface, and states for the method of front tracking. States contain
density, momentum, and energy density of the fluid, and references to the EOS model and other
parameters.
as PETSC [26] are used for solving the corresponding linear system of equations, and the electro-
magnetic terms are calculated. In the final stage of the time step, the states in the interior and on
interfaces are modified using electromagnetic terms.
3.2. Models for phase transition and cavitation
We are interested in the simulation of liquids in the presence of strong rarefaction waves. If the
liquid is not capable of sustaining tension, it breaks or cavitates [27]. The capability of liquids to
withstand tension depends on their purity. The homogeneous nucleation theory [27], applicable to
liquids without impurities that serve as initial cavitation centers, gives the following probability of
nucleation in a volume V during time interval τ :
Σ = 1− exp(−J0V t exp(−EC/(kbT ))) , (9)
where
EC =
16πS3
3△P 2
C
is the critical energy necessary to create the surface of a cavitation bubble of the critical radius RC,
and is related to the critical pressure (the critical strength of the tensile pressure in the liquid) as
Rc =
2S
∆Pc
,
and S is the surface tension. The coefficient J0 is
J0 = N
(
2S
πm
)1/2
,
where N is the number density of the liquid (molecules/volume) and m is the molecular mass. In
real fluids, the existence of impurities significantly reduces the values of the critical pressure and
energy compared to the values predicted by the homogeneous nucleation theory. Our DNS model of
cavitation is implemented as follows: when the pressure falls below the critical pressure, cavitation
bubbles are formed in some spatial distribution in the rarefaction wave by inserting numerically
tracked bubble surfaces and replacing liquid states with vapor states inside. The values of the
critical pressure and the initial number density of cavitation bubbles are selected individually for
each numerical experiment by estimating liquid conditions from available experimental data. In
most simulations, we assume that the amount of vapor in a bubble is constant and the expansion
43402-5
R. Samulyak et al.
of cavitation bubbles is caused only by pressure gradients. In other words, we neglect the phase
transition on the bubble surface. If the phase-transition-induced mass flux is essential for the
interface dynamics, the phase-transition algorithm developed in [28] can be used. Fast vaporization
or ablation is especially important in the pellet-fueling processes mentioned in the introduction.
Initial cavitation bubble sizes in real liquids (for instance, Rc = 1 micron for mercury at ∆Pc = 10
bar) is close to numerically resolved limits. However in practice, we frequently use larger initial
bubble sizes for coarser grid computations in large domains. This leads to some numerical errors
including the loss of mass and momentum conservation. While we estimate that such numerical
errors are insignificant in most of simulations, we are currently working on strict enforcement of the
conservative properties using a multi-scale coupling with a stand-alone, resolved, one-dimensional
simulation of a cavitation bubble in the ambient liquid. Similar cavitation algorithms have already
been used in 2D [20] and 3D simulations [29].
4. Applications to accelerator targets
In this section, we apply the front-tracking algorithms for multiphase MHD to the simulation of
the mercury target for the proposed Neutrino Factory / Muon Collider [18]. The target will contain
a series of 30-cm-long and 1-cm-diameter mercury jets entering a strong (∼ 15 Tesla) magnetic field
at a small angle to the solenoid axis. When each jet reaches the center of the solenoid, it interacts
with a powerful proton pulse penetrating the jet and depositing energy of the order of 100 J/g into
mercury. The purpose of our numerical simulations is to evaluate the states of the target before
and after the interaction with protons to optimize the target design. Preliminary simulations have
been reported in [19, 20, 22].
Figure 2. Schematic of the mercury jet target for Neutrino Factory / Muon Collider.
4.1. Simulation of the mercury jet entering a solenoid magnet
The nonuniform transverse component of the magnetic field with respect to the jet trajectory,
caused by a small angle between the jet and the magnetic solenoid axis, distorts the jet during the
motion toward the solenoid center. To evaluate the state of the jet target before the interaction
with the solenoid axis, we performed numerical simulations of the jet entering the solenoid using
the real profile of the magnetic field along the jet trajectory in the MERIT experiment. Results
are summarized in figure 3, which plots the transverse-distortion ratio of the jet or the maximum
radius in the transverse cross section normalized by the unperturbed initial radius of the jet. We
observe that the distortion strongly depends on the angle between the jet and the solenoid axis:
the maximum distortion ranges from 1.2 at the angle of 0.05 rad to 2.75 at the angle of 0.15 rad.
The latter distortion is unacceptable for the target: it transforms the jet into a thin sheet and
significantly reduces the effective cross section with the proton pulse and, as a result, the pion-
production rate. To reduce the amount of distortion, the angle of 0.033 rad was used in the MERIT
43402-6
Algorithms for multiphase magnetohydrodynamics and applications
experiment. Simulations of mercury jets in transverse magnetic fields similar to those presented in
figure 3 agree very well with theoretical calculations [10], other experiments [5], and simulations
using the HYMAG code [30]. However, the simulations underestimated the amount of distortion
for the MERIT experiment. Numerical simulation of the mercury jet entering the 15-T solenoid
at the angle of 0.033 rad predicted the distortion of 1.17 and the value of 1.8 was observed [21].
The reason for this disagreement is related to deviations from ideal simulation conditions which
are currently under investigation.
Figure 3. Normalized elliptic deformation of the mercury jet entering a 15-T solenoid at the
angle of 0.15 rad (top line), 0.1 rad (middle line), and 0.05 rad (bottom line) to the magnetic
axis.
4.2. Simulation of the mercury jet interacting with proton pulses.
When the mercury jet reaches the solenoid center, it interacts with a proton pulse, depositing
energy of the order of 100 J/g into mercury. Because of the short time scale of the interaction,
[]
(b)
Figure 4. Pressure distribution in mercury after the interaction with proton pulses. Left image:
cylindrical jet with a 24-GeV, 10-teraproton pluse, right image: elliptic jet with a 14-GeV, 10-
teraproton pluse.
43402-7
R. Samulyak et al.
we assume that the increase of the internal energy and pressure is an isochoric process. In com-
putations, we modify the internal energy and pressure states of the jet during the initial time
step according to predictions of atomistic Monte-Carlo simulations of the mercury-proton pulse
interaction using the MARS code [31]. In most of simulations, we assumed that the initial shape
of the jet is distorted by the transverse component of the magnetic field, and the cross section
of the jet is an ellipse with the long and short radii of 0.8 and 0.4 cm, correspondingly. But in
some simulations we also used the unperturbed cylindrical jet. The profiles of pressure after the
interaction with 24-GeV, 10-teraproton and 14-GeV, 10-teraproton pulses are shown in figure 4.
For the mercury jet with the elliptic cross section, we monitored the jet surface velocity in four
radial directions, or points A, B, C and D shown in figure 4(b).
Figure 5. Cavitation bubbles in the mercury jet at 20, 130, 200, and 250 µs after the interaction
of the cylindrical jet with the proton pulse.
After the energy deposition, the high-pressure wave propagates outward and reflects from the
mercury-air interface as a strong rarefaction wave. Since the mercury is incapable of sustaining
such a large tension, it cavitates and the growth of cavitation bubbles causes a reduction of tension,
rapid jet expansion and surface instabilities. A snapshot of cavitation bubbles in the jet is shown
in figure 5.
Figure 6. Filaments on the surface of initially
cylindrical jet at 150 µs after the interaction
with the proton pulse in magnetic fields rang-
ing from 0 to 15 Tesla.
If the jet cavitation and expansion occurs
in a longitudinal magnetic field, the radial mo-
tion of the fluid induces vortices of azimuthal
current. Then, the Lorentz force opposes the
fluid motion. Hence the magnetic field reduces
the amount of cavitation and tends to stabilize
the mercury jet. Snapshots of the jet surfaces
at 100 µs after the interaction with the proton
pulse in magnetic fields ranging from 0 to 15
Tesla are shown in figure 6.
Figures 7–10 show the evolution of veloci-
ties of elliptic-jet surface filaments in four radial
directions A, B, C, and D, as explained above.
The maximum velocity of filaments ejected in
the direction of the short axis reaches 35 m/s.
In all directions, the velocity decreases with
the increase of the magnetic field strength. We
would like to emphasize that the formation
and evolution of filaments cannot be attributed
solely to classical fluid interface instabilities
such as the Rayleigh-Taylor and Richtmyer-
43402-8
Algorithms for multiphase magnetohydrodynamics and applications
Meshkov instabilities. We have shown that the formation of mercury jet filaments critically depends
on the presence of cavitation.
Figure 7. Expansion velocity of jet surface fi-
laments at zero magnetic field.
Figure 8. Expansion velocity of jet surface fi-
laments in 5 Tesla magnetic field.
Figure 9. Expansion velocity of jet surface fi-
laments in 10 Tesla magnetic field.
Figure 10. Expansion velocity of jet surface
filaments in 15 Tesla magnetic field.
In figure 11, we compare the velocity of jet-surface filaments with experimental data, which is
sparse; it is available for different values of the beam intensity but not always for the 10-teraproton
beam which was used for numerical simulations. The simulated filament velocity at 5-T magnetic
field agrees well with the experiment but the experimental value must be assigned a large error
bar because of inconsistency of measurements for 10-, 15- and 20-teraproton beams. The simulated
velocity at 10-T field agrees well with the interpolated value using the experimental data points at
5 and 15 teraproton. There is only one measurement available for the 15-T magnetic field which
shows no velocity reduction compared to the 10-T field. The corresponding simulation indicates a
further reduction of the filament velocity compared to 5- and 10-T fields. Because of the sparseness
of experimental data and limited number of spikes used in the data analysis, it is difficult to make
a conclusion on the last case.
43402-9
R. Samulyak et al.
Figure 11. Comparison of experimental and simulated values of the jet surface filament velocity.
5. Conclusions and future work
We have described a numerical algorithm for the simulation of free-surface magnetohydrody-
namic flows at low magnetic Reynolds numbers. The corresponding governing equations constitute
a coupled hyperbolic-elliptic system in a geometrically complex and evolving domain. The numer-
ical algorithm includes the interface-tracking technique for the hyperbolic problem, a Riemann
problem for the material interface, discretization of elliptic equations in irregular domains with in-
terface constraints using the embedded-boundary method, and high-performance, parallel solvers
such as MUSCL-type schemes for hyperbolic problems and iterative solvers implemented in the
PETSc package. An extensive theoretical analysis of the method of front tracking for hyperbolic
systems of conservation laws has already been performed in earlier works, and the method has been
validated and tested on problems of Rayleigh-Taylor and Richtmyer-Meshkov surface instabilities.
The embedded-boundary technique for irregular geometric domains within the method of front
tracking was developed in [10] and generalized to include elliptic problems with discontinuities
along interfaces in [25], and shown to be second order accurate for both the electric potential and
its gradient.
The front-tracking MHD algorithm has been used for the simulation of tokamak fueling by the
ablation of cryogenic pellets [15, 16] and the dynamics of laser ablation plumes [17]. In this paper,
we applied the numerical MHD algorithms for the simulation of liquid mercury jet targets for
the Neutrino Factory/Muon Collider. Simulations of the mercury jet entrance into a non-uniform
magnetic field at small angles to the solenoid axis showed that the transverse component of the
magnetic field distorts the jet. The jet cross section becomes an elongated ellipse with the long
axis parallel to the direction of the transverse field. Simulations of the mercury jet interacting
with proton pulses in magnetic fields have also been performed. After the proton-pulse energy
deposition, the pressure in the center ranges from 10 to 30 kbar, depending on the proton-pulse
properties. Then the high-pressure wave propagates towards the jet boundary and causes tension
that leads to severe cavitation in mercury and the formation and growth of surface filaments. The
fastest filaments reach velocities of 35 m/s. The cavitation and surface filamentation is reduced by
longitudinal magnetic fields. Simulations are in reasonably good agreement with the results of the
mercury target experiment MERIT [21] conducted at CERN in the fall of 2007.
In our future work, we will develop a front-tracking MHD algorithm for multiphase or free-
surface incompressible fluids in the low-magnetic-Reynolds-number approximation. The algorithm
will enable the FronTier code to perform long-time-scale simulations of liquid metals in magnetic
fields for applications which do not involve strong waves. Such are liquid-metal blankets on the
interface with plasma that can be used for the protection of tokamak walls, mercury-handling
components in accelerator targets and other applications. We are also working on the free-surface
algorithm for the opposite MHD regime: the ideal MHD approximation for highly conducting
43402-10
Algorithms for multiphase magnetohydrodynamics and applications
plasmas. Such an algorithm will allow us to simulate the behavior of magnetized plasma targets
compressed by liners with the purpose of achieving nuclear fusion ignition. Hybrid magneto-inertial
fusion methods [32] have recently gained large attention due to potential capability of overcoming
difficulties of traditional approaches.
Acknowledgments
This manuscript has been authored in part by Brookhaven Science Associates, LLC, under
Contract No. DE–AC02–98CH10886 with the U.S. Department of Energy. The United States Gov-
ernment retains, and the publisher, by accepting the article for publication, acknowledges, a world-
wide license to publish or reproduce the published form of this manuscript, or allow others to do
so, for the United States Government purpose. This research utilized resources at the New York
Center for Computational Sciences at Stony Brook University/Brookhaven National Laboratory
which is supported by the U.S. Department of Energy under Contract No. DE–AC02–98CH10886
and by the State of New York.
References
1. ITER: International Thermonuclear Experimental Reactor, http://www.iter.org/
2. Freidberg J., Ideal Magnetohydrodynamics (Modern Perspectives in Energy). Springer, 1987.
3. Pegourie B., Plasma Phys. Control. Fusion, 2007, 49, R87.
4. Molokov S., Reed C.B., Review of free-surface MHD experiments and modeling, Technical Report
ANL/TD/TM99–08, Argonne National Laboratory, 1999.
5. Oshima S., Yamane R., Moshimaru Y., Matsuoka T., JCME Int. J., 1987, 30, 437.
6. Muller I., Punset C., Ammelt E., Purwins H.-G., Boeuf J.P., IEEE Trans. Plasma Sci., 1999, 27, 20.
7. Brauer I., Punset C., Purwins H.-G., Boeuf J.P., J. Appl. Phys., 1999, 85, 7569.
8. Datsko B., Semiconductors, 1997, 31, 146.
9. Hoover W., Smooth particle applied mechanics: the state of the art. World Scientific Publishing Co.,
Inc., 2006.
10. Samulyak R., Du J., Glimm J., Xu Z., J. Comp. Phys., 2007, 226, 1532.
11. Morley N., Smolentsev S., Munipalli R., Ni M.-J., Gao D., Abdou M., Fusion Eng. Design, 2004, 72,
3.
12. Glimm J., Grove J., Li X.L., Tan D.C., SIAM J. Sci. Comp., 2000, 21, 2240.
13. Du J., Fix B., GlimmJ., Jia X., Li X., Li Y., Wu L., J. Comp. Phys., 2006, 213, 613.
14. Fix B., Glimm J., Li X., Li Y., Liu X., Samulyak R., Xu Z., J. Phys.: Conf. Series, 2005, 16, 471.
15. Samulyak R., Lu T., Parks P., Nuclear Fusion, 2007, 47, 103.
16. Parks P., Lu T., Samulyak R., Physics of Plasmas, 2009, 16, 060705.
17. Du J., Lu T., Samulyak R., J. Nanoscience and Nanotechnology, 2008, 8, 3674.
18. The Neutrino Factory and Muon Collider collaboration, http://www.cap.bnl.gov/mumu/
19. Samulyak R., Prykarpatskyy Y., Mathematics and Computers in Simulations, 2004, 65, 431.
20. Samulyak R., Lu T., Prykarpatskyy Y., Glimm J., Xu Z., Kim M.N., Int. J. Multiscale Comp. Eng.,
2006, 4, 377.
21. McDonald K. et al., The MERIT high-power target experiment at the CERN PS, Proc. PAC09 (Van-
couver, Canada), TU4GRI03, 2009.
22. Samulyak R., Glimm J., Oh W., Kirk H., McDonald K., Lecture Notes in Comp. Sci., 2003, 2667,
558.
23. Menikoff R., Plohr B., Rev. Mod. Phys., 1989, 61, 75.
24. Johansen H., Colella P., J. Comput. Phys., 1998, 147, 60.
25. Wang S., Samulyak R., Guo T., Acta Mathematica Scientia, 2010, 30B, No. 2, 499.
26. Balay S., Gropp W., McInnes L., Smith B., PETSc Users Manual, ANL–95/11 – Revision 2.1.1,
Argonne National Laboratory, 2001.
27. Brennen C.E., Cavitation and Bubble dynamics. Oxford University Press, 1995.
28. Xu Z. L., Lu T., Samulyak R., Glimm J., Ji X.M., SIAM J. Sci. Computing, 2008, 30, No. 2, 895.
29. Bo W., Liu X., Glimm J., Li X., ASME Journal of Fluids Engineering, 2010 (submitted).
30. Morley N., private communication, 2005.
43402-11
R. Samulyak et al.
31. Striganov S., MARS15 simulations of the MERIT Mercury Target Experiment, Third High-Power
Targetry Workshop, Bad Zurzach, Switzerland, September 13, 2007.
32. Thio Y.C.F., Panarella E., Kirkpatrick R.C., Knupp C.E., Wysocki F., Parks P., Schmidt G., Magne-
tized target fusion in a spheroidal geometry with standoff drivers, in Current Trends in International
Fusion Research II, edited by Panarella E., National Research Council Canada, Ottawa, Canada, 1999
Розрахунковi алгоритми для багатофазної
магнiтогiдродинамiки i застосування до мiшеней
прискорювачiв
Р.В. Самуляк1,2, В. Бо2, Х.Лi 2, Г. Кiрк 3, K. К. МакДональд4
1 Факультет прикладної математики i статистики, Унiверситет Стонi Брук, Стонi Брук, США
2 Науковий обчислювальний центр, Нацiональна лабораторiя Брукхейвена, Уптон, США
3 Фiзичний факультет, Нацiональна лабораторiя Брукхейвена, Уптон, США
4 Лабораторiя Дж. Генрi, Прiнстонський унiверситет, Прiнстон, США
Описано iнтерфейс-простежуючий чисельний алгоритм для моделювання магнiтогiдродинамiчних
(MHD) потокiв мультифаза/вiльна поверхня у наближеннi низькомагнiтних чисел Рейнольдса
(Samulyak R., Du J., Glimm J., Xu Z., J. Comp. Phys., 2007, 226, 1532). Алгоритм застосовується
у мультифiзичному кодi FronTier i використовується для моделювання MHD процесiв у рiдинi
i слабо iонiзованiй плазмi. У цiй статтi здiйснено чисельне моделювання рiдкого ртутного
струменя, що входить у сильне i неоднорiдне магнiтне поле i взаємодiє з потужним протонним
iмпульсом, а також проведено порiвняння з експериментом. Такий ртутний струмiнь є прототипом
запропонованого майбутнього прискорювача частинок Мюонний коллайдер/Нейтринна фабрика.
Моделювання демонструє елiптичну дисторцiю ртутного струменя, коли вiн входить у магнiтне
поле при малому кутi до магнiтної осi, нестiйкостi струмiнь-поверхня (фiламентацiю), iндукованi
взаємодiєю з протонними iмпульсами i стабiлiзуючий ефект магнiтного поля.
Ключовi слова: MHD алгоритм, мультифазний потiк, ртутна мiшень
43402-12
Introduction
Governing equations
Numerical algorithms and implementation
Front-tracking algorithm for MHD equations
Models for phase transition and cavitation
Applications to accelerator targets
Simulation of the mercury jet entering a solenoid magnet
Simulation of the mercury jet interacting with proton pulses.
Conclusions and future work
|