# Magnetic field draping around clumpy high-velocity clouds in galactic halo

S. Lyla Jung1, Asger Grønnow2, N. M. McClure-Griffiths1

1 Research School of Astronomy & Astrophysics, The Australian National
University, Canberra ACT 2611, Australia
2 Kapteyn Astronomical Institute, University of Groningen, Landleven 12, 9747
AD Groningen, The Netherlands e-mail<EMAIL_ADDRESS>
(Last updated)

###### Abstract

Throughout the passage within the Galactic halo, high-velocity clouds (HVCs)
sweep up ambient magnetic fields and form stretched and draped configurations
of magnetic fields around them. This magnetized layer at the cloud-halo
interface can suppress the mixing between the cloud and the surrounding halo
gas. For simplicity, many earlier numerical studies adopt spherically
symmetric uniform-density clouds as initial conditions. However, observations
of HVCs indicate that the clouds are clumpy and turbulent. In this paper, we
perform 3D magnetohydrodynamic simulations to study the evolution of clouds
with more a realistic initial density distribution. In the parameter regime of
our models, simple spherical clouds only grow from the onset of the
simulations, while clumpy clouds lose cold gas at early times and then start
to grow unless they have very low metallicity. Some massive clumps composing a
clumpy cloud share broad similarities with single uniform-density clouds in
that they have clear head-tail morphology and magnetic fields draped around
themselves. Such similarities are marginal among clumps in the wake of the
cloud. With magnetic fields present, the growth of hydrodynamic instabilities
is suppressed along the direction of the cloud’s motion. Efficient radiative
cooling leads to a compact cloud with rich clumps and filaments condensed from
the surrounding gas. The slope of the initial density power spectrum is
closely related to the density contrast within a clumpy cloud. It, therefore,
determines how fast the material is pushed out of the initial cloud and mixed
with the halo gas.

###### keywords:

methods: numerical – magnetohydrodynamics (MHD) – instabilities – ISM: clouds
– galaxies: magnetic fields

††pubyear: 2022††pagerange: Magnetic field draping around clumpy high-velocity
clouds in galactic halo–Magnetic field draping around clumpy high-velocity
clouds in galactic halo

##  1 Introduction

High-velocity clouds (HVCs) are part of the multi-phase circumgalactic medium
(CGM) of the Milky Way. Though it is the hot gas ($T\sim 10^{6}\,\rm K$;
Henley & Shelton 2015) that fills much of the volume of the Galactic halo, a
significant portion of the CGM is in cooler gas phases, e.g., neutral hydrogen
(HI). The observed sky coverage of HI HVCs is $\sim 15\%$ with the column
density threshold of $N_{\rm HI}>2\times 10^{18}\,\rm cm^{-2}$ (Westmeier
2018) and this increases to $\sim 37\%$ with a lower threshold $N_{\rm
HI}>7\times 10^{17}\,\rm cm^{-2}$ (Murphy et al. 1995; Lockman et al. 2002).
The total mass estimate of the observed HVCs is largely dependent on the
accuracy of the distance measurements to the clouds. Putman et al. (2012)
estimate the total mass of HVCs not associated with the Magellanic System to
be $\sim 7.4\times 10^{7}\,\rm M_{\odot}$. When the Magellanic system (the
Magellanic Stream, Leading Arm, and Bridge) is included, the estimated HI mass
is an order of magnitude higher.

These cold phase gas clouds are a potential source of gas that could fuel the
future star formation activity of the Galaxy. The gas inflow and outflow rates
in the Milky Way halo are estimated from the observed kinematic properties of
HVCs. Considering HI HVCs including the Magellanic Stream, Richter (2012)
gives the inflow rate $\dot{M}_{\rm in}{\;\raise
1.29167pt\hbox{$>$\kern-7.5pt\raise-4.73611pt\hbox{$\sim$}}\;}0.7\,\rm
M_{\odot}yr^{-1}$. Similarly, Fox et al. (2019) estimate the inflow rate
$\dot{M}_{\rm in}{\;\raise
1.29167pt\hbox{$>$\kern-7.5pt\raise-4.73611pt\hbox{$\sim$}}\;}0.53\,\rm
M_{\odot}yr^{-1}$ and the outflow rate $\dot{M}_{\rm out}{\;\raise
1.29167pt\hbox{$>$\kern-7.5pt\raise-4.73611pt\hbox{$\sim$}}\;}0.16\,\rm
M_{\odot}yr^{-1}$ excluding the contributions of the Magellanic Stream and the
Fermi Bubbles (similar results are also presented in Clark et al. 2021). It is
clear that the Milky Way not only has a very apparent accretion of the
Magellanic system, but also the non-Magellanic origin HVCs are overall
inflowing into the Galaxy rather than outflowing. Such inflow might be the
dominant source of gas needed to sustain the star formation activity of the
Milky Way.

However, it is still an open question how cold clouds can survive the dynamic
environment of the CGM and deliver gas to the Galaxy. Surrounded by hot halo
gas, clouds falling towards the Galaxy are under severe ram pressure that can
strip gas away from the cloud. Also, hydrodynamic instabilities (typically,
Rayleigh-Taylor and Kelvin-Helmholtz instabilities) provoke mixing across the
cloud’s interface between the two media with different density and velocity.
Observed head-tail morphology (e.g., For et al. 2013) and complex small-scale
structures arising across the edges of HVCs (Barger et al. 2020) are direct
evidence of on-going ablation of the clouds moving at high-velocities.

On the other hand, the mixing of cold and hot gas produces gas at intermediate
density and temperature where radiative cooling can be efficient. Under
certain conditions where the cooling is faster than the mixing and stripping,
the total cold gas mass of the cloud system can grow with time (e.g.,
Marinacci et al. 2010; Armillotta et al. 2016; Gronke & Oh 2018; Li et al.
2020; Heitsch et al. 2021; Grønnow et al. 2022).

Introducing further complications, magnetic fields in the halo significantly
affects the overall evolution of the clouds. When a cloud is moving super-
Alfvenically, it sweeps up the magnetic field lines in the halo and develops a
strongly magnetized layer at the cloud-halo interface. The effect of magnetic
fields on the survivability of clouds moving at high velocities has been
actively studied from the numerical simulations’ perspective (Jones et al.
1996; Gregori et al. 1999, 2000; Santillán et al. 1999; Dursi & Pfrommer 2008;
Shin et al. 2008; Kwak et al. 2009; McCourt et al. 2015; Banda-Barragán et al.
2016; Grønnow et al. 2017; Banda-Barragán et al. 2018; Cottle et al. 2020;
Sparre et al. 2020; Grønnow et al. 2022). Earlier models clearly show that the
magnetized layer provides tension and suppresses hydrodynamic instabilities
across the field lines. Yet, the extent of stability provided by magnetic
fields to the clouds is not well constrained as it depends on, for example,
types of instabilities, the strength of the magnetic field, and the direction
of the field with respect to the clouds’ motion. Earlier work by Gregori et
al. (1999) found enhanced Rayleigh-Taylor instability along the direction
perpendicular to the cloud’s motion and the initial ambient magnetic field.
More recent works with higher resolution simulations give mixed results
especially when radiative cooling is included: for example, Sparre et al.
(2020) and Cottle et al. (2020) find only a small effect of magnetic fields on
the survival of the clouds while McCourt et al. (2015) and Grønnow et al.
(2018) predict a large effect. In either case, studies agree that the magnetic
field does not lead to faster destruction.

Many of the studies mentioned above, with notable exceptions of Banda-Barragán
et al. (2018); Banda-Barragán et al. (2019), assume a simple spherical cloud
with a uniform density core for initial conditions. In such a setting, the
size of the cloud is one of the important deciding factors of the survival
timescale of the clouds. For example, Li et al. (2020) present that the
lifetime of their simulated clouds closely follows a simple power-law
dependency on the cloud size and other properties (see their equation 19 for
the exact form of the fitting function). However, observations clearly show
that HVCs are in fact clumpy and turbulent, similar to much of the
interstellar medium (ISM); a large cloud complex may have multiple dense cores
and the spatial power spectrum of the column density follows a power law
distribution (e.g., Marchal et al. 2021). Hence, there is a need to examine to
what extent predictions from such spherical cloud simulations are valid when
it comes to clouds with more realistic structures. There have been a handful
of simulations probing the evolution of clouds with complex and realistic
density structures, however, many of them do not include magnetic fields
(e.g., Patnaude & Fesen 2005; Cooper et al. 2009; Schneider & Robertson 2017)
or are not in the parameter regime (for example, in terms of metallicity or
relative velocity) applicable to infalling HVCs (e.g., Banda-Barragán et al.
2018; Banda-Barragán et al. 2019).

In this paper, we investigate the properties and the evolution of clouds that
have an initial density distribution motivated from observations.We choose the
physical parameters of the fiducial simulation, e.g., gas metallicity, cloud
velocity, and density structure, to be in line with the observed HVCs. More
details on model settings are in Section 2. In Section 3, we compare the
evolution of the total cold gas mass of each model and evaluate the survival
of the clouds under various physical conditions. In Section 4, we discuss
processes that shape the evolution of the clouds by comparing structural
properties of the clouds in different models. Discussions on the resolution
effect, caveats and limitations in interpreting our models are in Section 5.1
and Section 5.2. We connect the results of the simulations to observed
properties of the Milky Way HVCs and make predictions for future observations
in Section 5.3. Section 6 is the summary and conclusion.

##  2 Methodology

###  2.1 Simulations overview

