##  1 Introduction

Coronal mass ejections (CMEs) are plasma and magnetic field structures that
usually erupt from highly magnetic, closed-field regions of the Sun called
active regions (ARs). CMEs can carry large kinetic and magnetic energy with
their erupted mass ranging between $10^{9}$ kg and $10^{13}$ kg [Vourlidas et
al., 2010], and speeds that may exceed 3000 km/s. An Earth-directed CME
carrying a negative out-of-the-ecliptic magnetic field component
($B_{\mathrm{z}}$) can cause extreme space weather events that can affect
space-borne and ground-based technological systems. Therefore, predicting CME
arrival times and magnetic field values at Earth by properly treating their
interplanetary propagation is of immense importance.

A number of CME propagation models have been developed over the years to
reproduce their travel through interplanetary space and forecast arrival times
at Earth. Such models include empirical models [e.g. Vandas et al., 1996,
Brueckner et al., 1998, Gopalswamy et al., 2001, 2005, Wang et al., 2002],
drag based models [Vršnak & Gopalswamy, 2002], and physics-based models, such
as the Shock Time of Arrival (STOA) [Dryer & Smart, 1984], STOA-2 [Moon et
al., 2002], and HAFv.2 [Fry et al., 2001]. Newer-generation prediction models
such as Wang-Sheeley-Arge (WSA)-ENLIL-Cone model [Odstrcil et al., 2004] use
magnetohydrodynamic (MHD) simulations in which CMEs are inserted into the
ambient flow as non-magnetized over-pressured plasma blobs. The accuracy of
this type of models is not much higher than in the previously-mentioned models
of lesser complexity and have a mean absolute error of more than 10 hr [Riley
et al., 2018, Amerstorfer et al., 2018]. Gopalswamy et al. [2013] have shown
that the performance of the ENLIL-cone model is similar to that of the
empirical shock arrival model [Gopalswamy et al., 2005]. The ENLIL-cone and
HAFv.2 models can calculate magnetic field in the ICME sheath region because a
preceeding shock distorts the interplanetary magnetic field. However, they are
unable to describe magnetic field in ICME clouds at Earth because of the
absence of an explicitly defined magnetic field inside them [Fry et al., 2001,
Odstrcil et al., 2004]. Attempts have been made to remedy this by adding a
simple magnetic field into the plasma blobs [see, e.g., Shen et al., 2011].

According to the standard model of CME eruption, also called the CSHKP model
[Carmichael, 1964, Sturrock, 1968, Hirayama, 1974, Kopp & Pneuman, 1976], pre-
eruptive states are characterized by the presence of an axial magnetic flux,
which may be either a sheared arcade or a flux rope located above the magnetic
polarity inversion line separating positive and negative radial magnetic
fluxes in an active region. Such flux ropes can rise due to different magnetic
force imbalances, and the overlying magnetic field lines add the majority of
poloidal flux to the flux rope during an eruption [see the review of Chen,
2011]. The interplanetary counterparts of CME flux ropes which are
characterized by a large, rotating magnetic field, as revealed by in-situ
measurements, are called the magnetic clouds [Burlaga et al., 1981].

Using in-situ observations of ICMEs between 1978–1982, Gosling [1990] showed
that $\sim 30\%$ of ICMEs sampled at 1 AU had magnetic clouds associated with
them. However, Marubashi [2000] later claimed that up to $80\%$ of all ICMEs
involve magnetic clouds. A recent study by Song et al. [2020] argues that all
CMEs observed at 1 AU have magnetic clouds, though some are not directly
observed because the observing spacecraft crosses the flanks of the ICME.
Studies, such as Gopalswamy et al. [2018] and Sarkar et al. [2020] show the
key importance of flux ropes in the development of forecast models capable of
predicting the CME magnetic properties at 1 AU. Therefore, the new generation
of CME forecast models is strongly associated with the advanced description of
flux rope structures [Shen et al., 2011, Shiota & Kataoka, 2016, Isavnin,
2016, Vandas & Romashets, 2017, Jin et al., 2017, Scolini et al., 2019, Singh
et al., 2020].

MHD simulations of the solar wind (SW) and flux-rope-based CMEs can either
start in the lower corona [Manchester et al., 2004, Jin et al., 2016, Singh et
al., 2018] or the inner heliosphere (IH) [Verbeke et al., 2019, Singh et al.,
2020]. If a simulation starts in the lower corona, one can allow a CME to be
automatically generated in the ambient solar corona background [Amari et al.,
2014]. Otherwise, an analytically-generated flux rope can be inserted [Titov &
Démoulin, 1999, Gibson & Low, 1998, Singh et al., 2020]. The latter approach
is less computationally expensive than the former and allows for better
constraining of magnetic flux inside flux-rope models [Singh et al., 2019].
The approach where CMEs are inserted in the IH is computationally less
expensive and easier to implement since the boundary condition treatment in
the superfast magnetosonic flow is mathematically trivial and no coronal
heating model is required. However, this approach can lead to errors
associated with the extrapolation of the SW and CMEs from the lower corona to
the inner boundary of the IH models. This approach is more common to the
current operational ICME predictions such as by NASA Goddard Space Flight
Center/Space Weather Research Center, NOAA/Space Weather Prediction Center
USA, Australian Bureau of Meteorology, Korean Space Weather Center, and UK Met
Office, which use the WSA-ENLIL-Cone model. This approach is also being
pursued by several other models, e.g., EUHFORIA [Verbeke et al., 2019] and the
one described by Shiota & Kataoka [2016].

In an IH model, a CME can either be inserted step by step at the inner
boundary [Verbeke et al., 2019] or in a single step above the inner boundary
[Singh et al., 2020]. If a CME flux rope is inserted by superimposing an
analytic model either in the corona or in the IH, its physical and magnetic
properties should be driven by observations for accurate forecasting. These
properties include the CME direction, tilt, shape, speed, mass, poloidal and
toroidal magnetic fluxes, and helicity sign. The CME tilt refers here to the
angle between the line connecting two footpoints of the CME flux rope and the
equatorial plane. This tilt can be seen as angle $\gamma$ in Figure 1 of
Thernisien et al. [2009]. Moreover, the ambient SW through which this CME will
propagate should also be based on data-driven models. CME observations used to
constrain any model are typically accompanied by uncertainties in different
measured quantities. It is crucial to understand how these uncertainties
propagate through the computational region and affect the model forecasts at
Earth. Ensemble modeling of CMEs is a viable approach to the study of the
propagation of these uncertainties and ensure the probabilistic forecasting
[Pizzo et al., 2015, Mays et al., 2015, Murray, 2018, Amerstorfer et al.,
2021].

In this paper, we use the constant-turn flux rope model to simulate a CME that
erupted on 12 July 2012 at 16:48. The geometry of this model is described by
the FRiED model of Isavnin [2016]. We use the analytic formulae governing the
constant turn flux rope, as described by Vandas & Romashets [2017], to specify
the magnetic field in this model. We estimate the 12 July 2012 CME direction,
tilt, half-angle, aspect ratio, and speed using the graduated cylindrical
shell (GCS) model and initialize our flux rope using these properties. The
magnetic flux associated with a CME in the flux-rope model is calculated using
the poloidal flux estimated from post-eruption arcades (PEAs).

We demonstrate that by using these observational initial conditions it is
possible to correctly reproduce the sign of all three magnetic field
components of the magnetic cloud of this ICME. The arrival time, density, and
speed are also reproduced with reasonably good accuracy as compared to
previous attempts to simulate this CME using different models [Shen et al.,
2014, Scolini et al., 2019, Singh et al., 2020]. We also perform an analysis
of differences in GCS fittings results reported in multiple publications and
online catalogs to quantify the uncertainties associated with such fittings.
We further investigate how the uncertainties in the GCS parameters, both
individually and collectively, affect the simulated CME properties at Earth by
using an ensemble modeling of the 12 July 2012 CME.

In Section 2, we describe the GCS model, SW model, and the flux-rope-based CME
model. We then show our CME simulation and GCS uncertainty analysis results in
Section 3, followed by the conclusions in Section 4.

##  2 Data and Models