Figure 1:  The top panel shows the cloud location and the distribution of AMR refinement levels at the initial snapshot of the MC model. The clumpy density structures are saturated due to the shared color scale that is more suitable for showing the cloud in the bottom panel. The bottom panel is the same information at the final snapshot ($t=44\,\rm Myr$). In both panels, we present the density distribution of gas with $C>0.25$ in order to compare it to the refinement levels. The density shown here is the maximum density along the projection (z-axis). We confirm that any structures at $C>0.25$ are covered by cells at the maximum refinement level (red) in both the snapshots.  Table 1: A list of simulations and input parameters. Column 1 is the name of the models, column 2 is the type of initial cloud presented in Table 2, and columns 3 to 6 are the gas velocity, density, metallicity, and magnetic field strength, respectively. Column 7 is the size of the smallest cell in the simulation domain, column 8 shows whether radiative cooling is implemented, and column 9 is the slope of the initial spatial power spectrum of the clouds.  | Cloud type | vhalo ($\rm kms^{-1}$) |  $\rho_{\rm halo}$ ($\rm gcm^{-3}$) | Metallicity (Z⊙) | B0 ($\rm\mu G$) | Max resolution (pc) | Cooling | Slope
---|---|---|---|---|---|---|---|---
MC | Clumpy-R32 | 200 | 10-27 | 0.3 | 0.1111This corresponds to plasma $\beta=P_{\rm thermal}/P_{\rm B}\approx 350$. | 3.125 | Yes | -3
HC | Clumpy-R32 | 200 | 10-27 | 0.3 | 0 | 3.125 | Yes | -3
MU | Uniform | 200 | 10-27 | 0.3 | 0.1 | 3.125 | Yes | -
HU | Uniform | 200 | 10-27 | 0.3 | 0 | 3.125 | Yes | -
MC-Z003 | Clumpy-R32 | 200 | 10-27 | 0.03 | 0.1 | 3.125 | Yes | -3
MC-NoCool | Clumpy-R32 | 200 | 10-27 | - | 0.1 | 3.125 | No | -3
MC-p2.7 | Clumpy-p2.7 | 200 | 10-27 | 0.3 | 0.1 | 3.125 | Yes | -2.7
MC-p3.5 | Clumpy-p3.5 | 200 | 10-27 | 0.3 | 0.1 | 3.125 | Yes | -3.5
MC-R16 | Clumpy-R16 | 200 | 10-27 | 0.3 | 0.1 | 6.25 | Yes | -3
MC-R64 | Clumpy-R64 | 200 | 10-27 | 0.3 | 0.1 | 1.5625 | Yes | -3
Table 2: A list of initial cloud type used in different models. See Fig. 2 for their 2D column density distributions. Cloud type |  $\chi$222$\chi\equiv\tilde{\rho}_{\rm cloud}/\rho_{\rm halo}$, where $\tilde{\rho}_{\rm cloud}$ is the median density of all gas cells with $C=1$. |  |  $M_{\rm cloud}$ [$\rm M_{\odot}$]
---
($T<2\times 10^{4}\,\rm K$)
|  $M_{\rm cloud}$ [$\rm M_{\odot}$]
---
($C$333See equation 3 and the related paragraph for the definition of
$C$.$=1$)
Uniform | 197 | $1.21\times 10^{4}$ | $1.12\times 10^{4}$
Clumpy-R16 | 149 | $1.23\times 10^{4}$ | $1.22\times 10^{4}$
Clumpy-R32 | 142 | $1.23\times 10^{4}$ | $1.23\times 10^{4}$
Clumpy-R64 | 133 | $1.23\times 10^{4}$ | $1.23\times 10^{4}$
Clumpy-p2.7 | 104 | $1.23\times 10^{4}$ | $1.21\times 10^{4}$
Clumpy-p3.5 | 181 | $1.23\times 10^{4}$ | $1.24\times 10^{4}$

We perform a set of 3D magnetohydrodynamic (MHD) simulations of moving clouds
using the PLUTO code version 4.1 (Mignone et al. 2012). Assuming infinite
conductivity, the fluid and magnetic fields in the simulation box evolve
following the ideal MHD equations. We use the HLLC Riemann solver (Toro et al.
1994) which has a good balance of accuracy and stability for all our
simulations. We have run our fiducial MHD simulation using HLLD solver as
well, because the HLLD solver is often considered to be a more accurate
extension of HLLC for MHD cases. We have verified that the differences between
the models with HLLC and HLLD solvers are not significant in our simulation
setting444Later in Fig. 4, we show the cold gas evolution in the HLLD run.. We
apply the hyperbolic divergence-cleaning algorithm introduced by Dedner et al.
(2002, see also ; ), that approximately enforces the solenoidal constraint of
magnetic fields ($\nabla\cdot\textbf{B}=0$).

The simulation domain is a rectangular box elongated along the $x$-axis, i.e.,
($L_{\rm x},L_{\rm y},L_{\rm z})=(10,2,2)\,\rm kpc$. Boundary conditions are
set to outflow, except at the gas injection boundary at $x=-1$. From the
injection boundary, hot halo gas is injected along the $x$-axis at a constant
positive velocity ($v_{\rm halo}$). Initially, we place an overdense structure
(the “cloud") centered at $1\,\rm kpc$ away from the injection boundary (i.e.,
$x=0\,\rm kpc$). See the middle panel of Fig. 1 for the initial setting of the
simulation domain. The rest of the volume is filled with halo gas with the
same properties as the medium injected at the boundary. The cloud is initially
static with respect to the simulation domain and in pressure equilibrium with
the halo gas.

We include radiative cooling via collisional ionization equilibrium based on
the tables of Sutherland & Dopita (1993), except for a model in that we turned
off the cooling for comparison study. The cooling rate of gas depends on the
metallicity and temperature of the gas. Cooling is switched off below a
temperature threshold of $10^{4}\,\rm K$. The metallicity is uniform
throughout the simulation domain.

We use the adaptive mesh refinement (AMR) technique to ensure that any dense
structures in the simulations are covered with high-resolution grids
throughout the time domain of the simulations. For the standard resolution, we
use a base grid of (200, 40, 40) cells along each axis and then allow
refinement by up to 4 levels based on the density gradient. Each level refines
cells to half the size of the previous level. Thus, the minimum AMR cell size
is $(10\,\rm{kpc}/200)\times 2^{-4}=3.125\,\rm pc$, equivalent to 32 cells per
cloud radius. See Fig. 1 for the distribution of cells in different refinement
levels at the initial (upper panel) and the final (lower panel) snapshot of
the fiducial simulation run.

Table 1 summaries the input parameters of the models in this study. In all
models, the density of the halo gas is $\rho_{\rm halo}=10^{-27}\,\rm
gcm^{-3}$, i.e., the particle number density $\sim 0.001\,\rm cm^{-3}$ and the
temperature $8.5\times 10^{5}\,\rm K$ that are representative of the inner
Milky Way corona (e.g., Henley & Shelton 2015; Miller & Bregman 2015). The
initial velocity of the clouds with respect to the surrounding halo gas is
$200\,\rm kms^{-1}$, that corresponds to a Mach number

$\mathcal{M_{\rm halo}}=\frac{v_{\rm halo}}{c_{\rm s}}=1.3,$ (1)

where the sound speed of the halo medium is $c_{\rm
s}=\sqrt{\frac{5P}{3\rho_{\rm halo}}}$ as we consider all gas in our models as
monatomic ideal gas and $P$ is the gas pressure. This choice of velocity is
motivated from the deviation velocity measurements of observed HVCs, i.e., the
deviation of observed clouds’ velocity from the Galactic rotation velocity
(Moss et al. 2013).

The halo gas has uniform magnetic field strength $B_{\rm 0}=0.1\,\rm\mu G$,
unless specified as non-magnetized models. This choice of the magnetic field
strength represents the halo gas $10\,\rm kpc$555This is a typical distance to
observed Milky Way HVCs (e.g., Wakker 2001; Thom et al. 2008). above the
Galactic disk in the solar neighbourhood (i.e., $8\,\rm kpc$ radius from the
Galactic center) according to the Galactic magnetic field model of Jansson &
Farrar (2012). Physically, this corresponds to the plasma $\beta=P_{\rm
thermal}/P_{\rm B}\approx 350$ and places the medium in the super-Alfvenic
regime666The Alfvenic Mach number $\mathcal{M_{\rm A}}=v_{\rm halo}/v_{\rm
A}\approx 22.4$, where the Alfven velocity is $v_{\rm A}=\sqrt{\frac{B_{\rm
0}^{2}}{4\pi\rho_{\rm halo}}}\approx 8.9\,\rm kms^{-1}$. The direction of the
halo magnetic field lines are perpendicular to the direction of the cloud
motion. Earlier studies by Banda-Barragán et al. (2016) and Grønnow et al.
(2018) show that only the transverse component drapes around the cloud and is
amplified even when there are magnetic field components both parallel and
transverse to the motion. As a result, the transverse component tends to
dominate the evolution for all field orientations in most cases, other than
when the halo magnetic fields are very close to parallel with negligible
transverse magnetic fields.

  * •

The MC model is a fiducial simulation representing typical observed HVCs in
the Milky Way halo. The metallicity is set to $0.3\,\rm Z_{\odot}$ which is
motivated by the typical range of metallicities of observed HVCs (e.g.,
Complex A $\sim 0.1\,\rm Z_{\odot}$, Complex C $0.1-0.3\,\rm Z_{\odot}$,
Magellanic Stream $0.1-0.3\,\rm Z_{\odot}$, Smith cloud $\sim 0.5\,\rm
Z_{\odot}$; from reviews by Putman et al. 2012 and Richter 2017, and
references therein). The cloud in this model has a clumpy density structure.
Details on how we generate the initially clumpy cloud will be presented
shortly in Section 2.2.

  * •

The HC model is a pure hydrodynamic (HD) counterpart of the fiducial model:
all the physical parameters are identical except that there are no magnetic
fields in the halo gas. See Section 4.2 for the detailed comparison of the
evolution of clouds in the two models.

  * •

The MU model and the HU model are MHD and HD simulations with a spherically
symmetric initial density cloud instead of the clumpy cloud.

  * •

For testing the effect of radiative cooling on the evolution of clouds, we
have the MC-Z003 model with a low metallicity ($0.03\,\rm Z_{\odot}$, i.e.,
less efficient cooling) throughout the simulation volume. This metallicity
roughly corresponds to the metallicity of HVCs with pristine intergalactic
origin. The MC-NoCool model has radiative cooling turned off completely.

  * •

We study the effect of the varying the slope of the initial density power
spectrum by comparing the MC-p2.7 model and the MC-p3.5 model which have the
power slope shallower (-2.7) and steeper (-3.5) than the MC model (-3),
respectively. See Section 2.2 on how the structure of a clumpy cloud changes
with the varying spectral slope.

  * •

We test the resolution effect using models with the maximum resolution lower
(MC-R16) and higher (MC-R64) than the fiducial MC model.

###  2.2 Initial cloud density distribution

For our “clumpy” cloud simulations in Table 1, we initialise the density
structure of the cloud to resemble the observed log-normal density
distribution of gas clouds. As a reference, we use a $21\,\rm cm$ HI map of
the Smith Cloud observed with the Robert C. Byrd Green Bank Telescope of the
National Radio Astronomy Observatory from Lockman et al. (2008). We
parameterise the properties of the cloud based on the spatial power spectrum
(SPS) and the column density probability distribution function (PDF). Only
pixels with the column density above 26 times the noise level (i.e.,
$N>10^{19.8}\,\rm cm^{-2}$) are considered in this process. Like many other
areas of the observed ISM (e.g., Pingel et al. 2018), the SPS of the HVC can
be readily approximated in a power-law form $P(k)=Ak^{\gamma}$, where $k$ is
the spatial frequency. The slope $\gamma$ represents the relative strength of
density fluctuations in different spatial frequencies. The amplitude $A$
determines the contrast of such fluctuation, which is directly related with
the width of the observed column density PDF ($\sigma_{\rm s}$) when $\gamma$
is determined. We use the SPS analysis toolkit Turbustat Python package (Koch
et al. 2019) to identify $\gamma$ from the HI image. We estimate $\sigma_{\rm
s}$ by fitting a log-normal distribution to the column density PDF. The
resulting $\gamma$ is -3.0 and $\sigma_{\rm s}$ is 0.24. Although these
parameters are direct measurements of the Smith cloud, we emphasize that the
purpose of this paper is not to conduct a case study of the Smith cloud.
Discussions in this paper are rather applicable to any clumpy clouds moving at
a high-velocity with respect to the surrounding gas.

We generate a 3D array of an overdensity ($\delta=\rho/\left<\rho\right>-1$)
that follows the statistics ($\gamma$ and $\sigma_{\rm s}$) obtained from the
above procedure777In addition to $\gamma=-3$, we also build clouds with
$\gamma=-2.7$ and -3.5 for testing the dependence of the initial slope of the
density power spectrum. The power spectrum amplitude is fixed in all cases..
We use the powerbox Python module (Murray 2018) which generates a periodic
homogeneous scalar field with a given power spectrum. We first generate a cube
at the highest resolution adopted in this work (i.e., 64 cells per cloud
radius; Clumpy-R64) and then degrade the 3D array to lower resolutions when
needed (e.g., Clumpy-R32 and Clumpy-R16). By doing so, we ensure that the
density distributions at different spatial resolutions are fundamentally the
same structures.

Next, we apply a spherical kernel to the overdensity cube to set a boundary of
the cloud that smoothly fades into the halo. The profile of the kernel is
described as a function of the distance from the kernel center (r):

$f(r)=\frac{1}{2}(n_{\rm cloud}-n_{\rm
halo})\left(1-\tanh\left[s\left(\frac{r}{r_{\rm
cloud}}-1\right)\right]\right),$ (2)

where $n_{\rm cloud}$ is the core number density of the cloud, and $n_{\rm
halo}$ is the gas number density in the halo. The steepness of the transition
between the cloud and the surrounding medium is controlled by the parameter
$s$. By definition, $r_{\rm cloud}$ ($=0.1\,\rm kpc$) is a radius where the
kernel density is halfway between $n_{\rm cloud}$ and $n_{\rm halo}$. We set a
relatively sharp kernel edge ($s=30$) that corresponds to a transition between
the cloud and the halo gas over $\sim 4$ cells at the standard resolution. We
move this spherical kernel around the overdensity cube until the center of the
mass of the structures within the kernel matches the kernel’s center within
the size of one cell. This step is to prevent any imbalance of the mass
distribution within the clumpy cloud. Once the kernel is multiplied to the
overdensity field, $n_{\rm halo}$ is added to all the cells. Finally, the cube
is plugged in to the simulation as an initial density field.

For the models with a spherically symmetric initial cloud density profile (MU
and HU in Table 1), we adopt a density profile that follows equation 2 with
the steepness parameter $s=10$. In this case, the cloud has a flat, uniform
density core and its edge smoothly fades into the surrounding halo gas.
Hereafter, we refer to this cloud structure as a “uniform density” cloud.

Figure 2:  The density structure of different initial clouds presented in
table 2. The color scale is fixed for easy comparison between the panels.
Figure 3:  The initial density PDF of the clouds shown in table 2.

In Fig. 2 we present the column density distribution of the clouds used in
this work. The color scale is fixed for all the panels for easy comparison
between the panels. The top three panels are the clouds with different power
spectrum slope $\gamma=-2.7,-3,$and $-3.5$ that we use for the MC-p2.7 model,
the MC model, and the MC-p3.5 model, from left to right. The choice of slopes
used in this study is made within the observed range of the slopes in the
Galactic ISM presented in Pingel et al. (2018) and references therein. A
shallower power spectrum slope (Clumpy-p2.7) leads to a clumpier cloud with
individual clumps highlighted and larger density contrast. In contrast, with a
steeper slope (Clumpy-p3.5), the density contrast within the cloud is smaller
and the clump structures are more diffuse. In the bottom panels of Fig. 2, we
show the clouds for lower and higher resolution models, and the uniform
density cloud.

Fig. 3 presents the density PDF of the initial clouds. The fiducial Clumpy-R32
cloud is shown in red solid line. The red dashed vertical line marks the core
density of the uniform density cloud. Depending on the power spectrum slope,
Clumpy-p2.7 cloud (green solid line) has a wider distribution and Clumpy-p3.5
cloud (green dashed line) has a narrower distribution which is in line with
the density contrast visible in Fig. 2. The PDFs of the clouds with different
resolutions are overall similar, with the distribution of Clumpy-R64 cloud
(blue dashed line) being slightly wider than that of the low-resolution cloud
(Clumpy-R16; blue solid line).

As mentioned earlier, the cloud is initially in pressure equilibrium with the
surrounding halo gas. As a result, the log-normal density distribution of the
cloud means that the temperature of the gas spans a certain range; a fraction
of initial cloud gas cells are in temperatures above the cooling threshold
($T>10^{4}\,\rm K$). With the radiative cooling implemented, this gas could
cool down to lower temperatures during the first few million years of the
simulation. We test the stability of the initial density field by letting the
cloud evolve with zero velocity relative to the surrounding halo gas. Because
the initial mass budget in the intermediate-temperature cells is low ($\sim
3.5\%$ of the initial mass), the change in the density structure is
negligible. The cloud stays stable throughout the time scale of the
simulations without any external influence.

In order to keep track of the advection of the material originally belonging
to the cloud, we introduce a passive scalar tracer $C$ that advects with the
flow:

$\frac{\partial(\rho C)}{\partial t}+\nabla\cdot(\rho C\bm{v})=0.$ (3)

Initially, $C$ is set to 1 within the radius $r_{\rm cloud}=0.1\,\rm kpc$ from
the cloud center and 0 elsewhere. In the following sections, we use $C$ as an
indicator of the degree of mixing between the cloud and the halo gas within a
certain cell.

Table 2 summarizes the properties of the initial clouds. The second column is
the density ratio between the cloud and the halo gas ($\chi$) calculated by
adopting the median density of all gas cells with $C=1$ as a representative
density of a cloud. In the third and the fourth column we present the mass of
the cloud calculated following two different definitions: (i) the mass of the
cold phase gas ($T<2\times 10^{4}\,\rm K$) and (ii) the mass within $r_{\rm
cloud}=0.1\,\rm kpc$ from the clouds’ center where the passive scalar tracer
$C$ is set to 1. As previously mentioned when we discuss the stability of the
initial cloud, the difference between the two masses is minimal since much of
the clouds’ mass within $r_{\rm cloud}$ exists in the dense cold phase.

###  2.3 Clump identification

In this paper, we frequently refer to “clumps” in a cloud. These clumps are
dense structures identified using the FellWalker algorithm (Berry 2015)
implemented in the CUPID package (Berry et al. 2007). The algorithm identifies
clumps based on the gradients in the 3D density distribution. We use any cells
with their density higher than twice the density of the halo gas for the clump
finding. If a dip connecting two clumps is smaller than 30 times the halo gas
density, the two clumps are considered as one. We reject clumps with smaller
than four cells along one of the axes in order to avoid clumps that are
marginally resolved given the spatial resolution of the simulations. We refer
interested readers to Berry (2015) for detailed descriptions of the algorithm.
The total number of individual clumps identified is on the order of dozens in
the fiducial model when it is reasonably evolved from the initial condition.
The number changes over time and models. We discuss these aspects in Section
4.

##  3 The evolution of total cold gas mass

Figure 4:  The evolution of the total cold gas mass in the simulations. The
gas mass is normalized by the mass at $t=1\,\rm Myr$ for an easy comparison
between the models.

The survivability of a cloud is closely related to whether or not a high-
velocity cloud is capable of delivering cold gas to the Galactic disk. If the
radiative cooling of hot gas is faster than the cloud disruption, the total
mass of the cold gas ($T<2\times 10^{4}\,\rm K$) in the system grows over
time. If not, the cloud eventually dissipates into the hot halo and no cold
phase gas remains. A characteristic timescale of high-speed relative motion
impacting a cloud in a halo is often parameterized using the cloud crushing
timescale ($t_{\rm cc}$); the time for the shock induced from the initial
collision between the cloud and the surrounding medium to travel through the
radius of the cloud

$t_{\rm cc}=\frac{r_{\rm cloud}}{v_{\rm halo}}\sqrt{\frac{\rho_{\rm
cloud}}{\rho_{\rm halo}}}.$ (4)

Much effort has been made in parameterizing the survival of the cloud using
scale-free parameters, e.g., $t_{\rm cool}/t_{\rm cc}$, where $t_{\rm
cool}=\frac{3nk_{\rm B}T}{2n_{\rm H}^{2}\Lambda}$ (Gronke & Oh 2018; Li et al.
2020; Sparre et al. 2020). However, in this work we stick to rather
qualitative interpretations. This is because the actual degree of impact of
the cloud-halo interaction changes significantly over time as the initial
clumpy cloud breaks into multiple smaller pieces. Therefore, while $t_{\rm
cc}$ can still be a useful guide for some models, we would not expect it to be
strictly applicable as a fundamental timescale of the evolution when it comes
to clumpy clouds as in that of uniform density clouds.

It is important to clarify that we discuss the survival of the clouds by the
total amount of cold gas in the system in this section. By doing so, we are
open to the possibility of the original cloud body being completely destroyed,
while the mass of cold gas can still grow due to the formation of cold clumps
condensed from the hot halo gas. This method is suitable for determining
whether high-velocity cloud systems are capable of delivering gas to central
galaxies, regardless of whether the gas is from the original cloud or not.

Fig. 4 shows the evolution of the total mass of cold gas in different
models888Due to the limited computing resources available, not all models are
covered up to the time when the upturn of the cold gas mass might take place.
For those models, discussions in this work are limited to the initial
disruption phase of clouds.. We first focus on the overall evolution of the
fiducial clumpy MC model (red solid line). Between 0 to $30\,\rm Myr$, there
is a slight (3%) decrease in the gas mass which is then followed by a dramatic
growth from $35\,\rm Myr$. We confirm that the initial decrease in the mass is
mainly due to the disruption of the cloud’s main body, especially, the
stripping of the lower-density interclump medium. The later growth is due to
the condensation of the intermediate-temperature gas that forms as a result of
mixing in the wake behind the cloud. Further discussions on this topic will be
presented in the following sections where we take a detailed look at the
morphology of the clouds.

Uniform density clouds (MU and HU; red and blue dashed lines) are always in
the growth regime throughout the time covered by the simulations regardless of
the presence of magnetic fields (MU: 104% and HU: 108% at $t\approx 28$). In
comparison, the mass of clumpy clouds (MC and HC; red and blue solid line)
decreases during the same time frame. We explain this discrepancy comes from
the different density contrast between the uniform/clumpy clouds and the halo
gas (see $\chi$ in Table 2)999 Note that we match the initial mass of the
uniform and the clumpy cloud to be identical. : the uniform cloud is overall
denser than the clumpy cloud, therefore, it is more resistant to the initial
disruption. The clumpy cloud has individual overdensities of much higher
densities compared to the uniform density cloud, but its low-density
interclump gas, which is susceptible to ram pressure stripping, is effectively
removed at its disruption phase.

Now we compare the evolution of the magnetized (MC and MU; red solid and
dashed lines) and non-magnetized (HC and HU; red solid and dashed lines)
models. In the uniform cloud cases, the cold gas mass of the HU cloud is
always higher than the MU cloud and the separation between the two curves
grows with time. These results agree with previous studies with similar
settings, e.g., Grønnow et al. (2018) and Kooij et al. 2021.The role of
magnetic fields in suppressing the growth of clouds is often discussed and
interpreted as an outcome of inhibited Kelvin-Helmholtz mixing at the cloud-
halo interface (Grønnow et al. 2018, 2022). We show evidence of suppressed
hydrodynamic instabilities in the presence of ambient magnetic fields in
Section 4.2. The difference between the HC model and the MC model is small
throughout the time domain ($<1.8\%$). However, the decrease of mass in the MC
model slows down and reaches a plateau by $t=28\,\rm Myr$, whereas the HC
model does not show such trend. However, we cannot know the fate of the HC
cloud at later times from this simulation; the cloud may keep loosing cold
gas, but it could also turn around later like most of the other models.