In this work, we use the level 0.5 image data collected on 12 July 2012
between 16:24 and 18:54 UT from the Sun-Earth Connection Coronal and
Heliospheric Investigation (SECCHI)/Cor2 [Howard et al., 2008] coronagraphs on
board the Solar Terrestrial Relations Observatory (STEREO) A & B [Kaiser et
al., 2008] spacecraft and process them to level 1 using the secchi_prep
program in the IDL SolarSoft library, which converts the units from data
numbers (DN) into Mean Solar Brightness (MSB). We also use level 1 image data
collected on 12 July 2012 between 16:24 and 17:18 UT from the Large Angle
Spectroscopic Coronagraph (LASCO)/C2/C3 [Brueckner et al., 1995] coronagraph
of the Solar and Heliospheric Observatory (SOHO) spacecraft. Level 1 data for
both these instruments are suitable for quantitative scientific analysis. The
coronagraph data are used to estimate the speed, direction, tilt, half-angle,
aspect ratio, and speed of the 12 July 2012 CME. The model used to estimate
these values is described in Section 2.1. We use the solar extreme ultraviolet
(EUV) images in 94, 131, and 193 Å wavelengths captured on 12 July 2012 at
22:30 UT and line of sight (LOS) magnetic field observations at 12 July 2012
16:10 UT from Solar Dynamics Observatory’s (SDO) Atmospheric Imaging Assembly
(AIA) [Pesnell et al., 2012] and Helioseismic and Magnetic Imager (HMI) [Schou
et al., 2012, Hoeksema et al., 2014], respectively. We use the aia_prep
program in SolarSoft to read in and calibrate AIA level 1 data to level 1.5.
Level 1 data include bad-pixel removal, de-spiking, and flat-fielding.
Additional roll-corrections applied to Level 1.5 data make the solar north
direction vertical in the images, re-scale the images to $0^{\prime\prime}.6$
pixels, and translate them to match the solar disk and image centers. The HMI
data are read using hmi_prep program in SolarSoft, which processes
magnetograms to ensure correct roll angles, and translates them to put the
solar disk center at a magnetogram center. These data are used to estimate the
magnetic flux and helicity sign of the 12 July 2012 CME. The methods used to
estimate these values are described in Section 3.1. We use 1 hour averaged
data provided by NASA/GSFC’s OMNI data set through OMNIWeb [King &
Papitashvili, 2005] to compare the in-situ properties at Earth with our
simulations of 12 July 2012 CME.

###  2.1 Graduated Cylindrical Shell model

The Graduated Cylindrical Shell (GCS) model uses three different views of a
CME from STEREO -A & -B, and SOHO coronagraphs to visually fit the observed
CME with a simplified structure consisting of conical legs and a curved front,
resembling a hollow croissant with a circular cross-section. The GCS method
has been implemented in IDL using the rtsccguicloud program [Thernisien et
al., 2006]. An example of GCS fitting using this tool is shown in Figure 1,
for the July 12, 2012 CME. The GCS model has six free parameters, namely, the
latitude, longitude, tilt, height, half-angle, and aspect ratio, which can be
modified to fit the 3D flux rope structure in all three images simultaneously.
The rtsccguicloud program allows for this modification by using sliders to
modify each parameter independently. The mathematical description of the
geometry of this model is described in detail by Thernisien [2011]. Notice
that the CME height, being one of the GCS parameters, can be converted to
speed by fitting a linear function to the height-time graph. We decided to use
linear fitting instead of quadratic fitting because the GCS uncertainty
analysis we perform later in the study is based on the catalogs that all
reported only the linear speeds. One cannot reliably derive GCS parameters
with a single coronagraph image, since the majority of structural information
of a CME is lost when its 3D shape is projected onto a 2D plane (plane of the
sky). The availability of multiple images made from different viewpoints is
critical for our ability to remove these projection effects. However, CMEs
typically have complicated shapes that are difficult to fit with the idealized
shapes assumed by the GCS model. Therefore, the CME parameters estimated by
the GCS model can have significant uncertainties. We will discuss them in more
detail in Section 3.2.

Figure 1: (Top panel, from left to right)The July 12, 2012 CME seen in STEREO
B Cor2, SOHO C2 and STEREO A Cor2 coronagraphs respectively. (Bottom panel)
The same CME fitted with the GCS model.

###  2.2 Magnetohydrodynamic models

In this study, we perform MHD modeling using the Multi-Scale Fluid-Kinetic
Simulation Suite (MS-FLUKSS), a collection of highly parallel modules capable
of performing adaptive mesh refinement (AMR) simulations of the SW in the
presence of neutral atoms, nonthermal ions, and turbulence, etc. [Pogorelov et
al., 2014]. We describe our SW and CME models in the following subsections.

####  2.2.1 Solar wind model

We simulate the ambient SW in the IH by solving the ideal MHD equations with
the finite-volume, total variation diminishing (TVD) approximation on a non-
uniform spherical grid with 150, 256, and 128 cells in the $r$, $\phi$, and
$\theta$ direction, respectively. We set the inner and outer boundaries of the
spherical domain at 0.1 and 1.5 AU, respectively. The IH model uses a time-
series of WSA maps as inner boundary conditions at 0.1 AU [e.g., Kim et al.,
2020]. We employ the Air Force Data Assimilative Photospheric Flux Transport
(ADAPT) synchronic maps using the NSO/GONG magnetograms [Arge et al., 2010,
2011, 2013, Hickmann et al., 2015] as input to the WSA model at the solar
surface. The WSA model extrapolates the photospheric magnetic field to a
spherical source surface at 2.5 $R_{\odot}$ using the potential field source
surface (PFSS) model, and then to the WSA outer boundary at 0.1 AU using the
Schatten current sheet model [Schatten, 1971]. The WSA model also calculates
the SW speed at its outer boundary as a function of the flux expansion factor
and distance to the nearest coronal hole boundary [Arge et al., 2003, 2005,
McGregor et al., 2011]. To generate the background SW for the flux-rope CME,
we use the best ADAPT-WSA realization (out of 12), based on their comparison
with near-Earth SW data.

We interpolate the WSA maps, both in time and space, from the original
$2.5^{\circ}\times 2.5^{\circ}$ resolution to the IH model base grid of $\sim
1.4^{\circ}\times 1.4^{\circ}$. We also scale the WSA magnetic field strengths
by a factor of 3 during this process to compensate for the systematic
underestimation of the open magnetic flux at 1 AU [Linker et al., 2016, 2017,
Wallace et al., 2019]. This factor has been selected to improve the match of
observations and simulations at Earth. We have used the same factor in Singh
et al. [2020]. We estimate the radial and azimuthal components of the magnetic
field at 0.1 AU from the WSA B values using the local SW speed to account for
the Sun’s rotation because the SW propagates radially outward from 1
$R_{\odot}$ to 0.1 AU [MacNeice et al., 2011]. We reduce the WSA speeds by 20%
to account for the difference in SW acceleration between the WSA and MS-FLUKSS
models [e.g., MacNeice et al., 2011, Kim et al., 2014]. Kim et al. [2014]
showed that using this reduction factor results in an improved matching of
observed and simulated parameters at Earth. To estimate the SW density and
temperature at 0.1 AU, we use the empirical correlations between the SW speed,
density, and temperature based on OMNI data [Elliott et al., 2016]. These
procedures give us all inner boundary conditions required to simulate SW in
the IH.

####  2.2.2 CME model

In this work, to simulate flux-rope-based CMEs, we use the geometry of the
FRiED model [Isavnin, 2016] and describe the initial magnetic field in it in
accordance with the uniform-twist analytic solution described by Vandas &
Romashets [2017]. The geometry of the FRiED model simplifies the CME shape to
a croissant-like structure with two legs rooted at the center of the Sun, as
shown in the left panel of Figure 2. This geometry is characterized by the
parameters $R_{\mathrm{t}}$, $R_{\mathrm{p}}$, and $\phi_{\mathrm{hw}}$, as
shown in the figure. The flux rope has a circular cross-section with a curved
axis shown by the dashed black line. The radius of this cross-section,
$R(\phi)$, is specified using the formula

$R(\phi)=\frac{R_{\mathrm{p}}}{R_{\mathrm{t}}}r(\phi),$ (1)

where $r(\phi)$ is the distance of the curved axis from the origin.

This distance is defined by the formula

$r(\phi)=R_{\mathrm{t}}\textup{cos}^{n}\left(\frac{\pi}{2}\frac{\phi}{\phi_{\mathrm{hw}}}\right),$
(2)

where, $n$ is a free parameter which can be used to adjust the flatness of the
flux rope shape.