With reduced or no radiative cooling, both the MC-Z003 model (bright purple
solid line) and the MC-NoCool model (bright purple dashed line) show
significant decrease in the cold gas (MC-Z003: 89% and MC-NoCool: 52% at
$t=28\,\rm Myr$). Only the initial $15\,\rm Myr$ of the evolution of the MC-
NoCool model is shown in Fig. 4 for visualization purposes. In Section 4.3, we
will show how the changed radiative cooling efficiency alters the density
structure and the flow pattern around the clouds. The clear discrepancy among
the MC model, the MC-Z0.03 model, and the MC-NoCool model hints for
possibilities of diverse evolutionary tracks of HVCs depending on their
metallicity (see Section 5.3.4 for further discussions).

Comparing the MC model to the models with shallower (MC-p2.7; green solid
line) and steeper (MC-p3.5; green dashed line) initial density power spectra,
we learn that varying the slope affects how significant the initial mass loss
is and how early the upturn of the mass we see in the evolution of the clumpy
cloud takes place. The MC-p2.7 model loses 3.5% of the initial cold gas mass
before it enters the growth phase at $t=26\,\rm Myr$. This is a faster gas
loss and earlier upturn of the mass compared to the evolution of the MC model
which loses 3% of the mass during the first $30\,\rm Myr$ and then starts to
grow in mass. The MC-p3.5 model, in comparison, shows only 1.5% initial mass
loss and has even later upturn that starts at around $t=36\,\rm Myr$. We
explain that this trend is determined by how fast the clouds’ main body breaks
up and mixes with the surrounding halo gas which leads to condensation of halo
gas in the wake. See Section 4.4 for further explanations.

##  4 Comparison of clouds evolved in various environments

Interesting insights come from comparing details between the models. By
comparing structural properties of clouds that started from the
identical/similar initial condition but evolved in different environments, we
can learn about what shapes the clouds’ structure throughout their evolution.
The evolution of the clouds are mainly governed by the following processes,
which will hence be the focus our analysis:

  1. 1.

Ram pressure pushes gas to the opposite direction of the cloud’s motion. Given
the high relative velocity between the halo and the cloud ($v_{\rm
halo}=200\,\rm kms^{-1}$), all our models are initially in the regime where
ram pressure stripping is highly effective. In particular, clump clouds in our
models span a range of densities, therefore, unlike a uniform density cloud in
which the ram pressure simply results in overall deceleration of the cloud,
the ram pressure can actually strip gas. As a result, the entire cloud system,
as well as individual clumps composing the clouds, develop a head-tail
structure. Detailed analysis and discussions are presented below in Section
4.1.

  2. 2.

Kelvin-Helmholtz and Rayleigh-Taylor instabilities develop across the cloud-
halo interface, leading to erosion of clouds from their edges. Such
instabilities grow differently with and without the presence of magnetic
fields, as magnetic tension provides extra stability to the mixing surface
(Chandrasekhar 1961). This topic will be discussed in Section 4.2 where we
perform thorough comparison between the MC model and the HC model.

  3. 3.

Radiative cooling takes place as cold cloud gas mixes with hot halo gas and
form unstable intermediate-temperature gas. As a result, dense clumps can form
in the wake behind the cloud main body (see. e.g., Marinacci et al. 2010;
Armillotta et al. 2016; Gronke & Oh 2020). The cooling rate ($\Lambda$) is
roughly proportional to the gas metallicity ($Z$). In Section 4.3, we examine
the effect of cooling on the cloud evolution by comparing the MC model
($Z=0.3Z_{\rm\odot}$) with its lower metallicity counterpart MC-Z0.03 model
($Z\approx 0.03Z_{\rm\odot}$) and the MC-NoCool model with cooling turned-off.
Also, the cooling rate is a strong function of gas temperature and density and
peaks at temperatures around $10^{5}\,\rm K$. Therefore, this intermediate
temperature gas is unstable and quickly cools down to lower temperatures $\sim
10^{4}\,\rm K$.

Figure 5:  The density (left column) and the magnetic field (right column)
distribution in different models. The model name is shown in the right corner
of each panel. The left column shows the maximum density along $z$-axis in
order to highlight dense clumps with small path lengths. In the right column,
magnetic field lines are coloured by the field strength along the lines.

In this section, we focus on the snapshot at $t\approx 28\,\rm Myr$, which
corresponds to approximately $4t_{\rm cc}$ since the onset of the simulations.
See equation 4 for the definition of $t_{\rm cc}$. That said, we do not impose
any further meaning on $t_{\rm cc}$ than a rough estimate of the model
environment; we do not attempt to, for example, utilize it as a quantitative
measure of the evolutionary stages of the clouds.

Fig. 5 displays all models performed in this study at a fixed time of
$t=28\,\rm Myr$. The left column (from panel a-1 to h-1) shows the density
distribution projected along the $z$-axis (i.e., the axis perpendicular to the
halo magnetic field and the direction of the cloud’s motion). The color shows
the maximum density along the projection axis in order to display dense
structures with small path lengths that are not evident if we plot the column
density distribution. The right column (from panel a-2 to h-2) is the magnetic
field lines colored by the field strength. Models that started from the clumpy
cloud initial condition (MC, HC, MC-Z003, MC-NoCool, MC-p2.7, and MC-p3.5)
share broad similarities in their morphology: the main cloud body has three
prominent “heads” and there are multiple filamentary tails extended along the
direction of the motion. In comparison, models with uniform density initial
cloud (MU and HU) have a razor-thin tail. Earlier studies have shown that such
an extremely thin plane-shaped tail is an artefact of the spherically
symmetric initial conditions for the cloud. When the perfect symmetry is
broken by, for example, random internal velocities (Kooij et al. 2021) or a
random internal magnetic field (McCourt et al. 2015), such a tail does not
form even behind uniform density clouds.

###  4.1 The fiducial model (MC)

Figure 6:  Panel (a): the density distribution of the MC model. The density
presented here is the maximum density along the projection axis. Panel (b),
(c), (d), and (e): the density and magnetic field (yellow arrows) slice along
different parts of the MC cloud. Each slice is taken along the centre of mass
$z$-coordinate of each clump. The three (panel b, d, and e) of the four panels
are showing the main heads of the cloud and panel (c) is showing the wake of
the cloud. The coordinates in the bottom right are the center of the mass
locations of the clumps shown in the panels.

We start by taking a close look at the cloud in the fiducial MC model. Fig. 6
shows the density and the magnetic field distribution in four selected regions
of the cloud. The main panel (a) in the middle is the projected maximum
density distribution of the cloud, the same as panel (a)-1 of Fig.5. Three of
the four zoom-in panels (panel b, d, and e) are slices of the each head of the
cloud and panel (c) takes a slice at the cloud’s wake where we find groups of
small cloudlets and filaments. The (x,y,z) coordinates in the bottom right
corner shows the center of mass location of the clumps presented in each
panel. The yellow arrows show the direction of magnetic field vectors on the
(x, y) slice and the size of the arrows is proportional to the field strength.

The slices along the heads (panel b, d, and e) share common features: there is
a clear head-tail density structure and the magnetic field lines are draped
around the overdensities. To a certain extent, one can say that each head acts
similarly to what earlier studies find using simulations of uniform density
spherical clouds. On the other hand, the density and magnetic field slice in
panel (c) show somewhat different features. Several clumps are threaded into
long filaments that extend along the direction of the cloud’s motion. The
magnetic field does not react to individual overdensities, but rather follows
a large-scale configuration that is affected by the magnetic field draping
around the entire cloud system. As a result, magnetic field vectors are
overall parallel to the filaments.

Figure 7: The properties of individual clumps in the MC model. Panel (a): 2D
projection of the clumps. The dashed line marks the center of mass location of
the cloud on the $x$-axis. Panel (b): density profile of individual clumps
along the $x$-axis. Panel (c): the relative velocity of the clumps as a
function of the median density. The size of the circles corresponds to the
mass average of the passive scalar tracer ($C$), i.e., the fraction of the
mass in each clump that was initially part of the cloud.

Using the clump finding algorithm described in Section 2.3, we identify 28
clumps in this snapshot. See panel (a) of Fig. 7 for their distribution on the
(x, y) projected plane. In this panel, as well as panel (b) and (c), clumps
are colour-coded differently for easy distinction. The vertical dashed line
marks the center of the mass of the cloud along the $x$-axis. The same line is
in panel (b) as well.

We present the density profiles of each clump along the $x$-axis in panel (b).
The solid line shows the median density at a given location $x$ in each clump
and the shaded area encloses $68\%$ of the distribution. Here, we determine
the head-tail structure as a density profile that decreases with increasing x.
Low density clumps are located further away from the cloud’s head which
demonstrates that the entire cloud system has a head-tail morphology When
focusing on clumps located close to the main body ($x{\;\raise
1.29167pt\hbox{$<$\kern-7.5pt\raise-4.73611pt\hbox{$\sim$}}\;}1.5\,\rm kpc$),
even individual clumps show the head-tail profile, whereas there is no such
clear trend among the clumps further away ($x{\;\raise
1.29167pt\hbox{$>$\kern-7.5pt\raise-4.73611pt\hbox{$\sim$}}\;}1.5\,\rm kpc$).
This is because the efficiency of ram pressure stripping as well as the growth
of Kelvin-Helmholtz instability depend on the relative velocity between a
clump and the surrounding halo gas.

In panel (c) of Fig. 7, we show the relative $x$-velocity between the clumps
and the surrounding gas (i.e., $v_{\rm rel}=v_{\rm halo}-v_{\rm x,clump}$,
where $v_{\rm x,clump}$ is the median $x$-direction velocity of a clump) as a
function of the median density of a clump. We only consider the $x$-direction
velocity here since the transverse velocity is almost negligible. Each circle
corresponds to a clump and the size of the circles is proportional to the
median passive scalar tracer ($0<=C<=1$) where $C=0$ means that all the cold
gas has condensed from hot gas while $C=1$ means that all the cold gas was
part of the initial cloud. See equation 3 for the definition of $C$. There is
a clear deceleration throughout the extent of the cloud with the velocity
gradient of about $130\,\rm kms^{-1}$ over $1.5\,\rm kpc$ across the
$x$-direction. The three head clumps have the relative velocities around
$140\,\rm kms^{-1}$. Compared to the initial velocity field where the halo gas
is moving at $v_{\rm halo}=200\,\rm kms^{-1}$ and the cloud is static with
respect to the simulation domain, these clumps are moving at $\sim 70\%$ of
their initial velocity and experiencing only half of the initial ram pressure
along the $x$-axis considering that the magnitude of ram pressure is
proportional to the square of the relative velocity ($P_{\rm ram}=\rho_{\rm
halo}v_{\rm rel}^{\rm 2}$).

In the very low velocity range (e.g., $v_{\rm rel}{\;\raise
1.29167pt\hbox{$<$\kern-7.5pt\raise-4.73611pt\hbox{$\sim$}}\;}20\,\rm
kms^{-1}$), there are some clumps with higher median density than the overall
density-velocity gradient. These clumps are located far behind the main cloud
($x>2\,\rm kpc$) and their passive scalar tracer is around 0.2 or smaller,
meaning that most of the mass comprising the clumps is from the surrounding
halo gas. We explain the formation of such clumps as follows: material
stripped from the main body acts as seeds for condensation of the halo gas and
decelerates as it accretes the low-momentum halo gas. The ram pressure
stripping and the growth of Kelvin-Helmholtz instability are negligible given
their small velocity. This explains why clumps at large $x$ do not reveal a
head-tail density profile.

In reality, Milky Way’s gravitational potential can further accelerate these
low-velocity clouds/clumps. Therefore, the external gravity field is required
for simulations to examine the evolution and possible accretion of these
condensed cold clumps in a more realistic sense. One example of such work is
Grønnow et al. (2022) where the authors show that even with the gravity
included many of the condensed clumps remain at low relative velocity with
respect to the halo gas as the magnetic tension decelerates the clumps. Such
effect has not been examined with clumpy clouds yet. We will postpone the
investigation to future studies as it is beyond the scope of this work.

###  4.2 Magnetic fields (MC vs HC)

In this section, we focus on the effect of magnetic fields by comparing the
fiducial MC model with its non-magnetized counterpart, the HC model. All other
physical and numerical settings are identical in both simulations. The density
distribution of the MC cloud and the HC cloud are shown in panel (a)-1, 2 and
(b)-1, 2 of Fig. 5. There are some similarities in their overall shape of the
clouds that arise from having the same initial condition, e.g., the presence
of the three main heads, however, clearly there are differences in their
detailed structures.

Figure 8:  The density of the 14 densest clumps in the MC model (top) and the
HC model (bottom).  Figure 9:  The effective size of MC (red $\times$), HC
(blue $\times$), MU (red circle), and HU (blue circle) clumps. Only the clumps
above $100\,\rm M_{\odot}$ are presented. Panel (a): the size along the two
directions perpendicular to the cloud’s motion ($r_{\rm z}$ vs $r_{\rm y}$).
Panel (b): the perpendicular size and the effective length along the direction
of motion ($r_{\rm yz}$ vs $r_{\rm x}$), where $r_{\rm yz}=(r_{\rm y}+r_{\rm
z})/2$. In both panels, the dashed line shows the one-to-one relation of the
two radii. The top and the right axis of the panels are the normalized
probability distribution of clumps in the MC model (red curve) and the HC
model (blue curve) and the location of the MU cloud (red arrow) and the HU
cloud (blue arrow).

As shown in the previous section, there 31 clumps in the MC model. With the
same clump finding algorithm, we find a significantly more clumps in the HC
model: in total 101 clumps with in general in the smaller mass range.
Furthermore, there is a notable difference in the overall shape of individual
clumps in the two models. In Fig. 8 we show the 14 densest clumps identified
in each model (top: MC, bottom: HC). Clumps in the MC model are in general
elongated along the $x$-axis, while HC clumps are fragmented into smaller
sizes. The difference in the size of the clumps along the $y$-axis is not as
significant, except for Clump 3 and 9 in the MC model that have low-density
tails notably extended along the $y$-direction.

For quantitative analysis, we measure the effective size of the clumps:

$r_{\rm i}=\sqrt{5\left(\left<X_{\rm i}^{2}\right>-\left<X_{\rm
i}\right>^{2}\right)},$ (5)

where $i$ stands for each of the x, y, and $z$-axis,

$\left<X_{\rm i}\right>=\frac{\int X_{\rm i}\rho dV}{\int\rho dV},$ (6)

and $C$ is the passive scalar tracer. In panel (a) of Fig. 9, we first compare
$r_{\rm y}$ and $r_{\rm z}$. Both $y$ and $z$-axis are perpendicular to the
direction of the cloud’s motion and the initial halo magnetic field is
parallel to the $y$-axis. In this figure, we limit the sample points to clumps
that are more massive than $100\,\rm M_{\odot}$. The red and blue crosses are
clumps from the MC model and the HC model, respectively. In addition, we
present the effective size of uniform density cloud models, MU (red open
circle) and HU (blue open circle). Both the uniform cloud models are
considered as one large clump contrary to the MC cloud and the HC cloud being
composed of a collection of smaller clumps. On the top and right side of the
panel, we show the normalized probability distributions of MC (red curve) and
HC clumps (blue curve). The red and blue arrows are the locations of the MU
and HU cloud on the $r_{\rm y}$ and $r_{\rm z}$ distribution. While the HU
cloud is perfectly axisymmetric along the transverse axes ($r_{\rm z}/r_{\rm
y}=1$), the MU cloud is clearly elongated along the $z$-axis, i.e., the
direction perpendicular to the halo magnetic field ($r_{\rm z}/r_{\rm
y}=2.23$). Such breach of the structural axisymmetry with the presence of the
ambient magnetic field is reported by earlier studies with magnetic fields
(e.g., Gregori et al. 1999; McCourt et al. 2015; Grønnow et al. 2017, 2018).
In the case of clumpy cloud models, we find that a couple of clumps in the MC
model are elongated along the $z$-axis as much as the MU cloud. Note that such
elongated shape is not as common among the smaller clumps and the distribution
of MC and HC clumps on the ($r_{\rm z}$, $r_{\rm y}$) plane is more or less
similar in this regime.

In panel (b) of Fig. 9, we compare the size of clumps along the $x$-axis
($r_{\rm x}$) to the size along the perpendicular axis. Again, we only present
clumps above $100\,\rm M_{\odot}$ here. We have already shown in the above
paragraph that $r_{\rm y}$ and $r_{\rm z}$ closely resemble each other for the
majority of MC and HC clumps, therefore, we take $r_{\rm yz}=(r_{\rm y}+r_{\rm
z})/2$ as a reasonable representation of the size of clumps along the plane
perpendicular to the motion. For the majority of the clumps, $r_{\rm x}$ is
larger than $r_{\rm yz}$, i.e., clumps are elongated along the $x$-axis. The
distribution of MC clumps are extended to larger $r_{\rm x}$ compared to that
of HC clumps (see the top axis for the normalized PDF) as we already pointed
out from the qualitative comparison in Fig. 8.

Figure 10:  The evolution of the spatial power spectrum in the MC model (red)
and the HC model (blue). The $x$-axis is the physical length scale ($1/k$) and
the $y$-axis is the ratio between the initial power and the power at a given
time $P(init)/P(t)$. The solid, dashed, and dotted lines show $t=28.3,19.6,$
and $9.7\,\rm Myr$, respectively. The gray vertical line shows the
observational beam size of the HI image from which the clump cloud’s initial
power spectrum was drawn. With time, $P(t)$ decreases significantly at small
scales, while there is no systematic evolution at the scale larger than the
beam size.

For further quantitative comparison between the MC model and the HC model, we
perform a spatial power spectrum (SPS) analysis. We obtain the 2D SPS of the
clouds by regriding all the AMR cells in the simulation domain to uniform
sized cells corresponding to the highest resolution and projecting the 3D
volume along the $z$-axis to get the column density distribution. We limit the
analysis to cells with density higher than twice the halo gas density
($>2\rho_{\rm halo}$). The 2D power-spectrum frequency scale is defined as
$k=\sqrt{k_{\rm x}^{2}+k_{\rm y}^{2}}$, where $k_{\rm x}$ and $k_{\rm y}$ are
the spatial frequencies along the $x$ and $y$-axes, respectively. In Fig. 10,
we show the evolution of the SPS for the MC model and the HC model. The
spatial frequency $k$ is converted to the physical length scale $1/k$ in units
of pc. The smallest scale shown in the figure is $\sim 10\,\rm pc$, which is
about three times the size of the maximum resolution cell of the simulation.
The vertical-axis of the figure is $P(init)/P(t)$, where $P(init)$ is the
median Fourier SPS amplitude at a given $k$ of the initial density
distribution and $P(t)$ is this amplitude at a given time $t$. The black
horizontal line is at $P(init)=P(t)$. The solid line shows the MC model (red)
and the HC model (blue) at $t=28\,\rm Myr$. Initially, both the MC model and
the HC model start from the same density power spectrum as described in
Section 2.2. At $t=28\,\rm Myr$, both the models have less power at smaller
length scales below $\approx 78\,\rm pc$ ($25\,\rm cells$) compared to the
initial spectrum, i.e., higher $P(init)/P(t)$. Furthermore, $P(init)/P(t)$ of
the MC model is always higher than that of the HC model at this length range
indicating that the growth of small-scale structures in the MC model is
considerably suppressed in comparison with the HC model. We confirm that these
trends are visible regardless of the projection axis, i.e., whether the 2D
column density is integrated along the initial magnetic field direction
(y-axis) or the axis perpendicular to it (z-axis).

For reference, we show $P(init)/P(t)$ at earlier snapshots in dotted
($t=9.7\,\rm Myr$) and dashed ($t=19.6\,\rm Myr$) lines. In both the MC and
the HC models, the difference between $P(init)$ and $P(t)$ grows over time to
the direction where the power at the small scales significantly decrease. This
is because smaller clumps are easier to be disrupted as seen from their small
$t_{\rm cc}$ (see equation 4) and possibly also due to the effect of numerical
diffusion. The gray vertical line shows the telescope beam size of the HI
observation that the reference SPS for the clumpy clouds is drawn from. The
angular beam size is converted to the physical size using the distance to the
Smith Cloud ($12.5\,\rm kpc$). It is important to point out that there is no
notable evolution of $P(t)$ at length scales larger than the beam. Therefore,
having an initial structure for the simulations built based on the image of
the HVC that has already been interacting with the halo gas does not conflict
with the SPS evolution we find in the models.

Meanwhile, $P(init)/P(t)$ of the MC model is always higher than that of the HC
model at small scales. The suppressed power at small scales in the MC model is
because of (i) a “tube” produced by magnetic field draping and (ii) the
suppression of hydrodynamic instabilities at the cloud-halo interface. As the
cloud material advects downstream clumps elongate along the cloud’s motion.
Magnetic fields drape around each clump and produce a tube along the $x$-axis
which pulls the cold gas along with it. This effect has been confirmed on the
largest scale where a single uniform cloud will develop a long thin tail
inside the draped tail field. In the HC model without magnetic fields, any
elongated structures can be easily perturbed as hydrodynamic instabilities
actively grow across the cloud-halo interface. This is why we see short
disturbed small clumps in this model. When there is external magnetic fields
in the halo, it provides extra tension at the cloud surface and stabilizes the
perturbation. Therefore, clumps can further extend along the advection, the
growth of small-scale structures is suppressed, and there is smaller power in
the SPS in smaller length scales.

###  4.3 Radiative cooling (MC vs MC-Z0.03 vs MC-NoCool)

We study the effect of radiative cooling on the cloud’s structural properties
by comparing three models with different cooling efficiency: the MC model with
$Z=0.3Z_{\odot}$, the MC-Z0.03 model with $Z=0.03Z_{\odot}$ that represents
the nearly pristine gas with less effective radiative cooling, and the MC-
NoCool model without the radiative cooling implementation. All other
parameters are set the same. The lower metallicity of the MC-Z0.03 model
compared to the MC model means less effective radiative cooling since the
cooling rate is roughly proportional to metallicity. Fig. 5 panel (a)-1, 2,
(e)-1, 2, and (f)-1, 2 show the density and magnetic field distributions of
the models at $t=28\,\rm Myr$. In this section, we will show that the more
efficient the radiative cooling is, the more clumps form in the wake and the
more compact the cloud main body.