We found that the flux rope geometry of the FRiED model can be matched very
well with the GCS model, designed with any suitable choice of its parameters,
by appropriately selecting the FRiED parameters. For example, the outer
boundary of a GCS-shaped flux rope with half-width
$hw_{\mathrm{GCS}}=33^{\circ}$, height $H_{\mathrm{GCS}}=70\,R_{\odot}$, and
aspect ratio $\kappa=0.5$ matches very well with the FRiED flux rope with
$R_{\mathrm{t}}=46.67\,R_{\odot}$, $R_{\mathrm{p}}=23.33\,R_{\odot}$, half-
width $\phi_{\mathrm{hw}}=40^{\circ}$, and $n=0.33$, as shown in the right
panel of Figure 2. The choice of parameters for these shapes is suitable for
the 12 July 2012 CME, which we are going to study below.

We use the uniform twist solution of Vandas & Romashets [2017] because some
recent studies [e.g., Hu et al., 2015] show that in situ signatures of the
ICME magnetic field are consistent with a constant turn configuration. Since
the FRiED model has a curved geometry, we cannot use constant turn cylindrical
flux rope models such as Gold–Hoyle model [Gold & Hoyle, 1960] to describe the
magnetic field in it. However, a magnetic field solution in a curved shape
such as a torus is more compatible with FRiED model geometry. Vandas &
Romashets [2017] give the analytic solution of a constant turn magnetic field
in a torus shape as:

$B_{\mathrm{r}}=0;\\\ B_{\phi}=\frac{B_{\mathrm{0}}}{1+b^{2}r^{2}};\\\
B_{\theta}=-\frac{B_{\mathrm{0}}R_{\mathrm{0}}br}{(1+b^{2}r^{2})(R_{\mathrm{0}}+r\cos\theta)}.$
(3)

Here $r$, $\theta$, and $\phi$ are the coordinates of the toroidally curved
cylindrical system, which is described in detail in Section 2 of Vandas &
Romashets [2017]. The parameters $B_{\mathrm{0}}$ and $b$ control the magnetic
flux and the number of magnetic field line turns, respectively.
$R_{\mathrm{0}}$ is the major radius of the torus. Since the FRiED geometry is
not exactly a torus, we define the magnetic field distribution in the
structure by assuming it to be locally toroidal. The parameter
$B_{\mathrm{0}}$ is modified inside the flux rope to conserve the inserted
axial (toroidal) flux. The value of $b$ can be estimated from the inserted
poloidal flux inside the flux rope. We fix $R_{\mathrm{0}}$ throughout the
flux rope so that the introduced positive and negative poloidal fluxes are
within 20% of each other. This means that the magnetic flux directed into the
plane containing the curved axis of the flux rope is within 20% of the
magnetic flux directed out of this plane.

In this model, we define the initial velocity at any location inside the flux
rope as a combination of the radial and expansion velocities, i.e.,
$\vec{V}=\vec{V}_{\mathrm{rad}}+\vec{V}_{\mathrm{exp}}$, where
$\vec{V}_{\mathrm{rad}}$ is in the radial direction away from the origin (the
Sun’s center) and $\vec{V}_{\mathrm{exp}}$ is in the direction pointing away
from the curved axis. The speed at the apex is considered to be the CME speed
$V_{\mathrm{CME}}$. To facilitate a self-similar expansion, we specify

$|V_{\mathrm{rad}}|=\frac{V_{\mathrm{CME}}}{1+R_{\mathrm{p}}/R_{\mathrm{t}}},\quad|V_{\mathrm{exp}}(r_{\mathrm{p}})|=\frac{r_{\mathrm{p}}}{R_{\mathrm{t}}}|V_{\mathrm{rad}}|,$

where $r_{\mathrm{p}}$ is the radial coordinate inside the local torus
structure. This specification can be readily derived from the self-similar
properties $R_{\mathrm{p}}\propto R_{\mathrm{t}}$ and $r_{\mathrm{p}}\propto
R(\phi)$. Initially, we assume a constant density inside the flux rope.

---
Figure 2: (Left panel) A cross-section of constant-turn flux rope by the plane
containing the curved axis of this flux rope (shown here by the dashed black
line). Blue lines represent the inner and outer edges of the flux rope in this
plane. The variable $r(\phi)$ is a polar coordinate describing the curved
axis. $R(\phi)$ is the radius of the out-of-plane circular cross section
(shown with orange color) centered at a point belonging to the curved axis,
while $r_{p}$ is the distance of a point on this circular cross section from
its center. $R(\phi)$ and $r(\phi)$ are given by Equations 1 and 2
respectively. (Right panel) The outer edge of the FRiED model with model
parameters $R_{\mathrm{t}}=46.67R_{\odot}$, $R_{\mathrm{p}}=23.33R_{\odot}$,
$\phi_{\mathrm{hw}}=40^{\circ}$, and $n=0.33$, is compared with the outer edge
of the GCS model with model parameters $\textup{half width}=33^{\circ}$,
$\textup{height}=70R_{\odot}$ and $\textup{aspect ratio}=0.5$.

We insert this flux rope into the ambient SW in such a way that the flux rope
is initially superimposed with the SW background. This superimposition is
described in detail in Singh et al. [2020]. It should be noted that this
approach is different from the other, more commonly used approach, where a
flux rope is inserted at the inner boundary gradually [e.g., Shiota & Kataoka,
2016, Verbeke et al., 2019]. Our method of flux-rope insertion into the SW
preserves the CME flux rope much more accurately, exhibiting a curved front
and two legs. On the other hand, the gradual, step-by-step insertion of flux
ropes at the inner boundary can be problematic because of their unnatural
initial expansion, which must then be compensated by introducing the flux rope
with initial speeds smaller than observed. We introduce the flux rope
parameters as follows:

  * •

$\vec{B}_{\mathrm{final}}=\vec{B}_{\mathrm{FR}}$,

  * •

$\rho_{\mathrm{final}}=\rho_{\mathrm{FR}}+\rho_{\mathrm{SW}}$,

  * •

$e_{\mathrm{final}}=e_{\mathrm{FR}}+e_{\mathrm{SW}}$.

Here $\vec{{B}}_{\mathrm{FR}}$ is given by Eq. (3), $\rho$ is the plasma
density and $e$ is the total energy density. We also note that

$e_{\mathrm{FR}}=\frac{|\vec{B}_{\mathrm{FR}}|^{2}}{8\pi},\quad
e_{\mathrm{SW}}=\frac{p_{\mathrm{SW}}}{\gamma-1}+\frac{|\vec{B}_{\mathrm{SW}}|^{2}}{8\pi}+\frac{\rho_{\mathrm{SW}}|\vec{v}_{\mathrm{SW}}|^{2}}{2},$

where $p$ and $\vec{v}$ are the thermal pressure and bulk velocity,
respectively. We have kept the adiabatic index $\gamma=1.5$ in this study.

In Figure 3, we show an example of such insertion. The flux rope has the
following FRiED parameters: $R_{\mathrm{t}}=46.67$, $R_{\mathrm{p}}=23.33$,
$\phi_{\mathrm{hw}}=40^{\circ}$, and $n=0.33$. The flux rope is assigned an
initial speed of 1100 $km/s$, a poloidal flux of $1.4\times 10^{22}$ Mx, and a
toroidal flux of $7.6\times 10^{21}$ Mx. It is initialized with zero latitude,
longitude, and tilt. The helicity sign is kept positive. A flux rope with
these properties is chosen because we will later simulate the 12 July 2012 CME
with the same parameters. The left panel of Figure 3 shows the velocity
distribution in the $z=0$ slice of this flux rope, along with the magnetic
field lines inside it. The right panel of Figure 3 shows the distribution of
$B_{\mathrm{z}}$ in the $z=0$ plane. The gray sphere represents the inner
boundary of our simulation domain at 0.1 AU. The SW background shown here is
on 12 July 2012 at 21:30 UT.

---
Figure 3: The initial structure of the constant turn flux rope model with the
apex speed of 1100 km/s, poloidal flux of $1.4\times 10^{22}$ Mx, and toroidal
flux of $7.6\times 10^{21}$ Mx when inserted into the inner heliosphere. The
apex of the flux rope is at $70R_{\odot}$, while its latitude, longitude, and
tilt zero. The $z=0$ semi-translucent slices of (left panel) radial velocity
and (right panel) $B_{\mathrm{z}}$ are shown along with the magnetic field
lines. The gray sphere represents the inner boundary of the IH model at 0.1
AU.