Figure 11:  The distribution of gas on the density - passive scalar tracer $C$
plane, where $C$ is the mass fraction of the gas that was part of the original
cloud. The color of each hexagon bin represent the mean $x$-direction relative
velocity as defined in Section 4.1. The black line encloses hexagon bins with
cold gas ($T<2\times 10^{4}\,\rm K$). In the MC model, there are dense gas
that are significantly decelerated (purple color) and mostly condensed from
the halo gas (bottom right corner of the panel). In comparison, gas in the MC-
NoCool show a smooth transition from the dense clouds (top right) to low-
density halo gas (bottom left). The MC-Z0.03 model shows characteristics
inbetween the MC and the MC-NoCool models.

In Fig. 11, we present the distribution of all gas cells in the MC model
(left), the MC-Z0.03 model (middle), and the MC-NoCool model (right) on the
density - scalar tracer $C$ plane, where $C$ is the fraction of the gas that
was part of the original cloud and, therefore, a measure of the cold gas
condensation. The color of each hexagon bin is the mean relative velocity
($v_{\rm rel}=v_{\rm halo}-v_{\rm x,clump}$) of gas cells within the bin. The
black line encloses the hexagon bins with cold gas ($T<2\times 10^{4}\,\rm
K$). These plots need to be interpreted with caution as the mass in different
hexagonal bins differs and the large majority of the mass is, in all three
cases, contained in the dense remnants of the original cloud, i.e., in the
upper-right bins. As pointed out earlier in Section 4.1, the MC model has
decelerated (purple color) high-density gas cells that are mostly condensed
from the halo (low $C$). In contrast, we do not see such a population of gas
cells in the MC-NoCool model. Gas cells in the models show a smooth transition
from the high-density high-$C$ cloud material to the low-density low-$C$ halo
gas. This is because the material stripped out of the original cloud smoothly
mixes and diffuses into the halo without any further cooling. In the MC-Z0.03
model, radiative cooling exists but is not as efficient as in the MC model.
Therefore, the gas cells’ properties lie between the MC model and the MC-
NoCool model. This explains why we do not see many dense clumps and filaments
in the wake behind the clouds in the MC-Z0.03 model and the MC-NoCool model.

Figure 12:  The column density distribution of the MC model, the MC-Z0.03
model, and the MC-NoCool model viewed from the head-on direction, i.e.,
projected along the $x$-axis. The dotted and the dashed circles have the
radius enclosing 50% and 90% of the total cloud mass ($R_{\rm 50}$ and $R_{\rm
90}$). The ratio $R_{\rm 90}/R_{\rm 50}$ for the MC model, the MC-Z0.03 model,
and the MC-NoCool model is 1.53, 1.56, and 1.7, respectively.

Furthermore, we see the effect of radiative cooling not only in the region
behind the cloud, but also in the main cloud body itself. The change in the
compactness of the cloud is best presented in Fig. 12 which shows the column
density distribution of the clouds viewed from the axis parallel to the
cloud’s motion. The inner dotted (outer dashed) circle is $R_{\rm 50}$
($R_{\rm 90}$), the radius enclosing 50% (90%) of the total mass of the dense
gas ($>2\rho_{\rm halo}$). As the cooling becomes inefficient, both $R_{\rm
50}$ and $R_{\rm 90}$ increase, i.e., the mass is more sparsely distributed.

The distribution of magnetic field lines also follows a similar trend. In
panel a-2, e-2, and f-2 of Fig. 5, we find that the field lines are much
tighter, denser, and better aligned in the tail of the MC model than in the
MC-Z0.03 model and the MC-NoCool model. Consequently, the field strength is
also overall higher in the MC model as it is given by the density of field
lines. The tighter field is due to the field lines being draped over the
smaller surface of the more compact MC clumps.

Figure 13: The average projected pressure (left) and $y$-velocity, i.e., the
velocity component along the magnetic fields which is perpendicular to the
cloud’s main motion (right) around the clouds in different models. The white
contour shows the density distribution.

The above results indicate that radiative cooling can change the overall flow
pattern in and around the clouds. Fig. 13 presents the gas pressure (left
column) and the $y$-velocity (right column) averaged along the $z$-axis. Note
that our choice to show the average properties along the projection axis
(z-axis) is to present the overall patterns of the pressure and the
$y$-velocity. When the properties across a slice (e.g., $z=0$ slice) is
considered, it is not as straightforward to compare the patterns between the
models because the clouds in different models have different structures. The
white contours show the density distribution as shown earlier in Fig. 5.

The pressure distribution in the MU model (panel a-1) shows a clear pressure
gradient in the front and behind the cloud as a result of the shock
compression due to the cloud’s high-velocity motion. The $y$-velocity field
(panel a-2) around the MU cloud has a flow pattern close to the quadrupolar
distribution: the flow of halo gas at $y>0$ and $y<0$ separate at the cloud’s
head and converge behind the cloud. The pressure and the $y$-velocity patterns
in the MC model (panel b-1 and b-2) are largely alike to that of the MU model,
but with a slight distortion due to the breached symmetry. Such similarity
indicates that the main cloud body of the MC model acts in a similar way to
the one in the MU model in determining the flow around the cloud. On the
contrary, we find noticeable differences in the flow pattern in models with
no/less cooling (panel c-1, 2, and panel d-1, 2). The diverging flow along the
$y$-axis at the cloud head is still noticeable, but the converging flow behind
the cloud is not as coherent as in the fiducial metallicity ($Z=0.3\,\rm
Z_{\odot}$) models. The presence of the low-pressure region behind the cloud
is also not as obvious.

At the initial snapshot, all three models (MC, MC-Z0.03, and MC-NoCool) start
with the same clumpy cloud density distribution; see the Clumpy-R32 cloud in
Fig. 2 and Fig. 3. Dense cores of the clumpy cloud are embedded in the inter-
clump gas that is not as dense as the cores and denser than the halo gas.
These intermediate-density materials are more susceptible to the ram pressure
stripping and therefore the halo gas can easily penetrate the cloud along the
narrow gap between the clumps. Such internal flow reduces the pressure
gradient and easily produces a disturbed flow pattern that continues to the
wake behind the clouds. In contrast, when radiative cooling is efficient,
inter-clump gas can easily cool down to a higher-density phase and act like
bonds between the existing cores. In this case, the halo gas faces stronger
resistance in flowing through the cloud, so much of the flow goes around the
cloud’s outer surface instead.

In summary, the role of radiative cooling in the cloud evolution is to (i)
produce dense clumps and filaments in the wake behind the cloud and (ii) keep
the main body more compact.

###  4.4 Density power spectrum (MC vs MC-p2.7 vs MC-p3.5)

The slope of the density spatial power spectrum ($\gamma$) describes the
relative strength of density fluctuations in different spatial scales. The
density and magnetic field distributions of the models at $t=28\,\rm Myr$ are
in panel (a), (g), (h)-1,2 of Fig. 5. Regardless of the power spectrum slope,
all the models develop the three leading head structures. However, the extent
to which the heads are exposed from the main cloud body varies: in the MC-p3.5
model all the heads are still tightly attached to the trailing main cloud
while in the MC-p2.7 model they are well separated from each other. Also,
another clear difference we find with the change of the power spectrum slope
is the amount of dense filamentary structures we see in the wake behind the
cloud.

The differences between the models come from the different density contrast
within the initial clumpy clouds. Earlier in Section 2.2, we compare the
structures of Clumpy-R32, Clumpy-p2.7, and Clumpy-p3.5 clouds that enter as
initial conditions of the MC model, the MC-p2.7 model, and the MC-p3.5 model,
respectively. In short, the shallower the slope of the power spectrum is, the
higher the density contrast within a cloud, i.e., the cloud has denser cores
and less dense inter-clump material. With time, such low-density inter-clump
gas can easily be pushed away from the cloud body. Therefore, we see dense
heads well separated from each other in the MC-p2.7 model compared to the
other models. Furthermore, the displaced inter-clump gas mixes with the hot
halo gas and condenses in the wake behind the cloud which explains why we find
more dense filaments and clumps in the tail region of the MC-p2.7 model
compared to the other models. This effect is also noticeable in the cold gas
mass evolution of the MC-p2.7 model (see Section 3 for detailed discussion).

##  5 Discussion

###  5.1 Resolution effect

Figure 14: The same format as Fig. 4, showing the evolution of models with
different spatial resolutions.

To assess the effect of the maximum resolution, we ran simulations identical
to the MC model, but with one less and one additional AMR level. This halves
(MC-R16) and doubles (MC-R64) the maximum resolution, respectively. As can be
seen in Fig. 14, the lower resolution simulation loses less cold gas and
enters the growth regime earlier than at higher resolution. It is overall more
similar to the uniform cloud case MU than at higher resolution. This is
expected because the lower resolution smoothes out some of the density
structures making the distribution slightly more uniform. We could only run
the higher resolution case until $t\sim 15\,\rm Myr$ due to its high
computational cost, however during this time it quite closely follows our
standard resolution MC model. It loses cold gas slightly faster than for the
MC model as expected given that the MC model also loses cold gas faster than
the lowest resolution case. Given this, we expect MC-R64 to also eventually
enter the growth regime but a bit later than the MC model. This is in line
with the findings for the resolution dependence of cold gas growth for
magnetized uniform clouds (Grønnow et al. 2018).

###  5.2 Caveats

Observational evidence suggests that the ISM as well as the surrounding hot
halo gas are turbulent (Jansson & Farrar 2012; Beck 2016). In this work, we
simulate a cloud with a realistic clumpy density structure motivated by
observed HI clouds. But in fact, not only is the clouds’ density distribution
clumpy, but their velocity and magnetic field structures are also turbulent. A
number of studies to date have taken into account the turbulent nature of the
cloud-halo interaction to some degree, e.g., internal turbulent velocity
fields of clouds (Heitsch et al. 2021), tangled magnetic fields within clouds
(McCourt et al. 2015; Gronke & Oh 2020) and turbulent magnetic fields in the
halo gas (Sparre et al. 2020), or all of the above (Banda-Barragán et al.
2018; Banda-Barragán et al. 2019). There is a clear indication that turbulent
halo magnetic fields act somewhat differently to uniform fields when draping a
cloud: they suppress Kelvin-Helmholtz instabilities in both the directions
that are perpendicular to the cloud’s motion (Sparre et al. 2020). That said,
at the very vicinity of the cloud-halo interaction front, random magnetic
fields are stretched and flattened out and therefore act like a transverse,
ordered field (Asai et al. 2007) similarly to the uniform halo magnetic field
adopted in the models in this study.

We do not include such non-uniform velocities or magnetic fields in models in
this study in order to isolate the effect of the clumpy density structure from
that of velocities or magnetic fields as well as to keep the number of free
parameters in the simulations manageable. In addition, the turbulent velocity
and tangled magnetic fields of HVCs are not as well identified by observations
as their density structure since these are more complicated and hard to
observe.