After the initial insertion, the model flux rope propagates through the IH as
an ICME. In Figure 4, We show the model ICME after 30 hrs of propagation. The
apex speed of the ICME has dropped to about 700 $km/s$ by this time and the
$B_{\mathrm{z}}$ values have also reduced due to the ICME expansion. This
example shows that our model can be used to simulate ICMEs in the IH.

---
Figure 4: The flux rope shown in Figure 3 30 hours after the initial insertion
is propagating as an ICME. The $z=0$ semi-translucent slices of (left panel)
radial velocity and (right panel) $B_{\mathrm{z}}$ are shown along with the
magnetic field lines. The gray sphere represents the inner boundary of the IH
model at 0.1 AU.

##  3 Results

###  3.1 12 July 2012 CME

In this section, we present the simulation results for the 12 July 2012 CME.
We have found the linear speed, latitude, longitude, tilt, half-angle, and
aspect ratio of this CME to be 1265 km/s, $-10^{\circ}$, $1^{\circ}$,
$52^{\circ}$, $33^{\circ}$, and 0.5, respectively, at $15\,R_{\odot}$. Since
we introduce our model CME into the SW background in a single step, we have to
introduce the CME when its apex has already reached a height of
$70\,R_{\odot}$, as we show in Figure 3. We assumed a self-similar CME
expansion from $15R_{\odot}$ to $70R_{\odot}$ to propagate the model CME from
$15\,R_{\odot}$ to $70\,R_{\odot}$. Using the drag-based model (DBM) [Vršnak &
Zic, 2007], we estimated that the CME should reach the apex height of
$70R_{\odot}$ on 13-Jul-2012 03:33 with a speed of about 1097 km/s. We used
the drag parameter of $0.1\times 10^{-7}$ and the asymptotic SW speed equal to
450 km/s here. Vršnak et al. [2014] found that this combination of drag
parameters ensures roughly the same arrival time accuracy for the DBM and WSA-
ENLIL-Cone models. We also found that the poloidal flux carried by this CME is
$1.4\times 10^{22}$ Mx. This value is found from the reconnected flux under
the PEA of the source active region using the method described by Gopalswamy
et al. [2018]. We used this value in the empirical relation between the ICME
poloidal and toroidal fluxes [Qiu et al., 2007] to calculate the toroidal flux
of this CME to be $7.6\times 10^{21}$ Mx. The helicity sign of this CME was
found to be positive [Singh et al., 2019]. Using the CME brightness in the
coronagraph images, we found that the CME mass is $1.65\times 10^{16}$ g. The
brightness is due to the Thomson scattering of photospheric light by the
plasma electrons [Billings, 1966]. By integrating over the CME area and
removing the projection effects with the help of multiple coronagraph images,
we can calculate the true mass of the CME [Colaninno & Vourlidas, 2009], which
includes the contribution from the CME sheath, cavity, and core. However, we
insert only the flux-rope part of a CME into the ambient SW, which resides in
the cavity region and has a density much lower than the sheath and core
densities [Riley et al., 2008]. Therefore, we added only $1.65\times 10^{13}$
g of uniformly distributed extra mass inside the model flux-rope to the
ambient (background) SW. As shown below, introducing this fraction of the
total mass ensures a reasonable agreement with the observed ICME density at
Earth. However, it needs to be investigated further if this fraction of total
mass is appropriate for other CMEs.

Figure 5 shows our comparison between the observations and simulations at
Earth. We use the 1 hour averaged SW data provided by NASA/GSFC through
OMNIWeb [King & Papitashvili, 2005]. Visual inspection shows that the ICME
arrives at Earth with a turbulent sheath, followed by a magnetic cloud with a
smoothly rotating magnetic field. The distribution of magnetic field in the
cloud has peaks in the negative direction for all $RTN$ components. Please
note that by peaks, we mean maximum excursion, and this definition is valid
for the rest of the paper. Comparing our simulation results to observations,
we find that there is a reasonable agreement in the density, speed, and
magnetic field values. Our simulations were able to reproduce the signs of the
peaks in all three components of the magnetic field inside the ICME cloud. The
peak values of the $RTN$ components of the magnetic field inside the observed
cloud were (-12,-9,-18) nT. The corresponding values in the simulated ICME
cloud were (-11, -14, -14) nT. This agreement in the magnetic field seems to
be better than in the previously reported simulation studies of the same CME
[Shen et al., 2014, Scolini et al., 2019, Singh et al., 2020]. The simulated
ICME arrived just 1 hour behind the observed ICME. The simulated density, both
in the sheath region and the magnetic cloud, matches well with the
observations, although we inserted only 0.1% of the observed mass in the
initial flux rope. However, we do notice that the density in the sheath region
is much smoother in our simulation as compared to the observation. We also
note that the simulated speed at the ICME front agrees with observations.
Inside the ICME, our simulation shows speeds moderately lower than in the
observations. Worth noticing are some differences between the model background
SW and OMNI data ahead of the interplanetary shock arrival on 14 July 2012,
which can be a possible source of errors in the sheath region. Firstly, the
background magnetic field polarity is different in the modeled and observed
SW, which contributes to the errors in the simulated sheath magnetic field.
Secondly, the model SW is slightly denser and faster than observed.

In our previous study, we simulated this CME using the modified spheromak
model [Singh et al., 2020], which did not reproduce the negative value of
$B_{\mathrm{r}}$. In future studies, we will compare the modified spheromak
and constant turn models in more depth to determine whether the constant turn
model performs consistently better than the modified spheromak model. This
ICME had been modeled previously by Shen et al. [2014] and Scolini et al.
[2019] with the magnetized plasma blob and spheromak models, respectively. The
magnetic field components obtained with our constant turn model seem to agree
with observations much better than in the above-mentioned studies.

The kinematic and magnetic properties of CMEs derived from observations
contain some level of uncertainty. In the next subsections, we will address
the uncertainties associated with GCS fittings and analyze their effects on
the CME simulation results at Earth. We only focus on the uncertainties due to
GCS fitting. The impact of the uncertainties in the magnetic field flux
estimates on the CME properties at Earth will be considered elsewhere.

---
Figure 5: Comparison of the 1 hr averaged OMNI SW data (blue) and simulation
results at Earth (black). The simulated ICME arrived 1 hour after the observed
ICME. The signs of the peaks in the magnetic field distribution inside the
ICME cloud agree well for all three components. The $N$-component of any
vector in the $RTN$ coordinate system is comparable to the $z$-component in
the GSE coordinate system.

###  3.2 GCS uncertainty analysis

The GCS model is a useful tool to derive the 3D properties of CMEs using
multiple viewpoint images. However, the derived properties can have large
uncertainties, due to the idealization of the CME shape to be a symmetrical,
two-conical-legs structure with a curved front that is propagating radially
outward. In reality, CMEs can be highly asymmetrical, skewed, and/or slanted
with respect to the radial direction. Moreover, SOHO and STEREO A&B roughly
belong to the same ecliptic plane. The direction of the CME out of this
ecliptic plane is roughly its latitude value, and it can be easily fitted with
the GCS model because all three viewpoints will show a similar extent of the
CME out of this plane. The direction of the CME in the ecliptic plane is
roughly its longitude value. Estimating this direction with the GCS model is
challenging compared to estimating the latitude because all three viewpoints
will show different extents of the CME in this direction. This means that GCS
fitting based on SOHO and STEREO data can result in larger ambiguities in the
longitude estimates as compared to latitude estimates.

Apparently, the estimates of GCS parameters can vary depending on the fitting
approach. Here we try to estimate the average subjective uncertainties for the
same CMEs when fitted by different GCS users. Using GCS fitting results
reported in different publications, as well as in online catalogs, we compiled
6 lists of CMEs. For example, Gopalswamy et al. [2014] fitted 77 CMEs that had
erupted between Feb 2010 and Jan 2014 and reported their latitudes,
longitudes, and speeds. Thernisien et al. [2009] reported the latitudes,
longitudes, tilts, and half-angles of 26 CMEs that erupted between Nov 2007
and Aug 2008. Shi et al. [2015] reported latitudes, longitudes, and speeds of
21 CMEs that erupted between Dec 2008 and Oct 2012. CME Kinematic Database
(KINCAT) [Millward et al., 2013] reports GCS latitudes, longitudes, tilts,
speeds, half-angles, and aspect ratios of 122 CMEs that erupted between May
2007 and Oct 2013. The full halo CME (FHCME) catalog [Shen et al., 2013] has
GCS fitting results for the latitudes, longitudes, and speeds of 39 front-side
halo CMEs between Dec 2009 and May 2012. Finally, we also used our own
unpublished catalog of GCS fittings for 49 CMEs that erupted between Feb 2010
and July 2012. All these GCS fittings apply linear regression to the height-
time data to obtain CME speed.

Using the above-mentioned CME lists, we identified the CMEs that appeared in
at least two of the lists and published the GCS parameters of those CMEs (see
https://doi.org/10.5281/zenodo.5515483), which allowed us to compare the
estimated latitude and longitude values of 56 CMEs, tilt values of 22 CMEs,
and speed values of 18 CMEs. We define the subjective uncertainty in a GCS
parameter for a CME as the range of measurements i.e. difference between
maximum and minimum values of this parameter reported in the lists that
contained the same CME. For example, if a CME appeared in three of the lists
that had the latitude estimates of $1^{\circ}$, $5^{\circ}$, and $10^{\circ}$,
the subjective uncertainty in latitude for that CME would be $9^{\circ}$. The
speed uncertainty can be better represented by the speed fractional
uncertainty because speed uncertainty itself can be higher for faster CMEs
than for slower ones. The speed fractional uncertainty can be specified as

$\textup{Speed Fractional Uncertainty}=\frac{\textup{Speed
Uncertainty}}{\textup{Average speed in separate fits}}$

The distributions of subjective uncertainties for GCS estimates of the
considered CMEs are shown in Figure 6. We further found the average of these
uncertainties. We found that the GCS estimates of latitude, longitude, and
tilt show an average uncertainty of $5.7^{\circ}$, $11.2^{\circ}$, and
$24.7^{\circ}$, respectively. We can see that the latitude uncertainties are
significantly smaller than the longitude uncertainties, which is because the
observing spacecraft are all roughly in the ecliptic plane. The tilt
uncertainties are much higher than the latitude and longitude uncertainties.
We can see that 4 of the 22 CMEs had tilt uncertainties of more than
$40^{\circ}$. We also found that, on average, CMEs can have uncertainties in
their speed of about 11.4%.

Thernisien et al. [2009] used the sensitivity analysis to estimate the average
uncertainties in GCS parameters of 26 CMEs they considered in their study.
They found the average uncertainties in latitude, longitude, and tilt to be
$1.8^{\circ}$, $4.3^{\circ}$, and $22^{\circ}$, respectively, which are
smaller than the values we found by comparing the fittings provided by
different users of the GCS model. The difference in these values is due to the
choice of CMEs used in the studies and/or due to the different methods of
finding these uncertainties. Thernisien et al. [2009] used a theoretical
approach to estimate the uncertainties whereas we used a practical approach in
which we directly compared the fitting results of different users of the GCS
model. The average uncertainties in the CME tilt in Thernisien et al. [2009]
and in our analysis are comparable.

---
Figure 6: Histograms showing the distribution of uncertainties in GCS
estimates of (top-left) latitude, (top-right) longitude, (bottom-left) tilt,
and (bottom-right) speed fraction.

###  3.3 Ensemble modeling

Pizzo et al. [2015] showed that MHD modeling of CMEs does not depend on the
changes in their initial properties chaotically. This means that ensemble
modeling of CMEs can be performed to get estimates of the CME forecast
uncertainties by creating ensemble members according to the uncertainties in
the initial CME model parameters. To study the propagation of uncertainties in
the GCS parameters to the CME forecasting uncertainties, we have performed an
ensemble modeling for the 12 July 2012 CME using our constant-turn model. We
created the ensemble members by taking into consideration the average
uncertainties in latitude, longitude, tilt, and speed. Out of the five
external lists of reported GCS parameters considered in our study, only the
KINCAT database reports half-angles and aspect ratios. Therefore, we could not
find the average uncertainties of aspect ratio and half-angle by comparing
them across the lists. Instead, we used the uncertainty values reported in
Thernisien et al. [2009]. They were $+0.07/-0.04$ and $+13^{\circ}/-7^{\circ}$
for the aspect ratios and half-angles, respectively.

We created 77 ensemble members for our CME. These members are described by
different GCS parameters used to initiate the CME model. The first member is
the one with our own estimated values of the 12 July 2012 CME parameters. An
ensemble member with the GCS properties of latitude $\theta$, longitude
$\phi$, tilt $\gamma$, speed $V_{\mathrm{CME}}$, half-angle $\alpha$, and
aspect ratio $\kappa$ can be represented by a set
$\\{\theta,\phi,\gamma,V_{\mathrm{CME}},\alpha,\kappa\\}$. The other 76
ensemble members can be represented by the following sets:

  1. 1.

$\\{\theta\pm 5.7^{\circ},\phi\pm 11.2^{\circ},\gamma\pm
24.7^{\circ},V_{\mathrm{CME}}\pm
11.4\%,\alpha_{-7^{\circ}}^{+13^{\circ}},\kappa_{-0.04}^{+0.07}\\}$, (64
members)

  2. 2.

$\\{\theta\pm 5.7^{\circ},\phi,\gamma,V_{\mathrm{CME}},\alpha,\kappa\\}$ (2
members)

  3. 3.

$\\{\theta,\phi\pm 11.2^{\circ},\gamma,V_{\mathrm{CME}},\alpha,\kappa\\}$ (2
members)

  4. 4.

$\\{\theta,\phi,\gamma\pm 24.7^{\circ},V_{\mathrm{CME}},\alpha,\kappa\\}$ (2
members)

  5. 5.

$\\{\theta,\phi,\gamma,V_{\mathrm{CME}}\pm 11.4\%,\alpha,\kappa\\}$ (2
members)

  6. 6.

$\\{\theta,\phi,\gamma,V_{\mathrm{CME}},\alpha_{-7^{\circ}}^{+13^{\circ}},\kappa\\}$
(2 members)

  7. 7.

$\\{\theta,\phi,\gamma,V_{\mathrm{CME}},\alpha,\kappa_{-0.04}^{+0.07}\\}$ (2
members)

The 64 members in set 1 represent all possible combinations of GCS parameter
extremes, where the GCS parameters we estimated for the 12 July 2012 CME
differ according to the average uncertainties in these parameters. The
simulation of these ensemble members should give us the upper and lower bounds
of simulation uncertainties due to subjective uncertainties in the GCS
fitting. The 12 members in sets 2–7 are used in this work to study the
individual contribution of uncertainties in each parameter to the simulation
results at Earth.

In the left panel of Figure 7, we compare all ensemble members with OMNI data
at Earth. The members show a 12 hr wide arrival window bounded by grey
vertical lines. The peak speed at the time of arrival ranges between 480 and
610 km/s. The peak density of the front of ICMEs ranges between 18 and 26
cm-3. In the right panel of Figure 7, we show the ensemble members by shifting
them in time to match the individual arrival times to that of the first
ensemble member. This is done for an easier qualitative comparison of the
plasma and magnetic field values in the ensemble members. While 55% of the
ensemble members correctly reproduced the peak signs for all three components
of the cloud magnetic field, a significant fraction did not correctly
reproduce the sign of the peak in the radial magnetic field component.

---
Figure 7: (Left panel) Plasma and magnetic field in the 1 hr averaged OMNI SW
data (blue) and simulation results probed at Earth (black and green). The
black line represents the first ensemble member with GCS parameters found by
us. The green lines represent the rest of the 76 ensemble members defined in
Section 3.3. The grey vertical lines show the 12 hour wide arrival time window
of all the ensemble members. (Right panel) Same as the left panel, with the
arrival time of all ensemble members matched with the first member.

To investigate the effects of uncertainty in each GCS parameter on the
simulation results, we used sets 2–7, where each of the GCS parameters was
modified according to its average uncertainties, keeping the rest of the
parameters unchanged. This exercise provides a rather good insight into the
role of uncertainties in the individual GCS parameters in producing the CME
model uncertainties at Earth.

The left panel of Figure 8 shows the effect of uncertainties in latitude on
CME forecasting at Earth. The latitude uncertainty of $\pm 5.74^{\circ}$
widens the arrival time window by 2 hours. This uncertainty, however, does not
change the signs of the peaks in the cloud magnetic field components, which
all remain negative. As mentioned in Section 3.1, the peak $RTN$ components of
the ICME magnetic cloud of our simulated CME were (-11, -14, -14) nT. By
varying the latitude by $\pm 5.74^{\circ}$, we obtain the peak values of the
$RTN$ components equal to (-9, -18, -19) and (-8, -11, -11) nT.

The effect of the longitude uncertainty on magnetic field forecasts in the
cloud is much more pronounced, as can be seen in the right panel of Figure 8.
The longitude uncertainty of $\pm 11.23^{\circ}$ gives the range of the peaks
in the $RTN$ components equal to (-7, -7, -7) and (4, -20, -16) nT. We can see
that a $-11.23^{\circ}$ change in latitude results in the wrong sign of
$B_{\mathrm{r}}$. Besides, a $\pm 11.23^{\circ}$ uncertainty in longitude
widens the arrival time window by 3 hours.

---
Figure 8: The distributions of plasma velocity and density, and magnetic field
components in the 1 hr averaged OMNI data (blue) and simulation results at
Earth (black, red, and green). The black lines represent the first ensemble
member with our derived GCS parameters. The red and green lines represent the
two ensemble members defined in set 2 (Left panel) and 3 (Right panel) in
Section 3.3.

The tilt uncertainty of $\pm 24.71^{\circ}$ can significantly alter the
magnetic field profile in the magnetic cloud of ICME at Earth, as shown in the
left panel of Figure 9. Such tilt range can change the peak values of the
$RTN$ components of the cloud magnetic field from (-11, -14, -14) nT to (-5,
-17, -20) or (-8, -10, -10) nT. The tilt uncertainty, however, does not change
the signs of the peak magnetic field components from being negative. The
arrival time window is widened by 1 hour due to the tilt uncertainty.

In the right panel of Figure 9, we show the effect of a $\pm 11.4\%$
uncertainty in speed on the simulated signatures of ICMEs at Earth. One can
see that the speed uncertainty does not change the magnetic field profile in
the ICME cloud. This uncertainty, however, widens the arrival time window by 6
hours. We also found that the peak density stays unaffected by the ICME speed
change. This result is in line with the recent observational analysis by
Temmer et al. [2021]. It should also be noted that the arrival time window
with the full set of ensemble members had a width of 12 hours. This shows that
the arrival time window can be widened by as much as 100% by the combination
of uncertainties in CME parameters, as compared to uncertainties in initial
speed only.

---
Figure 9: Plasma speed and density, and magnetic field component distributions
in the 1 hr averaged OMNI SW data (blue) and simulation results at Earth
(black, red, and green). The black line represents the first ensemble member
with our derived GCS parameters. The red and green lines represent the two
ensemble members defined in set 4 (Left panel) and 5 (Right panel) in Section
3.3.

Figure 10 shows the effect of $+13^{\circ}/-7^{\circ}$ half-width
$+0.07/-0.04$ aspect ratio uncertainties (the left and right panels,
respectively). Both of these uncertainties result in very small changes in the
magnetic field and plasma properties at Earth. The arrival time is also
largely unaffected by these uncertainties. The ensemble member with the larger
half-width has a slightly larger density peak of the front (by 2 cm-3), which
might be due to a larger mass pile-up at the front during its interplanetary
travel. This has been shown also in the observations of Temmer et al. [2021],
which show that wider ICMEs tend to have stronger mass pile-up at the front.
Since the 12 July, 2012 CME was Earth-directed, we did not expect the half-
width to play a major role in controlling the 1 AU signature of this ICME.
However, off the Sun-Earth line CMEs can be affected by the uncertainties in
half-width to a larger extent. This is because such uncertainties can create
an ICME with magnetic clouds either passing over Earth or missing it.

---
Figure 10: Plasma speed and density, and magnetic field components in the 1 hr
averaged OMNI SW data (blue) and simulation results at Earth (black, red, and
green). The black line represents the first ensemble member with our derived
GCS parameters. The red and green lines represent the two ensemble members
defined in set 6 (Left panel) and 7 (Right panel) in Section 3.3.

##  4 Conclusion

In this study, we have implemented a flux rope model with a constant-turn
magnetic field [Vandas & Romashets, 2017] to perform MHD simulations of the 12
July 2012 CME propagation through a time-dependent, data-driven ambient SW.
The CME shape used in our model is based on the FRiED approach described by
Isavnin [2016]. It can be matched well with the shape approximation used in
the GCS model, which is frequently used to derive 3D properties of CMEs from
the analysis of multi-viewpoint images made by coronagraphs. We constrain the
poloidal magnetic flux in this model with the measured reconnected flux in the
PEAs of the source active region. The toroidal magnetic flux is estimated from
the poloidal flux using the empirical relation given by Qiu et al. [2007]. We
have shown that by constraining the flux rope model with GCS and PEA
observations of the 12 July 2012 CME, the ICME features at Earth are
reproduced with reasonable accuracy. Our simulated CME agrees with the
observations of speed, density, and magnetic field components signs in the
ICME magnetic cloud. The arrival time of our simulated ICME was off by just 1
hour.

We compared the GCS parameters calculated independently by multiple users to
derive the average subjective uncertainties in the CME latitude, longitude,
tilt, and speed, which turned out to be $5.74^{\circ}$, $11.23^{\circ}$,
$24.71^{\circ}$, and 11.4% respectively. Using these uncertainties, as well as
those in the CME half-angle and aspect ratio ($+13^{\circ}/-7^{\circ}$ and
$+0.07/-0.04$, respectively) from Thernisien et al. [2009], we have created 77
ensemble members for the 12 July 2012 CME. The results of our ensemble
modeling are summarized as follows:

  1. 1.

The arrival time window for the full ensemble is 12 hours wide. This value is
similar to the typical arrival time errors obtained with different models.
This 12-hour width of the full ensemble is due to the combined non-linear
impact of uncertainties in all GCS model parameters.

  2. 2.

55% of the ensemble members correctly reproduced the signs of all $RTN$
components of the peak magnetic field in the ICME magnetic cloud, while 28% of
the ensemble members did not produce the correct sign of $B_{\mathrm{r}}$.

  3. 3.

Subjective uncertainties of $\pm 11.23^{\circ}$ in longitude and $\pm
24.71^{\circ}$ in the CME tilt are found to be the main sources of uncertainty
in the simulated magnetic field. Subjective uncertainty of $\pm 5.74^{\circ}$
in latitude had a lesser impact on the simulated magnetic field values. This
highlights the importance of out-of-ecliptic coronagraph observations of CMEs,
which would be very important for better constraining the CME longitude.

  4. 4.

Uncertainties in the half-width and the aspect ratio are hardly affecting the
ICME features at 1 AU. This is not surprising, since our CME was Earth-
directed. The half-angle of a model CME may play an important role for the
CMEs that are not aligned with the Sun-Earth line.

It takes about 2 hours of computation for MS-FLUKSS to simulate about 6 days
of physical time on a $150\times 128\times 256$ grid using 128 CPUs. Since we
can run all ensemble members in parallel, the ensemble modeling we have
performed here is feasible for real-time space weather forecasting. It is of
importance to note that we have not considered the uncertainties in the input
magnetic fluxes. Incorporating those uncertainties into the analysis may
significantly increase our ensemble size. However, since the half-width and
aspect ratio of the CME were shown to have minimal effect on the simulation
results at Earth, we could still create ensembles of reasonable size by
excluding those members.

ACKNOWLEDGEMENT

The authors acknowledge support from NASA/NSF SWQU grant 2028154. TKK
acknowledges support from AFOSR grant FA9550-19-1-0027. NP was also supported,
in part, by NSF-BSF grant 2010450 and NASA grants 80NSSC19K0075 and
80NSSC21K0004. Supercomputer allocations were provided on SGI Pleiades by NASA
High-End Computing Program award SMD-21-44038581 and also on TACC Stampede2
and SDSC Expanse by NSF XSEDE project MCA07S033. We acknowledge the NASA/GSFC
Space Physics Data Facility’s OMNIWeb for the SW and IMF data used in this
study. We also used SOHO and STEREO coronagraph data from https://lasco-
www.nrl.navy.mil and stereo-ssc.nascom.nasa.gov respectively. SDO EUV and
magnetogram data have been taken from
http://jsoc.stanford.edu/ajax/exportdata.html. This work utilizes data
produced collaboratively between the Air Force Research Laboratory (AFRL) and
the National Solar Observatory.

## References

  * Amari et al. [2014] Amari, T., Canou, A., & Aly, J.-J. 2014, Nature, 514, 465, doi: 10.1038/nature13815
  * Amerstorfer et al. [2018] Amerstorfer, T., Möstl, C., Hess, P., et al. 2018, Space Weather, 16, 784, doi: https://doi.org/10.1029/2017SW001786
  * Amerstorfer et al. [2021] Amerstorfer, T., Hinterreiter, J., Reiss, M. A., et al. 2021, Space Weather, 19, e2020SW002553, doi: https://doi.org/10.1029/2020SW002553
  * Arge et al. [2005] Arge, C. N., de Toma, G., & Luhmann, J. G. 2005, in Astronomical Society of the Pacific Conference Series, Vol. 346, Large-scale Structures and their Role in Solar Activity, ed. K. Sankarasubramanian, M. Penn, & A. Pevtsov, 371
  * Arge et al. [2013] Arge, C. N., Henney, C. J., Hernandez, I. G., et al. 2013, AIP Conference Proceedings, 1539, 11, doi: 10.1063/1.4810977
  * Arge et al. [2010] Arge, C. N., Henney, C. J., Koller, J., et al. 2010, AIP Conference Proceedings, 1216, 343, doi: 10.1063/1.3395870
  * Arge et al. [2011] Arge, C. N., Henney, C. J., Koller, J., et al. 2011, in Astronomical Society of the Pacific Conference Series, Vol. 444, 5th International Conference of Numerical Modeling of Space Plasma Flows (ASTRONUM 2010), ed. N. V. Pogorelov, E. Audit, & G. P. Zank, 99
  * Arge et al. [2003] Arge, C. N., Odstrcil, D., Pizzo, V. J., & Mayer, L. R. 2003, AIP Conference Proceedings, 679, 190, doi: 10.1063/1.1618574
  * Billings [1966] Billings, D. E. 1966, A guide to the solar corona
  * Brueckner et al. [1995] Brueckner, G. E., Howard, R. A., Koomen, M. J., et al. 1995, Solar Physics, 162, 357, doi: 10.1007/BF00733434
  * Brueckner et al. [1998] Brueckner, G. E., Delaboudiniere, J. P., Howard, R. A., et al. 1998, grl, 25, 3019, doi: 10.1029/98GL00704
  * Burlaga et al. [1981] Burlaga, L., Sittler, E., Mariani, F., & Schwenn, R. 1981, Journal of Geophysical Research: Space Physics, 86, 6673, doi: 10.1029/JA086iA08p06673
  * Carmichael [1964] Carmichael, H. 1964, A Process for Flares, Vol. 50, 451
  * Chen [2011] Chen, P. F. 2011, Living Reviews in Solar Physics, 8, 1, doi: 10.12942/lrsp-2011-1
  * Colaninno & Vourlidas [2009] Colaninno, R. C., & Vourlidas, A. 2009, The Astrophysical Journal, 698, 852, doi: 10.1088/0004-637X/698/1/852
  * Dryer & Smart [1984] Dryer, M., & Smart, D. 1984, Advances in Space Research, 4, 291 , doi: https://doi.org/10.1016/0273-1177(84)90200-X
  * Elliott et al. [2016] Elliott, H. A., McComas, D. J., & DeForest, C. E. 2016, The Astrophysical Journal, 832, 66, doi: 10.3847/0004-637x/832/1/66
  * Fry et al. [2001] Fry, C. D., Sun, W., Deehr, C. S., et al. 2001, jgr, 106, 20985, doi: 10.1029/2000JA000220
  * Gibson & Low [1998] Gibson, S. E., & Low, B. C. 1998, apj, 493, 460, doi: 10.1086/305107
  * Gold & Hoyle [1960] Gold, T., & Hoyle, F. 1960, Monthly Notices of the Royal Astronomical Society, 120, 89, doi: 10.1093/mnras/120.2.89
  * Gopalswamy et al. [2018] Gopalswamy, N., Akiyama, S., Yashiro, S., & Xie, H. 2018, in IAU Symposium, Vol. 335, Space Weather of the Heliosphere: Processes and Forecasts, ed. C. Foullon & O. E. Malandraki, 258–262, doi: 10.1017/S1743921317011048
  * Gopalswamy et al. [2018] Gopalswamy, N., Akiyama, S., Yashiro, S., & Xie, H. 2018, Journal of Atmospheric and Solar-Terrestrial Physics, 180, 35 , doi: https://doi.org/10.1016/j.jastp.2017.06.004
  * Gopalswamy et al. [2005] Gopalswamy, N., Lara, A., Manoharan, P. K., & Howard, R. A. 2005, Advances in Space Research, 36, 2289, doi: 10.1016/j.asr.2004.07.014
  * Gopalswamy et al. [2001] Gopalswamy, N., Lara, A., Yashiro, S., Kaiser, M. L., & Howard, R. A. 2001, jgr, 106, 29207, doi: 10.1029/2001JA000177
  * Gopalswamy et al. [2013] Gopalswamy, N., Mäkelä, P., Xie, H., & Yashiro, S. 2013, Space Weather, 11, 661, doi: 10.1002/2013SW000945
  * Gopalswamy et al. [2014] Gopalswamy, N., Xie, H., Akiyama, S., Mäkelä, P. A., & Yashiro, S. 2014, Earth, Planets and Space, 66, 104, doi: 10.1186/1880-5981-66-104
  * Gosling [1990] Gosling, J. T. 1990, Washington DC American Geophysical Union Geophysical Monograph Series, 58, 343, doi: 10.1029/GM058p0343
  * Hickmann et al. [2015] Hickmann, K. S., Godinez, H. C., Henney, C. J., & Arge, C. N. 2015, Solar Physics, 290, 1105
  * Hirayama [1974] Hirayama, T. 1974, solphys, 34, 323, doi: 10.1007/BF00153671
  * Hoeksema et al. [2014] Hoeksema, J. T., Liu, Y., Hayashi, K., et al. 2014, Solar Physics, 289, 3483, doi: 10.1007/s11207-014-0516-8
  * Howard et al. [2008] Howard, R. A., Moses, J. D., Vourlidas, A., et al. 2008, Space Science Reviews, 136, 67, doi: 10.1007/s11214-008-9341-4
  * Hu et al. [2015] Hu, Q., Qiu, J., & Krucker, S. 2015, Journal of Geophysical Research (Space Physics), 120, 5266, doi: 10.1002/2015JA021133
  * Isavnin [2016] Isavnin, A. 2016, apj, 833, 267, doi: 10.3847/1538-4357/833/2/267
  * Jin et al. [2017] Jin, M., Manchester, W. B., van der Holst, B., et al. 2017, The Astrophysical Journal, 834, 172, doi: 10.3847/1538-4357/834/2/172
  * Jin et al. [2016] Jin, M., Schrijver, C. J., Cheung, M. C. M., et al. 2016, The Astrophysical Journal, 820, 16, doi: 10.3847/0004-637X/820/1/16
  * Kaiser et al. [2008] Kaiser, M. L., Kucera, T. A., Davila, J. M., et al. 2008, Space Science Reviews, 136, 5, doi: 10.1007/s11214-007-9277-0
  * Kim et al. [2014] Kim, T. K., Pogorelov, N. V., Borovikov, S. N., et al. 2014, Journal of Geophysical Research (Space Physics), 119, 7981, doi: 10.1002/2013JA019755
  * Kim et al. [2020] Kim, T. K., Pogorelov, N. V., Arge, C. N., et al. 2020, apjs, 246, 40, doi: 10.3847/1538-4365/ab58c9
  * King & Papitashvili [2005] King, J. H., & Papitashvili, N. E. 2005, Journal of Geophysical Research: Space Physics, 110, doi: 10.1029/2004JA010649
  * Kopp & Pneuman [1976] Kopp, R. A., & Pneuman, G. W. 1976, solphys, 50, 85, doi: 10.1007/BF00206193
  * Linker et al. [2016] Linker, J. A., Caplan, R. M., Downs, C., et al. 2016, Journal of Physics: Conference Series, 719, 012012, doi: 10.1088/1742-6596/719/1/012012
  * Linker et al. [2017] —. 2017, The Astrophysical Journal, 848, 70, doi: 10.3847/1538-4357/aa8a70
  * MacNeice et al. [2011] MacNeice, P., Elliott, B., & Acebal, A. 2011, Space Weather, 9, doi: 10.1029/2011SW000665
  * Manchester et al. [2004] Manchester, W. B., Gombosi, T. I., Roussev, I., et al. 2004, Journal of Geophysical Research (Space Physics), 109, A02107, doi: 10.1029/2003JA010150
  * Marubashi [2000] Marubashi, K. 2000, Advances in Space Research, 26, 55 , doi: https://doi.org/10.1016/S0273-1177(99)01026-1
  * Mays et al. [2015] Mays, M. L., Taktakishvili, A., Pulkkinen, A., et al. 2015, solphys, 290, 1775, doi: 10.1007/s11207-015-0692-1
  * McGregor et al. [2011] McGregor, S. L., Hughes, W. J., Arge, C. N., Owens, M. J., & Odstrcil, D. 2011, Journal of Geophysical Research (Space Physics), 116, A03101, doi: 10.1029/2010JA015881
  * Millward et al. [2013] Millward, G., Biesecker, D., Pizzo, V., & de Koning, C. A. 2013, Space Weather, 11, 57, doi: https://doi.org/10.1002/swe.20024
  * Moon et al. [2002] Moon, Y. J., Dryer, M., Smith, Z., Park, Y. D., & Cho, K. S. 2002, grl, 29, 1390, doi: 10.1029/2002GL014865
  * Murray [2018] Murray, S. A. 2018, Space Weather, 16, 777, doi: https://doi.org/10.1029/2018SW001861
  * Odstrcil et al. [2004] Odstrcil, D., Riley, P., & Zhao, X. P. 2004, Journal of Geophysical Research (Space Physics), 109, A02116, doi: 10.1029/2003JA010135
  * Pesnell et al. [2012] Pesnell, W. D., Thompson, B. J., & Chamberlin, P. C. 2012, Solar Physics, 275, 3, doi: 10.1007/s11207-011-9841-3
  * Pizzo et al. [2015] Pizzo, V. J., de Koning, C., Cash, M., et al. 2015, Space Weather, 13, 676, doi: https://doi.org/10.1002/2015SW001221
  * Pogorelov et al. [2014] Pogorelov, N., Borovikov, S., Heerikhuisen, J., et al. 2014, in Proceedings of the 2014 Annual Conference on Extreme Science and Engineering Discovery Environment, XSEDE ’14 (New York, NY, USA: ACM), 22:1–22:8, doi: 10.1145/2616498.2616499
  * Qiu et al. [2007] Qiu, J., Hu, Q., Howard, T. A., & Yurchyshyn, V. B. 2007, The Astrophysical Journal, 659, 758, doi: 10.1086/512060
  * Riley et al. [2008] Riley, P., Lionello, R., Mikić, Z., & Linker, J. 2008, The Astrophysical Journal, 672, 1221, doi: 10.1086/523893
  * Riley et al. [2018] Riley, P., Mays, M. L., Andries, J., et al. 2018, Space Weather, 16, 1245, doi: 10.1029/2018SW001962
  * Sarkar et al. [2020] Sarkar, R., Gopalswamy, N., & Srivastava, N. 2020, apj, 888, 121, doi: 10.3847/1538-4357/ab5fd7
  * Schatten [1971] Schatten, K. H. 1971, Cosmic Electrodynamics, 2, 232
  * Schou et al. [2012] Schou, J., Scherrer, P. H., Bush, R. I., et al. 2012, Solar Physics, 275, 229, doi: 10.1007/s11207-011-9842-2
  * Scolini et al. [2019] Scolini, C., Rodriguez, L., Mierla, M., Pomoell, J., & Poedts, S. 2019, A&A, 626, A122, doi: 10.1051/0004-6361/201935053
  * Shen et al. [2013] Shen, C., Wang, Y., Pan, Z., et al. 2013, Journal of Geophysical Research: Space Physics, 118, 6858, doi: https://doi.org/10.1002/2013JA018872
  * Shen et al. [2011] Shen, F., Feng, X. S., Wu, S. T., Xiang, C. Q., & Song, W. B. 2011, Journal of Geophysical Research: Space Physics, 116, doi: https://doi.org/10.1029/2010JA015809
  * Shen et al. [2014] Shen, F., Shen, C., Zhang, J., et al. 2014, Journal of Geophysical Research: Space Physics, 119, 7128, doi: 10.1002/2014JA020365
  * Shi et al. [2015] Shi, T., Wang, Y., Wan, L., et al. 2015, The Astrophysical Journal, 806, 271, doi: 10.1088/0004-637x/806/2/271
  * Shiota & Kataoka [2016] Shiota, D., & Kataoka, R. 2016, Space Weather, 14, 56, doi: 10.1002/2015SW001308
  * Singh et al. [2020] Singh, T., 1Yalim, M. S., Pogorelov, N. V., & Gopalswamy, N. 2020, The Astrophysical Journal, Submitted
  * Singh et al. [2020] Singh, T., Kim, T. K., Pogorelov, N. V., & Arge, C. N. 2020, Space Weather, 18, e02405, doi: 10.1029/2019SW002405
  * Singh et al. [2018] Singh, T., Yalim, M. S., & Pogorelov, N. V. 2018, The Astrophysical Journal, 864, 18, doi: 10.3847/1538-4357/aad3b4
  * Singh et al. [2019] Singh, T., Yalim, M. S., Pogorelov, N. V., & Gopalswamy, N. 2019, The Astrophysical Journal, 875, L17, doi: 10.3847/2041-8213/ab14e9
  * Song et al. [2020] Song, H. Q., Zhang, J., Cheng, X., et al. 2020, The Astrophysical Journal, 901, L21, doi: 10.3847/2041-8213/abb6ec
  * Sturrock [1968] Sturrock, P. A. 1968, in IAU Symposium, Vol. 35, Structure and Development of Solar Active Regions, ed. K. O. Kiepenheuer, 471
  * Temmer et al. [2021] Temmer, M., Holzknecht, L., Dumbović, M., et al. 2021, Journal of Geophysical Research: Space Physics, 126, e2020JA028380, doi: https://doi.org/10.1029/2020JA028380
  * Thernisien [2011] Thernisien, A. 2011, apjs, 194, 33, doi: 10.1088/0067-0049/194/2/33
  * Thernisien et al. [2009] Thernisien, A., Vourlidas, A., & Howard, R. A. 2009, solphys, 256, 111, doi: 10.1007/s11207-009-9346-5
  * Thernisien et al. [2006] Thernisien, A. F. R., Howard, R. A., & Vourlidas, A. 2006, The Astrophysical Journal, 652, 763, doi: 10.1086/508254
  * Titov & Démoulin [1999] Titov, V. S., & Démoulin, P. 1999, aap, 351, 707
  * Vandas et al. [1996] Vandas, M., Fischer, S., Dryer, M., Smith, Z., & Detman, T. 1996, jgr, 101, 2505, doi: 10.1029/95JA02446
  * Vandas & Romashets [2017] Vandas, M., & Romashets, E. 2017, aap, 608, A118, doi: 10.1051/0004-6361/201731412
  * Verbeke et al. [2019] Verbeke, C., Pomoell, J., & Poedts, S. 2019, aap, 627, A111, doi: 10.1051/0004-6361/201834702
  * Vourlidas et al. [2010] Vourlidas, A., Howard, R. A., Esfandiari, E., et al. 2010, The Astrophysical Journal, 722, 1522, doi: 10.1088/0004-637x/722/2/1522
  * Vršnak & Zic [2007] Vršnak, B., & Zic, T. 2007, A&A, 472, 937, doi: 10.1051/0004-6361:20077499
  * Vršnak et al. [2014] Vršnak, B., Temmer, M., Žic, T., et al. 2014, The Astrophysical Journal Supplement Series, 213, 21, doi: 10.1088/0067-0049/213/2/21
  * Vršnak & Gopalswamy [2002] Vršnak, B., & Gopalswamy, N. 2002, Journal of Geophysical Research (Space Physics), 107, 1019, doi: 10.1029/2001JA000120
  * Wallace et al. [2019] Wallace, S., Arge, C. N., Pattichis, M., Hock-Mysliwiec, R. A., & Henney, C. J. 2019, Solar Physics, 294, 1
  * Wang et al. [2002] Wang, Y. M., Ye, P. Z., Wang, S., Zhou, G. P., & Wang, J. X. 2002, Journal of Geophysical Research: Space Physics, 107, SSH 2, doi: 10.1029/2002JA009244