Also it is worth noting that many of the studies mentioned above typically aim
to simulate a cloud being launched from a galactic disc, which is clearly a
turbulent environment. The parameter space of our models represents HVCs at
larger distances infalling from the intergalactic space. Thus far, it is not
well studied to what extent such clouds are affected by the turbulent
environment.

In this work, we do not simulate the metallicity difference between the cloud
and the halo gas. We instead set a uniform metallicity throughout the
simulation domain as previously mentioned in Section 2. Miller & Bregman
(2013) estimated a lower limit of the hot-phase halo gas metallicity to be
$\sim 0.2\,\rm Z_{\odot}$ by comparing the electron density estimated from
their OVII K$\alpha$ absorption line measurements and the dispersion measure
of pulsars in the Large Magellanic Cloud (Anderson & Bregman 2010). Having
said that, using the same metallicity for the halo and HVC is well motivated
by the fact that the observed metallicities of the halo gas ($\sim 0.2\,\rm
Z_{\odot}$) and HVCs ($0.1-0.3\,\rm Z_{\odot}$) largely overlap.

Another physical process that may affect the evolution of cold clouds but
omitted in this work is thermal conduction. Hydrodynamic simulations in
general show that thermal conduction suppresses the growth of hydrodynamical
instaibilities at the cloud-halo interface (Vieser & Hensler 2007a, b;
Armillotta et al. 2016, 2017; Brüggen & Scannapieco 2016). With the presence
of magnetic fields, electrons, which are the main mediators of thermal
conduction, travel along the magnetic field lines that are not necessarily in
the direction that reduces the temperature gradient efficiently. Thus, thermal
evaporation of clouds and cloudlets is not as significant in our case, though
not negligible (see, e.g., Kooij et al. 2021).

Also, we do not attempt to directly model the heating and ionization of gas
from UV photons emitted by galaxies. However, the effect of heating is
indirectly roughly taken into account at low temperature by the $T\sim
10^{4}\,\rm K$ cooling floor as mentioned in Section 2. At temperatures around
$T\sim 10^{5}\,\rm K$, the inclusion of UV photons will reduce the cooling
efficiency and as a result, the formation of cold clumps condensed from the
cloud-halo mixing will be not as efficient as we see in models discussed in
this paper.

###  5.3 Milky Way HVCs

In this section, we apply our findings from the simulated clouds to explain
and predict observable properties of Milky Way HVCs.

####  5.3.1 Detectability on the RM sky

So far searches for magnetized HVCs has been mainly focused on identifying
coherent structure on the Faraday rotation measure (RM) grid that spatially
overlaps with the distribution of observed HI clouds on the sky-plane. For a
magnetized cloud to be able to leave a clear spatially resolved RM structure,
it is favorable to have the coherent magnetic field over a large-scale. The
more complex geometry of magnetic fields expected for clumpy clouds could
significantly restrict the observational detectability of magnetized HVCs.

Instead, the imprint of magnetized HVCs could appear as the small-scale
variation of RM. As the large-scale Galactic magnetic fields vary smoothly
over the sky plane, magnetized HVCs with complex field geometry would appear
as a small-scale distortion in the RM grid. The angular scale of such RM
variation would depend on various conditions, such as the distance to the
HVCs. However, the RM grid density provided by currently available
polarization data is not sufficient for such a statistical approach to the
small-scale RM variation in HVCs. We expect upcoming polarization surveys, for
example, the Polarisation Sky Survey of the Universe’s Magnetism
(POSSUM101010https://askap.org/possum/), will be capable of performing further
investigations with significantly improved RM source density.

####  5.3.2 Decelerated structures

We find significant deceleration of clumps behind the main body of the clouds
in our models. Our models do not include the effect of Milky Way’s
gravitational potential, but as discussed earlier in Section 3, Grønnow et al.
(2022) show that the most of the clumps remain decelerated due to the magnetic
tension even when the gravity is considered. The presence of the low-velocity
clumps suggests that, observationally, HVCs are not limited to what we find in
the high-velocity range but likely to be associated with lower-velocity
components. Indeed, Lockman et al. (2008) and Hill et al. (2013) identify
several decelerated clumps in HI and H$\alpha$ emission respectively
associated with the Smith cloud. evident in the position-velocity space. We
note that the models we use in this paper do not involve physics to
differentiate various ISM phases (e.g., ionization) but rather treat the gas
as a single phase fluid with different temperature. Therefore, the
investigation on the spatial distribution of neutral and ionized medium is
beyond the scope of this paper.

####  5.3.3 Growth of instabilities and small-scale structures

Earlier in Fig. 5, we have shown that clouds in our clumpy cloud models share
roughly common morphological features. They all develop multiple leading heads
attached to the main body and trailing filamentary tails regardless of the
environment or physical conditions. On top of that, the presence of the halo’s
magnetic field significantly suppresses the growth of small-scale structures
and leads to the development of elongated filamentary clumps. It is worth
noting that it is the individual clumps that we find are elongated, not the
entire cloud. The mass profile of the entire clouds are rather more extended
along the direction of the cloud’s motion in the HC model than in the MC
model.

Observationally, such multi-scale analysis requires targeted high-quality
observation data as the smallest detectable scale is restricted by the
resolution and the signal-to-noise ratio. Thus, the search for small-scale
structures across the cloud-wind interface has mainly targeted HVC complexes
with a large sky coverage. The high-resolution Herschel image of the Draco
nebula presented by Miville-Deschênes et al. (2017) shows multiple heads (or
“fingers” as they refer to them in the paper) linked to each other and
pointing toward the Galactic plane. Similarly to what we see in our models,
the separation between the dense heads is periodic as a result of the RT
instability. Barger et al. (2020) identified an array of HI clumps hanging
from the edges of Complex A. Their morphological and kinematical analysis
explains that these structures are a snapshot of active growth of hydrodynamic
instabilities, both KH and RT. Marchal et al. (2021) examined a small area of
HVC Complex C and found multi-scale and multi-phase structures in the region.
The identified clumps are predominantly elongated (the aspect ratio of
$1.5-6$) and have preferred position angle orientations. Such highly elongated
substructures favour the magnetized cloud in our models.

On the face of it, it is tempting to conclude that Complex A is not magnetized
and Complex C is magnetized. However, we caution against drawing a firm
conclusion here for a number of reasons: (i) it is difficult to separately
evaluate whether the observed small-scale structures of the clouds are
outcomes of RT and KH instabilities. (ii) None of Complex A and Complex C
reveal a sign of magnetic fields in the RM grid. We require studies on a
larger sample of HVCs for a comprehensive discussion on how the hydrodynamic
instabilities develop in observed clouds in comparison to our models.

####  5.3.4 The origin and fate of clouds

Metallicity is often considered as an indicator of diverse origins of the
HVCs. Clouds formed from recycled materials of the Galactic fountain are
expected to have close-to-solar metallicities (e.g., the Smith Cloud,
intermediate-velocity clouds) and ones with extragalactic origin have the
metallicity of primodial intergalactic medium if not enriched by external
galaxies (e.g., the Magellanic system). Models in this study show significant
metallicity dependency on the clouds’ morphology, the magnetic field strength,
and the evolution of the clouds’ cold gas mass. Our results imply that HVCs
with diverse origins would evolve differently under the same environment.
Fountain-originated clouds are more efficient in radiative cooling than
extragalactic-originated clouds which increases the chances of the clouds
eventually reaching the Galactic gas disk and fueling the star formation
activity. A caveat here is that the metallicity of a cloud varies with time.
Clouds and the surrounding halo gas can have different metallicity; usually
halo gas is expected to be more metal-poor than fountain-driven clouds. In
Heitsch et al. (2021), the authors show steady replacement of cloud materials
by the ambient halo gas when the mixing and the condensation is efficient, and
therefore, the metallicity composing a cloud constantly changes with time.
This implies that the metal mixing at the cloud-halo boundary can alter the
cooling efficiency over the course of the cloud evolution.

##  6 Conclusion

Understanding the growth and disruption of cold clouds moving at high-velocity
provides insight into whether they are capable of delivering cold gas to the
galaxy and fuel the star formation. Many earlier numerical simulations that
studied the survival of moving clouds in hot magnetized halo assume simple
spherically symmetric uniform density distributions as an initial condition,
whereas observed high-velocity clouds in similar environments clearly have
non-uniform clumpy structures. Motivated by this, we focus on performing
simulations of initially clumpy clouds and examining whether the change of the
cloud structure alters our understanding of the evolution of clouds in
magnetized halos.

The total cold gas mass in the uniform density clouds grow from the beginning
as expected from earlier studies, while the clumpy clouds lose cold gas early
on instead. Still, later on even clumpy clouds experience a significant growth
of the mass. Hence, the finding of earlier works that HVCs might be able to
grow and feed the Galaxy still holds to clumpy clouds.

In the fiducial model, some massive clumps to some degree have similar
characteristics with uniform density clouds. They have a head-tail morphology
and magnetic field lines draped around each of them. In comparison, these
characteristics are lost among clumps trailing behind the cloud. Instead, the
configuration of magnetic fields is dominated by a large-scale pattern that
drapes around the entire cloud system, not individual clumps. Clumps in this
regime are mostly condensed from the halo gas and significantly decelerated
from the high relative velocity of the main cloud.

By comparing morphologies of clouds evolved under different environments, we
study the significance of ram pressure stripping, growth of hydrodynamic
instabilities, and radiative cooling. Specifically, we compare models (i) with
and without halo magnetic fields, (ii) with varying metallicity and the
cooling efficiency, and (iii) with different slopes of the initial density
power spectrum.

  1. 1.

With magnetic fields present, there is an indication of suppressed
hydrodynamic instabilities at the cloud-halo interface along the direction of
the cloud’s motion. The magnetized cloud (MC model) has clumps that in general
have large effective size along the direction of motion. In comparison, its
non-magnetized counterpart (HC model) is composed of a larger number of
smaller clumps because clumps break into smaller pieces rather than become
elongated. Likewise, the 2D spatial power spectrum (SPS) of magnetized clouds
has less power at small length scales compared to that of non-magnetized
clouds throughout the time span covered in this work.

  2. 2.

Radiative cooling efficiency and the metallicity of the cloud, has two
apparent effects on the cloud structure. When cooling is efficient
(metallicity is higher), first, a greater number of dense clumps condense from
the halo gas in the wake behind a cloud. These clumps are substantially slower
with respect to the halo compared to the initial high-velocity motion of the
cloud. Second, the distributions of density and magnetic field lines are more
compact. This is because the intermediate density inter-clump material, which
would have been easily pushed away by the ram pressure when cooling was not
efficient, cools down to lower temperature (i.e., higher density) and acts
like a bond between the clumps. We confirm that this effect leads to a change
in the overall flow pattern around a clumpy cloud.

  3. 3.

The slope of the initial density spatial power spectrum defines the relative
strength between small- and large-scale structures in a clumpy cloud. A cloud
with a shallower power spectrum slope (MC-p2.7) has a wider density PDF, i.e.,
high density contrast between clumps and inter-clump gas compared to its
steeper slope counterpart (MC-p3.5). Early in the simulation, the total cold
gas mass in the cloud with the shallower slope decreases faster because the
ram pressure can easily push the low-density inter-clump gas. On the other
hand in the later stage of cloud evolution, such stripped intermediate density
gas mixes with the halo gas and cools down to form dense filamentary tails
behind the cloud resulting in significant growth of the cold gas mass. In
comparison, the cloud with a steeper power spectrum shows less significant
mass loss at the beginning followed by a more delayed upturn of mass.

## Data availability

The data directly related this article will be shared on reasonable request to
the corresponding author.

## Acknowledgements

We acknowledge high-performance computing resources provided by the Australian
National Computational Infrastructure (project ix81) in the framework of the
ANU Merit Allocation Scheme. N. M. Mc-G acknowledges funding from the
Australian Research Council in the form of DP190101571 and FL210100039.

Our analysis was performed using the Python programming language (Python
Software Foundation, https://www.python.org). The following packages were used
throughout the analysis: numpy (Harris et al. 2020), SciPy (Virtanen et al.
2020), and matplotlib (Hunter 2007). This research also made use of the
publicly available packages and softwares including VisIt (Childs et al.,
2012), yt project (Turk et al., 2011), powerbox (Murray, 2018), Turbustat
(Koch et al., 2019), CUPID (Berry et al., 2007).

S.L.J. and N. M. Mc-G acknowledge the Ngunnawal and Ngambri people as the
traditional owners and ongoing custodians of the land on which the Research
School of Astronomy & Astrophysics is sited at Mt Stromlo.

## References

  * Anderson & Bregman (2010) Anderson M. E., Bregman J. N., 2010, ApJ, 714, 320
  * Armillotta et al. (2016) Armillotta L., Fraternali F., Marinacci F., 2016, MNRAS, 462, 4157
  * Armillotta et al. (2017) Armillotta L., Fraternali F., Werk J. K., Prochaska J. X., Marinacci F., 2017, MNRAS, 470, 114
  * Asai et al. (2007) Asai N., Fukuda N., Matsumoto R., 2007, ApJ, 663, 816
  * Banda-Barragán et al. (2016) Banda-Barragán W. E., Parkin E. R., Federrath C., Crocker R. M., Bicknell G. V., 2016, MNRAS, 455, 1309
  * Banda-Barragán et al. (2018) Banda-Barragán W. E., Federrath C., Crocker R. M., Bicknell G. V., 2018, MNRAS, 473, 3454
  * Banda-Barragán et al. (2019) Banda-Barragán W. E., Zertuche F. J., Federrath C., García Del Valle J., Brüggen M., Wagner A. Y., 2019, MNRAS, 486, 4526
  * Barger et al. (2020) Barger K. A., et al., 2020, ApJ, 902, 154
  * Beck (2016) Beck R., 2016, A&ARv, 24, 7
  * Berry (2015) Berry D. S., 2015, Astronomy and Computing, 10, 22
  * Berry et al. (2007) Berry D. S., Reinhold K., Jenness T., Economou F., 2007, in Shaw R. A., Hill F., Bell D. J., eds, Astronomical Society of the Pacific Conference Series Vol. 376, Astronomical Data Analysis Software and Systems XVI. p. 425
  * Brüggen & Scannapieco (2016) Brüggen M., Scannapieco E., 2016, ApJ, 822, 31
  * Chandrasekhar (1961) Chandrasekhar S., 1961, Hydrodynamic and hydromagnetic stability
  * Childs et al. (2012) Childs H., et al., 2012, in , High Performance Visualization–Enabling Extreme-Scale Scientific Insight. pp 357–372, doi:10.1201/b12985
  * Clark et al. (2021) Clark S., Bordoloi R., Fox A. J., 2021, arXiv e-prints, p. arXiv:2108.04295
  * Cooper et al. (2009) Cooper J. L., Bicknell G. V., Sutherland R. S., Bland-Hawthorn J., 2009, ApJ, 703, 330
  * Cottle et al. (2020) Cottle J., Scannapieco E., Brüggen M., Banda-Barragán W., Federrath C., 2020, ApJ, 892, 59
  * Dedner et al. (2002) Dedner A., Kemm F., Kröner D., Munz C. D., Schnitzer T., Wesenberg M., 2002, Journal of Computational Physics, 175, 645
  * Dursi & Pfrommer (2008) Dursi L. J., Pfrommer C., 2008, ApJ, 677, 993
  * For et al. (2013) For B.-Q., Staveley-Smith L., McClure-Griffiths N. M., 2013, ApJ, 764, 74
  * Fox et al. (2019) Fox A. J., Richter P., Ashley T., Heckman T. M., Lehner N., Werk J. K., Bordoloi R., Peeples M. S., 2019, ApJ, 884, 53
  * Gregori et al. (1999) Gregori G., Miniati F., Ryu D., Jones T. W., 1999, ApJ, 527, L113
  * Gregori et al. (2000) Gregori G., Miniati F., Ryu D., Jones T. W., 2000, ApJ, 543, 775
  * Gronke & Oh (2018) Gronke M., Oh S. P., 2018, MNRAS, 480, L111
  * Gronke & Oh (2020) Gronke M., Oh S. P., 2020, MNRAS, 492, 1970
  * Grønnow et al. (2017) Grønnow A., Tepper-García T., Bland -Hawthorn J., McClure-Griffiths N. M., 2017, ApJ, 845, 69
  * Grønnow et al. (2018) Grønnow A., Tepper-García T., Bland -Hawthorn J., 2018, ApJ, 865, 64
  * Grønnow et al. (2022) Grønnow A., Tepper-García T., Bland-Hawthorn J., Fraternali F., 2022, MNRAS, 509, 5756
  * Harris et al. (2020) Harris C. R., et al., 2020, Nature, 585, 357–362
  * Heitsch et al. (2021) Heitsch F., Marchal A., Miville-Deschênes M. A., Shull J. M., Fox A. J., 2021, MNRAS,
  * Henley & Shelton (2015) Henley D. B., Shelton R. L., 2015, ApJ, 808, 22
  * Hill et al. (2013) Hill A. S., Mao S. A., Benjamin R. A., Lockman F. J., McClure-Griffiths N. M., 2013, ApJ, 777, 55
  * Hunter (2007) Hunter J. D., 2007, Computing in Science and Engineering, 9, 90
  * Jansson & Farrar (2012) Jansson R., Farrar G. R., 2012, ApJ, 761, L11
  * Jones et al. (1996) Jones T. W., Ryu D., Tregillis I. L., 1996, ApJ, 473, 365
  * Koch et al. (2019) Koch E. W., Rosolowsky E. W., Boyden R. D., Burkhart B., Ginsburg A., Loeppky J. L., Offner S. S. R., 2019, TurbuStat: Turbulence statistics in spectral-line data cubes (ascl:1907.015)
  * Kooij et al. (2021) Kooij R., Grønnow A., Fraternali F., 2021, MNRAS, 502, 1263
  * Kwak et al. (2009) Kwak K., Shelton R. L., Raley E. A., 2009, ApJ, 699, 1775
  * Li et al. (2020) Li Z., Hopkins P. F., Squire J., Hummels C., 2020, MNRAS, 492, 1841
  * Lockman et al. (2002) Lockman F. J., Murphy E. M., Petty-Powell S., Urick V. J., 2002, ApJS, 140, 331
  * Lockman et al. (2008) Lockman F. J., Benjamin R. A., Heroux A. J., Langston G. I., 2008, ApJ, 679, L21
  * Marchal et al. (2021) Marchal A., Martin P. G., Gong M., 2021, arXiv e-prints, p. arXiv:2106.12683
  * Marinacci et al. (2010) Marinacci F., Binney J., Fraternali F., Nipoti C., Ciotti L., Londrillo P., 2010, MNRAS, 404, 1464
  * McCourt et al. (2015) McCourt M., O’Leary R. M., Madigan A.-M., Quataert E., 2015, MNRAS, 449, 2
  * Mignone & Tzeferacos (2010) Mignone A., Tzeferacos P., 2010, Journal of Computational Physics, 229, 2117
  * Mignone et al. (2010) Mignone A., Tzeferacos P., Bodo G., 2010, Journal of Computational Physics, 229, 5896
  * Mignone et al. (2012) Mignone A., Zanni C., Tzeferacos P., van Straalen B., Colella P., Bodo G., 2012, ApJS, 198, 7
  * Miller & Bregman (2013) Miller M. J., Bregman J. N., 2013, ApJ, 770, 118
  * Miller & Bregman (2015) Miller M. J., Bregman J. N., 2015, ApJ, 800, 14
  * Miville-Deschênes et al. (2017) Miville-Deschênes M. A., et al., 2017, A&A, 599, A109
  * Moss et al. (2013) Moss V. A., McClure-Griffiths N. M., Murphy T., Pisano D. J., Kummerfeld J. K., Curran J. R., 2013, ApJS, 209, 12
  * Murphy et al. (1995) Murphy E. M., Lockman F. J., Savage B. D., 1995, ApJ, 447, 642
  * Murray (2018) Murray S. G., 2018, The Journal of Open Source Software, 3, 850
  * Patnaude & Fesen (2005) Patnaude D. J., Fesen R. A., 2005, ApJ, 633, 240
  * Pingel et al. (2018) Pingel N. M., Lee M.-Y., Burkhart B., Stanimirović S., 2018, ApJ, 856, 136
  * Putman et al. (2012) Putman M. E., Peek J. E. G., Joung M. R., 2012, ARA&A, 50, 491
  * Richter (2012) Richter P., 2012, ApJ, 750, 165
  * Richter (2017) Richter P., 2017, in Fox A., Davé R., eds, Astrophysics and Space Science Library Vol. 430, Gas Accretion onto Galaxies. p. 15 (arXiv:1612.00449), doi:10.1007/978-3-319-52512-9_2
  * Santillán et al. (1999) Santillán A., Franco J., Martos M., Kim J., 1999, ApJ, 515, 657
  * Schneider & Robertson (2017) Schneider E. E., Robertson B. E., 2017, ApJ, 834, 144
  * Shin et al. (2008) Shin M.-S., Stone J. M., Snyder G. F., 2008, ApJ, 680, 336
  * Sparre et al. (2020) Sparre M., Pfrommer C., Ehlert K., 2020, MNRAS, 499, 4261
  * Sutherland & Dopita (1993) Sutherland R. S., Dopita M. A., 1993, ApJS, 88, 253
  * Thom et al. (2008) Thom C., Peek J. E. G., Putman M. E., Heiles C., Peek K. M. G., Wilhelm R., 2008, ApJ, 684, 364
  * Toro et al. (1994) Toro E. F., Spruce M., Speares W., 1994, Shock Waves, 4, 25
  * Turk et al. (2011) Turk M. J., Smith B. D., Oishi J. S., Skory S., Skillman S. W., Abel T., Norman M. L., 2011, The Astrophysical Journal Supplement Series, 192, 9
  * Vieser & Hensler (2007a) Vieser W., Hensler G., 2007a, A&A, 472, 141
  * Vieser & Hensler (2007b) Vieser W., Hensler G., 2007b, A&A, 475, 251
  * Virtanen et al. (2020) Virtanen P., et al., 2020, Nature Methods, 17, 261
  * Wakker (2001) Wakker B. P., 2001, ApJS, 136, 463
  * Westmeier (2018) Westmeier T., 2018, MNRAS, 474, 289