A Series of Advances in Analytic Interplanetary CME Modeling

Coronal mass ejections (CMEs) and high speed streams (HSSs) are large‐scale transient structures that routinely propagate away from the Sun. Individually, they can cause space weather effects at the Earth, or elsewhere in space, but many of the largest events occur when these structures interact during their interplanetary propagation. We present the initial coupling of Open Solar Physics Rapid Ensemble Information (OSPREI), a model for CME evolution, with Mostly Empirical Operational Wind with a High Speed Stream, a time‐dependent HSS model that can serve as a background for the OSPREI CME. We present several improvements made to OSPREI in order to take advantage of the new time‐dependent, higher‐dimension background. This includes an update in the drag calculation and the ability to determine the rotation of a yaw‐like angle. We present several theoretical case studies, describing the difference in the CME behavior between a HSS background and a quiescent one. This behavior includes interplanetary CME propagation, expansion, deformation, and rotation, as well as the formation of a CME‐driven sheath. We also determine how the CME behavior changes with the HSS size and initial front distance. Generally, for a fast CME, we see that the drag is greatly reduced within the HSS, leading to faster CMEs and shorter travel times. The drag reappears stronger if the CME reaches the stream interaction region or upstream solar wind, leading to a stronger shock with more compression until the CME sufficiently decelerates. We model a CME–HSS interaction event observed by Parker Solar Probe in January 2022. The model improvements create a better match to the observed in situ profiles.

• We combine two existing models to simulate the interplanetary interaction of a Coronal mass ejection (CME) with a time-dependent high speed stream (HSS) in the solar wind background • A HSS reduces the drag on a CME and causes rotation if only one flank is embedded in it; if the CME reaches the stream interaction region drag suddenly increases • The combined models can reproduce the in situ profiles for the interaction of a CME and HSS observed by Parker Solar Probe in January 2022 Supporting Information: Supporting Information may be found in the online version of this article.

OSPREI
OSPREI combines three different models to form a full Sun-to-Earth (or any other relevant endpoint) simulation.The first model is Forecasting a CME's Altered Trajectory (ForeCAT; Kay et al., 2013Kay et al., , 2015)), which simulates the nonradial motion of a CME in the corona based on the solar magnetic forces.Within ForeCAT, and the rest of OSPREI, a CME is represented by a rigid torus.Both the toroidal axis and the cross section perpendicular to the axis can take on elliptical shapes.The torus can expand in size and change in both aspect ratios (axial and cross-sectional), but it is not arbitrarily deformable.
In ForeCAT, the radial motion and expansion are prescribed by the user, which can be as simple as setting the maximum speed and angular width of the CME in the corona.ForeCAT simulates the deflection of CME from the magnetic pressure gradients and magnetic tension from the background solar magnetic field, causing a change in the latitude and longitude of the CME.Any net torque from these magnetic forces causes a rotation about the radial vector extending through the nose of the CME torus.In most cases, these magnetic forces become negligible by 5 R ⊙ (Kay & Opher, 2015).Typically, we run ForeCAT until 21.5 R ⊙ (0.1 au) to ensure that we capture any coronal deflection or rotation and because it is the distance most commonly used as the inner boundary of the heliospheric domain in interplanetary simulations (e.g., Odstrcil, 2003;Verbeke et al., 2019).
At 0.1 au, OSPREI transitions into the ANother Type of Ensemble Arrival Time Results (ANTEATR; Kay & Gopalswamy, 2018;Kay, Mays, & Collado-Vega, 2022;Kay, Mays, & Verbeke, 2020;Kay & Nieves-Chinchilla, 2021a, 2021b) component.ANTEATR started out as a simple, 1D drag model.In ANTEATR, the drag force, F D is calculated in the same manner as hydrodynamic drag where C D is the drag coefficient, A CME is the cross-sectional area of the CME in the direction of motion, ρ SW is the SW density, and v CME and v SW are the CME and SW speeds, respectively.
ANTEATR originally represented the background SW with a constant velocity and the density inversely proportional to radial distance squared, with the parameters scaled using their 1 au values.A single drag force was calculated using the CME speed and the expected SW speed at the distance of the toroidal axis at the CME nose.At this time, the CME's angular size and shape were assumed to remain constant.While the drag force and evolution of the CME velocity were 1D, ANTEATR still used the full 3D CME shape to determine the time of first impact at a desired location.Kay, Mays, and Verbeke (2020) showed that using an ensemble of this first, most-simplified, version of ANTEATR was still able to reproduce arrival time results within about 6 hr.
A series of improvements have been made to ANTEATR in recent years.Kay and Nieves-Chinchilla (2021a) introduced the Physics-driven Approach to Realistic Axial Deformation and Expansion (PARADE) version.With PARADE, magnetic, thermal, and drag forces are calculated near the CME nose and flanks and used to simulate the expansion and deformation (change in aspect ratio) of a CME, rather than assuming the size and shape remain fixed.Calculating these forces requires a flux rope model and the internal temperature of the CME.As the CME has an elliptical cross section, PARADE incorporated the elliptic-cylindrical (EC) flux rope model of Nieves-Chinchilla et al. (2018), and assumed a constant internal temperature.PARADE also requires the magnetic field and temperature of the background SW, which were set to follow a Parker spiral magnetic field and a power-law dependence for the temperature.PARADE reproduces the observed tendency of CMEs to "pancake" or flatten in the radial direction relative to the perpendicular direction (e.g., Riley & Crooker, 2004).PARADE generates CMEs with average internal density, magnetic field strength, and temperature generally similar to those observed near 1 au.Kay, Nieves-Chinchilla, et al. (2022) added the Pile Up Procedure (PUP) to ANTEATR, which simulates the development of a CME-driven sheath during interplanetary propagation.Using the Rankine-Hugoniot jump conditions, PUP determines the compression at the shock (or discontinuity).This is then used to determine the amount of SW swept up into the CME-driven sheath at each time step, as well as the properties within the sheath.The current version of OSPREI includes both the PARADE and PUP updates to ANTEATR, but also allows the user to turn them off if desired.
The final component of OSPREI is the ForeCAT In situ Data Observer (FIDO; Kay et al., 2017), which generates in situ profiles given a synthetic satellite location.This initially was limited to magnetic field profiles simply using a flux rope model and the relative location of a satellite as the CME passes over it.As ANTEATR has evolved to simulate more properties, FIDO has also expanded to include the CME velocity, density, and temperature, as well as the sheath profile.Initially, ANTEATR and FIDO were run sequentially with ANTEATR ceasing at the time of first contact and FIDO assuming that the CME continued propagation with the same bulk and expansion velocities.This, however, excluded any evolution that continued beyond the time of first contact.While mostly minimal in previous versions, we expect there can be more significant continued evolution when the CME is actively interacting with a HSS.As such, ANTEATR and FIDO now run in parallel within OSPREI so that the CME continues evolving as synthetic in situ measurements are being generated.For further information on OSPREI see Kay, Mays, and Collado-Vega (2022).We also include in Figure S1 a revised version of the schematic within Kay, Mays, and Collado-Vega (2022) that shows the updated version of OSPREI used in this work.

MEOW-HiSS
MEOW-HiSS is essentially a set of coefficients, based on multilayered regressions, that can be used to scale a predetermined set of empirical functions to mimic radial profiles of a HSS from an MHD simulation.It is based upon a set of European Heliospheric Forecasting Information Asset (EUHFORIA; Pomoell & Poedts, 2018) simulations for an idealized HSS emanating from a circular coronal hole (CH) at the inner boundary of the heliospheric domain (set at 21.5 R ⊙ or 0.1 au).For each CH area, radial profiles were extracted at multiple simulation times, which correspond to different HSS front distances along a constant radial direction.The profiles for individual plasma parameters can be broken down into small regions that are well-described by simple mathematical functions (e.g., a straight line or an exponential).
MEOW-HISS uses three regions to represent the SIR, one for the HSS plateau, and another for the tail, in addition to the upstream and downstream regions (ahead of and behind the HSS, respectively).Kay et al. (2023) first found the location of the region boundaries within each profile.For each CH, a regression was performed to get each boundary location as a function of time.In most cases a first order polynomial is sufficient, suggesting that each boundary moves with nearly constant velocity in the MHD simulation.These time-dependent polynomial coefficients are then fit as a function of CH area to get a new set of regression coefficients.In practice, MEOW-HiSS takes an input CH area and the HSS front distance at the start of the simulation.It uses the CH area to determine the appropriate time-dependent coefficients for this HSS.The front distance then relates MEOW-HiSS simulation time to the corresponding time in the original MHD simulation.The boundary locations can then be determined at any time relative to the start of the MEOW-HiSS simulation.
The same approach is used for the actual HSS properties within MEOW-HiSS.For a single plasma parameter (radial speed, density, radial magnetic field, longitudinal magnetic field, or temperature), the same type of mathematical function is used to create a segment within a specific region across all times and CH areas.Each segment is constrained by one to three critical values, depending on the specific function, which are all determined directly from the MHD simulation.For example, for all CH areas and simulation times, a linear profile is used to represent the radial velocity in the first two regions of the SIR.This linear segment is constrained by critical values at two points.Kay et al. (2023) found the critical values for all combinations of CH areas and HSS front distances (front distance being a more convenient measure than MHD time).A first set of regressions relates the value at each point to the front distance, then a second regression relates these distance coefficients to the CH area.In operation, MEOW-HiSS runs in reverse of the development with the CH area determining the distance coefficients, then the distance coefficients determining the critical values, which then establish the segment within a region.Kay et al. (2023) contains a more detailed description of the development and use of MEOW-HiSS.
Figure 1 demonstrates MEOW-HiSS results in the equatorial plane for a 8 × 10 10 km 2 area CH when the HSS front is at 1 au.Panel (a) shows the regions in the empirical HSS in alternating colors of blue and orange.The front three regions correspond to the SIR, followed by an orange plateau region, and a blue tail.The black region corresponds to quiescent SW.We also label regions of "higher" and "lower" longitude.We use a heliocentric coordinate system (typically Stonyhurst) where the longitude increases toward solar west, which we will refer to as higher longitude throughout this work.At the higher longitude edge of the MEOW-HiSS HSS we see the frontmost region of the SIR extend to far radial distances but sharply drop off in longitude.This is where the empirical relations begin to break down as the HSS disappears from the domain.The effect is only this drastic in the contours of distinct regions, the relations still create a gradual decay in plasma properties.Figure 1b shows the modeled velocity for the same HSS configuration, which has a much smoother boundary at higher longitudes.
We note that this is the first iteration of MEOW-HiSS, which is based exclusively on MHD simulations of a HSS generated from a circular CH at the equator.In real cases, we expect there to be latitudinal variations and asymmetry in the CH shape.These aspects will be incorporated into future versions of MEOW-HiSS and their effects on CME propagation analyzed at that time.

Improvements in OSPREI
All of the regressions were performed on HSS MHD profiles normalized by a quiescent MHD profile.Accordingly, MEOW-HiSS returns the relative change in a parameter, which then needs to be scaled back to physical units using a quiescent profile.In OSPREI, we use the same quiescent profile that was previously used as the ANTEATR background-constant speed, density that falls with distance-squared, Parker spiral magnetic field, and power-law temperature, all of which are normalized by the values at 1 au (or the appropriate distance of interest for cases in which the satellite is not near Earth).OSPREI has a very modular design so coupling with MEOW-HiSS is essentially just an extra call to the MEOW-HiSS module after the call to the quiescent SW module.No other modifications to the code are required, however, we have made several improvements, specifically to ANTEATR, in light of having a significantly improved SW background.

Drag Calculation
The first improvement is in the manner in which we calculate the drag force, which we illustrate in the top panels of Figure 2. The gray torus represents an edge-on view of the simulated CME, and different colored arrows represent different forces.Previously, ANTEATR calculated the drag at a single point in the center of the cross section using the bulk CME speed and the expected SW velocity at this point.This is a simplification as an expanding (or contracting) CME will have a gradient in the local velocity as one passes through the CME.However, particularly in the case of a uniform background, the single-point approximation is not likely that different from the average drag in a precise calculation.
MEOW-HiSS introduces the possibility of having very different environments ahead of and behind the CME.
For example, if a HSS begins overtaking a slow CME, then we expect that the drag on the back side of the CME would tend to accelerate it.The CME front, however, may still be interacting with slow, quiescent wind, so there the drag force would want to decelerate the CME.Calculating the expected drag in the center might produce something near the average of these two forces, but it is an oversimplification that neglects the possibility that competing forces at the front and back could create a compression or squashing of the CME.We note that this is an idealized conceptualization for a simplified model, in a real event the regions directly upstream and downstream will not be pristine HSS, but rather interaction regions influenced by both the HSS and the CME.
We have converted ANTEATR to use a two-point drag calculation where we calculate the forces at the front and back of the cross section at the CME nose.The illustration in Figure 2 shows a less extreme example than the previous hypothetical example.Here, we have a slightly stronger drag force at the front (dark blue) than at the back (light blue), but both correspond to decelerating the CME.This is representative of what one might see in an expanding CME propagating faster than a uniform background SW.The pictorial vector math illustrates how the front and back drag forces are averaged to determine the bulk drag (purple arrow) but the difference determines an inward squashing force (maroon arrows).
While the example is only shown at the nose, we note that we actually calculate the drag force at both the nose and the flanks as the PARADE version of ANTEATR uses the difference in the forces to determine the deformation of the toroidal axis.The new version does the two-point calculation at both the nose and flanks.The difference in the bulk drag determines the axial deformation.We use the average squashing from the nose and flanks as we currently can only simulate a cross section that is uniform along the full toroidal axis.

Yaw Rotation
The second improvement is the addition of interplanetary rotation about an axis running through the CME nose and perpendicular to the toroidal axis.The rotation is determined from the differential drag forces at the flanks of the CME.In the case of a uniform background SW and a symmetric CME, there is no net torque, so we would not have expected any rotation with the previous version of OSPREI.MEOW-HiSS, however, breaks the symmetry.
The bottom half of Figure 2 illustrates this rotation.On the left, the black dot at the CME nose shows the axis of rotation (perpendicular to the page) and the dashed lines represent the lever arms for determining the net torque.
The blue arrows represent different drag forces, with the bottom flank having a stronger decelerating force in this example.The right shows the corresponding expected rotation.We refer to this as a rotation in the "yaw," though we note that a true yaw rotation is through the center of mass whereas we determine a rotation about the nose, primarily for the computational simplicity of keeping the nose pointed in the radial direction.
Simulating the yaw rotation requires calculating the net torque on the CME and using the moment of inertia to convert this into an angular acceleration.We first determine the radial drag force at each flank.For the CME speed, we determine the full velocity vector (bulk motion and expansion) at the outer edge of each flank (see the location of the blue arrows in the bottom left of Figure 2) and take the component in the radial direction.MEOW-HiSS only provides a 2D SW, so we project the outer flank locations onto the equatorial plane and use the corresponding position to determine the local SW properties.This ignores any latitudinal variation in the HSS, which certainly exists.This oversimplification will be the worst for a highly inclined CME.A highly inclined CME, however, will have minimal longitudinal separation so our approach will underestimate any variation between the flanks and therefore underestimate any yaw rotation.Future work will expand MEOW-HiSS to account for latitudinal variations and at that point we will refine this drag force calculation.
We then determine the lever arms, L, from the flank edges to the nose (dashed lines in bottom left panel of Figure 2) and use these, in combination with the forces, to determine the torque, τ on the CME.
The moment of inertia, I CME can then be used to determine the angular acceleration, α.To keep the problem tractable, we derived I CME for an elliptic torus with a circular cross section as where ρ CME is the density of the CME, R Ax is the half-width of the toroidal axis in the perpendicular direction, R CS is the radius of the cross section, δ Ax is the aspect ratio of the axis (radial divided by perpendicular), and K and E are complete elliptic integrals of the first and second kinds.We approximate R CS as the cross-sectional radius in the perpendicular direction for CMEs with an elliptical cross section.
We determine α at each time step, which determines a change in the angular velocity, ω, of the CME.We conserve angular momentum, I CME ω, and adjust ω at each time step based on any changes in I CME before including any additional acceleration.If α becomes negligible then the yaw rotation will slow down as the CME grows in size, which increases I CME .
This form of solid body rotation is something that can be easily calculated in an analytic code such as OSPREI, but is clearly an oversimplification of the physical processes at play.Similar approximations are often made in simple models of CME behavior, such as the frequent use of one-dimensional drag in arrival time models.CMEs are not rigid, coherent structures on global scales (Owens et al., 2017), so their overall behavior is the aggregate of a series of localized distortions.Still, we believe that the solid body rotation is a useful mathematical construct that moves OSPREI closer toward reproducing the real world behavior of CMEs.

Theoretical Case Studies
We first present several purely theoretical cases (without specific observational counterparts) and closely examine the interaction between the CME and HSS in detail before expanding to parameter space analysis and reproducing actual observations.All CMEs are initiated at 0° latitude, 0° longitude, and with their toroidal axis lying in the equatorial plane (0° inclination).Table 1 lists the remaining input parameters for these cases.For the initial properties of the CME itself, this includes the velocity of the CME front, v F , the mass, M CME , the half face-on width, AW, the half edge-on width, AW ⊥ , the yaw, the aspect ratio of the axis, δ Ax , the aspect ratio of the cross section, δ CS , the internal temperature, T CME , the magnetic field at the center of the flux rope, B, the corresponding flux rope parameter B 0 , the EC flux rope model parameters τ and C, the adiabatic index, γ, an expansion parameter that varies between initial self-similar expansion or convective expansion, f Exp , the drag parameter, C d , and the initial distance.For the SW, this includes the number density, n SW , radial velocity, v SW , total magnetic field strength, B SW , and temperature T SW .For a more thorough description of these inputs see Kay, Mays, and Collado-Vega (2022).
We use a fast and an average CME as an example of the interaction.We note that the speed of the average CME is slightly above average, but we refer to it as average for simplicity.The fast and average CME properties are the same as used in Kay, Nieves-Chinchilla, et al. (2022).For the quiescent SW at 1 au, we use the default values from Kay, Mays, and Collado-Vega (2022), which differ slightly from those in Kay, Nieves-Chinchilla, et al. (2022) where they were tuned to match the MHD simulation.
To start, we run only the ANTEATR component of OSPREI.The third row shows the edge-on angular width of the CME (black) and the axial aspect ratio (maroon).The fourth row shows the yaw angle of the CME (black) and the rotational speed (maroon).The fifth row shows the Alfvénic Mach number of the shock (black) and the number density within the sheath (maroon).Finally, the bottom panel shows the signed logarithm of the acceleration from various forces: the total radial acceleration of the CME (black), the bulk drag (purple), the pancaking from drag (maroon), the axial magnetic forces (dashed light blue), the magnetic pressure (dashed green), and the thermal pressure (dashed yellow).The fine dashed, horizontal line indicates the zero value.
shaded).The second row shows the speed of the shock (maroon) and the CME front (black), and the background SW speed (blue dashed line).The third row shows AW (black) and δ Ax (maroon).In general, the behavior of the axial and cross-sectional width and aspect ratio are similar, which is why we only show the axial values.The fourth row shows the yaw (black) and the angular velocity (maroon).The fifth row shows the Alfvénic Mach number of the shock (M A , black), and the number density in the sheath (n Sh , maroon).Finally, the bottom row shows accelerations (calculated from the forces) acting upon the CME.This includes the total acceleration in the radial direction (black), the bulk drag (purple), the pancaking or squashing drag (maroon), the axial magnetic forces (dashed light blue), the magnetic pressure (dashed green), and the thermal pressure (yellow).The accelerations are shown in a "signed logarithm" format such that for some arbitrary a and any values of |a| less than 1 are set to a.This unscaled region for absolute values less than one is shaded light gray.For the signed logarithm, positive values indicate an outward acceleration in the radial direction, either in the sense of a bulk motion or an outward expansion of the cross section.
For cases including a HSS, we show the location of the HSS within the standardized profile figure sing the same alternating bands of light blue and orange as in Figure 1.Within the top panel, we shade the background to show the radial extent of the different regions at each time step/CME front distance.For a static HSS, the regions will appear as perfectly horizontal lines within the panel.When the HSS propagates outward the regions become inclined.The three regions at the farthest distance (along the y-axis) correspond to the SIR, the next closest to the Sun the plateau, with the tail being closest to the Sun.All five regions may not be visible at all distances if the HSS extends beyond the range of the panel.Within the other five panels, the background is shaded along the x-axis of the panel according to the region of the CME front at that time/distance.

Fast CME
Figure 3a shows results for the fast CME in a quiescent background.Figure 3b shows results including a static HSS with the front at 0.8 au, generated from an 8 × 10 10 km 2 CH.Figures 3c and 3d show results for a time-dependent version of the same size HSS, with the CME either entering (c) or exiting (d) the HSS at the same radial distance as the static case.
Figure 3a shows that in the quiescent background, the fast CME drives a shock and, throughout its propagation, the CME gradually decreases in velocity, slightly expands in size, pancakes, and accumulates a sheath.We see no yaw rotation due to the longitudinally-uniform background SW.Looking at the forces, we see that the drag force acts to decelerate the CME, as expected given that it propagates faster than the background SW.The axial magnetic forces also slow down the radial motion.The thermal pressure continually causes the CME to expand outward throughout propagation.The magnetic pressure oscillates between expansion and contraction as the CME magnetic pressure changes relative to that of the background SW.The CME conserves magnetic flux so its magnetic pressure changes as the CME volume changes.This case is expanding in angular width but pancaking in aspect ratio, producing competing effects that cause the initial variations in the direction of this force.The drag forces remain high (absolute value > 100 cm/s 2 ) through the duration of the simulation but the internal forces drop below 1 cm/s 2 by about 0.6 au.
In Figure 3b, the CME mostly behaves the same as in the quiescent background until the nose reaches the back of the tail of the static HSS at about 0.18 au.As the CME front enters the HSS tail, the drag between the CME and the background SW decreases, which can be seen in the acceleration panel.This causes there to be less decelera tion in the bulk velocity, less expansion and pancaking, and a weaker shock and less accumulation of material into the sheath.The CME continues experiencing weaker drag throughout the plateau (0.3-0.73 au) until it emerges into the SIR.At this point, the drag rapidly increases because the background velocity begins approaching the quiescent value while the background density is high in this interaction region.The CME begins decelerating, expanding, and pancaking at a faster rate while driving a stronger shock and accumulating more SW material into the sheath.
We do see a slight rotation in the yaw before the front enters the tail at 0.18 au.At close distances, one flank of the CME experiences the tail of the HSS while the other is in the quiescent SW.The side of the CME at higher longitudes first experiences the higher speed HSS background, so it experiences less drag on this side, causing it to rotate in a negative or clockwise direction when observed from above from solar north.As the CME pushes into the HSS it reaches a point when both flanks are embedded in the fast wind and the torque disappears.The continued radial motion eventually brings the higher longitude flank beyond the SIR so that it is embedded in the quiescent SW while the other flank is still in the core of the HSS.At this point, the rotational velocity becomes positive and begins slowing down the clockwise rotation.We only see a few degrees of total yaw rotation.Intuition, based on general remote-sensing CME observations and more sophisticated simulations, suggests that we certainly do not expect there to be tens of degrees of yaw rotation, but, in later sections, we will explore whether these relatively small rotations are at all meaningful for in situ predictions.
Using a static HSS versus a time-dependent background does not fundamentally change any of the ways in which the CME interacts with the different regions of the HSS, it only changes the distance at which those interactions occur and the duration within each region.For a fast CME, the drag will decrease as the CME enters the tail, be at a minimum within the plateau, then rapidly increase as the CME reaches the SIR and eventually the upstream SW.
When the CME nose is behind the HSS the flank at higher longitude experiences less drag and the yaw rotates clockwise.As the CME moves toward and beyond the HSS front, the flank at lower longitudes experiences less drag and the yaw rotates counterclockwise.
In Figure 3c, the CME starts behind the HSS, propagates into the tail then remains in the plateau for the remainder of the propagation to 1 au.This causes it to spend more time in regions with weak drag and never reach the SIR where the properties begin "recovering" toward the quiescent values.In contrast, for the case in Figure 3d, the CME front begins in the plateau and eventually moves to within the SIR as the CME approaches 1 au.The CME experiences less drag for the first 0.6 au of propagation so it maintains a velocity well above 1,000 km/s during this time and drives a stronger shock (relative to the quiescent case) at farther distances.The results for the interaction of the CME with a time-dependent HSS are analogous to those seen with a static HSS in Kay, Nieves-Chinchilla, et al. ( 2022), the time-dependence just changes the timing of the different behaviors.In Section 5 we quantify how parameter space variations in the HSS size and initial distance correspond to timing differences that affect the final properties at 1 au.

Average CME
Figure 4 shows the same as Figure 3, but for the average CME.This CME moves slightly faster than the quiescent background so we see a slight deceleration due to drag in Figure 4a, but the magnitude is much less than for the fast case.The CME exhibits expansion and pancaking and is able to form a weak shock and accumulate some SW material.The sheath duration is actually comparable to that of the fast case, but the average case has much less compression within this region.
For the static HSS (Figure 4b), the CME moves slower than the background SW within most of the HSS.The shock speed drops to zero midway through the HSS tail when the HSS speed exceeds that of the CME front.At this point, the CME continues with whatever sheath it has accumulated but no additional material is gained.The bulk drag force becomes positive and works to accelerate the CME.The CME velocity does increase by about 10 km/s, but this is not noticeable on the scale of Figure 4b.Once the CME exits the HSS, which it can only do because this is a static case, the background SW speed decreases and the drag force returns to decelerating the CME.The average CME exhibits the same clockwise rotation followed by counterclockwise motion as the fast CME, but the magnitude is less than a degree.
For the time-dependent background cases, the initial CME speed is barely faster than the propagation speed of this HSS.In Figure 4c the CME enters the tail of the HSS and remains there for the duration of the propagation.At 0.58 au, it eventually reaches the point where the HSS speed exceeds that of the CME and the shock disappears.The velocity differential is minimal so we do not get a significant acceleration from the bulk drag.In Figure 4d, the CME begins near the outer boundary of the plateau so no shock initially forms and the drag acts to accelerate the CME.The CME reaches a peak speed of 660 km/s before starting to decelerate as it enters the SIR.Once the SW speed drops sufficiently the shock forms and the CME sheath begins to grow.In this case, we consistently see a large pancaking force.Typically, pancaking is expected from a fast CME running into the slower upstream SW, we see that it also can occur when a slower CME is pushed from behind by fast SW.As with the fast CME, these time-dependent results largely reproduce the behavior found for the static case in Kay, Nieves-Chinchilla, et al. ( 2022), although we see some added effects when a CME is pushed from behind due to the addition of the two-point drag sampling the background both upstream and downstream of the CME.

An Extreme Illustrative Case
We present another purely theoretical case where we have tuned the simulation to maximize the effects of a HSS pushing a slow CME from behind.The parameters for this case are listed in the "slow" column of Table 1.It uses the same quiescent SW parameters as the previous theoretical cases, but we initiate MEOW-HiSS with a HSS front at 0.05 au and originating from a 500 × 10 10 km 2 CH.This embeds the CME toward the back of the SIR at the start of the simulation.Figure 5 shows the results for this case in the same format as Figures 3 and 4.
The velocity panel shows that the CME speed is initially lower than that of the HSS and no shock forms.The accelerations show a strong forward push from the bulk drag, but the differential drag causes a strong compression of the CME.The CME speed slowly increases toward that of the background SW, finally reaching it around 1 au.As the CME speed increases, both drag accelerations slowly decrease.
We see that the CME AW initially increases due to the magnetic and thermal pressure.As the aspect ratio changes acceleration the magnetic pressure changes directions several times below 0.3 au.The thermal pressure continues to correspond to acceleration, but the magnetic pressure is much more sensitive to the aspect ratio of the cross section, which is continually decreasing at this point.Beyond 0.3 au we see a slow decrease in the AW.The CME increases radial speed but the expansion velocities remain lower than the values corresponding to self-similar expansion, so both the AW and the aspect ratio decrease.
We see a small counterclockwise rotation in the yaw.The rotational velocity is largest near the Sun, where there is the most imbalance of forces with one edge embedded in the CME and the other in the quiescent SW.As the CME propagates the HSS moves toward the nose of the CME and the force imbalance decreases and disappears.The rotational velocity decreases over time to conserve angular momentum as the CME expands.
To better visualize the interaction between this slow CME and the HSS we show contours of the total velocity within the equatorial plane in Figure 6.The CME is oriented so that the toroidal axis of the CME lies within the equatorial plane.The white line outlines the outer edge of the CME flux rope.We emphasize that while our results can be visualized similarly to the common method for many hydrodynamic and MHD simulations, OSPREI is not a full 3D simulation where the plasma properties are self-consistently evolved across grid cells.MEOW-HiSS does provide plasma properties for the time-dependent background SW in the equatorial plane, but the CME is simply added on top within these images.OSPREI only calculates forces at a few specific points and uses them to change the speed, size, and shape of a torus, then uses conservation laws to evolve the internal properties.We then just use these to replace the portion of the equatorial plane corresponding to the CME with the appropriate speed.
The top row shows the evolution of the slow CME interacting with the HSS and the bottom row shows simulation results using the same CME inputs without any background HSS.The pale blue dots represent the location of four synthetic satellites at 0.12, 0.35, 0.70, and 0.99 au (25, 75, 150, and 213 R ⊙ ).We use FIDO to simulate profiles at these points, which we discuss later on.The three panels within each row correspond to the first contact at each of the three outermost satellites, with the corresponding simulation time listed above each panel.We see that the CME interacting with the HSS arrives at each satellite earlier than the CME in the quiescent background.The difference in the arrival time increases as the CME propagates farther, reaching a difference of 22.8 hr at the farthest satellite close to 1 au.
Figure 6 also highlights the differing evolution in the CME speed and shape.Even within the leftmost panel, there is a noticeable difference in speed between the two cases, with the HSS case being faster, and this difference continues to increase over the temporal domain of the simulation.We note that both cases do have the CME accelerating, as the slow CME is slower than even the quiescent background.For both cases, we see that the CME is faster near the edges than at the nose, which corresponds to the pancaking of the CME.We also can see the difference in the width of the CME cross section in the radial direction.At 0.99 au, the CME is twice as wide in the quiescent simulation than in that with the HSS.Combining this with the difference in speed leads to a much shorter duration when the CME is pushed by the HSS.The yaw rotation is minimal for this case, so no significant rotation can be seen in Figure 6.
Figure 7 shows the in situ profiles corresponding to the synthetic satellites shown in Figure 6.The red family of colors (pink, red, maroon, and brown) corresponds to the case with the HSS, and the blue family of colors (light blue, blue, navy, and dark blue) to the quiescent case.The thin portion of each line represents the ambient SW and the thick portion is the CME flux rope.All profiles are shown on the same scale to facilitate comparison.
From top to bottom, Figure 7 shows the normalized total magnetic field, B R , B T , and B N field components, velocity, temperature, and number density.Since most plasma properties vary by orders of magnitude between the Sun and Earth, at each satellite we show the profiles normalized by the maximum value within the flux rope portion of the profile for the case with the HSS.For example, for the first satellite, we find the maximum in the thick portion of the pink profile and use this to normalize both the pink and light blue profiles.The same approach is used at each subsequent satellite.The velocity is the only parameter shown without normalization.
At the first satellite (0.13 au, pink and light blue), the profiles are nearly identical as minimal evolution has occurred.The most obvious difference is in the SW speed and temperature following the CME, which shows the HSS in the pink case.The pink CME does show a slightly shorter duration and temperature by this distance.By the second satellite (0.35 au, red and blue) we already see significant differences in the synthetic profiles.
The red HSS-driven case is shifting toward an earlier arrival.The HSS accelerates the CME to higher speeds and compresses it.This compression enhances the magnetic field and temperature, though at this distance the density remains similar.All three properties vary with the CME size, but the density is based on the total volume whereas the magnetic field and temperature depend more strongly on the cross-sectional width in the model.We note that the magnetic field increases in intensity, but the rotation in the individual vector components remains the same.
These effects are further magnified by the third (0.70 au, maroon and navy) and fourth (0.99 au, brown and dark blue) satellites.The arrival times continue to separate, the CME in the quiescent background actually arrives at the third satellite around the same time as the HSS-driven CME reaches the fourth satellite.The HSS-driven case continues to shorten (relatively) in duration and enhance in internal properties, including density.These profiles highlight the extreme effects that one may expect from the interaction between a CME and HSS by the time it reaches 1 au.

Parameter Space Results
As done in Kay, Nieves-Chinchilla, et al. ( 2022), we perform a set of observations sampling parameter space to better understand the effect of different HSS properties on the CME-HSS interactions.We vary the CH area between 200 × 10 8 km 2 and 2,500 × 10 8 km 2 and the initial HSS front between 0 au and 1.2 au.We use a resolution of 200 × 10 8 km 2 in the area and 0.1 au in the distance.For each combination of CH area and HSS distance, we simulate both a fast and an average CME using the same CME properties as in Section 3.
Before looking at the parameter space results, we consider what this variation in CH area and HSS front distance corresponds to with our time-dependent simulations.Figure 8 shows contours of different properties versus the CH area and initial HSS front distance.Figure 8a shows the region of the HSS in which the CME front is located at the start of the simulation.Region 0 corresponds to the upstream SW, 1-3 to the three SIR regions, 4 to the plateau, 5 to the tail, and 6 to the downstream SW.The difference in HSS size can be seen as the range of front distances corresponding to regions 1-5 increases with CH area.We note that the step-like features at the region boundaries are due to the coarse resolution of the grid.
Figure 8b shows the longitude at which the HSS front is initially at 1 au given each pair of HSS input parameters.
The HSS follows a Parker spiral due to the solar rotation so, at a single point in time, the distance of the HSS front decreases as one moves toward higher longitudes (using Carrington-like or Stonyhurst-like coordinates).We find minimal variation with the CH area as the curvature of the spiral is essentially driven by the solar rotation.The corresponding longitudes of the 1 au fronts vary between roughly 50° east and 10° west of the source location.
The bottom row shows the region in which the CME nose is embedded when it reaches 1 au, analogous to panel (a). Figure 8c shows results for the fast CME and Figure 8d for the average CME.For the fast CME, we find that the CME will have reached the upstream quiescent SW if the initial HSS front distance is below about 0.3 au.It only remains in the downstream quiescent field for the smallest CH areas and farthest initial distances.In contrast, almost none of the average CME cases reach the upstream SW.The average CME moves slightly faster than the bulk motion of the HSS, but if the HSS starts sufficiently far out we do see a section of parameter space where the average CME cannot catch up to it and remains in the downstream region.

Comparison With Previous Results
Table 2 lists the results of the parameter space study for both the fast and average cases.The left column lists the different outputs from OSPREI for both CME and sheath properties (divided by the horizontal line).The two vertical lines separate results for the fast and average cases.Within each region the first column lists the final value at 1 au for the quiescent case, followed by the most negative and most positive change in that value over the HSS parameter space exploration.We do not include a detailed analysis of each row within this table but the values are included for completeness.The general results are the same as in the analogous study within Kay, Nieves-Chinchilla, et al. (2022), just with small changes in magnitude due to the new time-dependent nature of the background HSS.
We briefly look at contours of the variation in parameter space for several outputs.Again, the results are fundamentally the same as Kay, Nieves-Chinchilla, et al.From top to bottom, the rows show the normalized total magnetic field, B R , B T , B N , velocity, temperature, and number density.The pink, red, maroon, and brown profiles correspond to results using an HSS background, and respectively to increasing satellite distances (25 R ⊙ , 75 R ⊙ , 150 R ⊙ , and 213 R ⊙ ).The corresponding profiles for the ambient background are respectively the light blue, blue, navy, and dark blue lines.The thicker portion of each profile corresponds to the CME and the rest to solar wind (either ambient or HSS).For each satellite distance, both the ambient and HSS profiles have been normalized using the maximum value within the CME for the HSS case (e.g., both red and blue have been normalized by the maximum in the bold part of red).For the magnetic field components, we use the maximum in the total magnetic field strength for normalization.
that work), but we include a few examples for illustrative purposes.Figure 9 shows several contour plots for the fast case (left column) and the average (case).The dashed black lines indicate the final region in which the CME front is located at 1 au (the same as Figures 8c and 8d).From top to bottom, Figure 9 includes the transit time, axial aspect ratio, magnetic field strength in the sheath, and the yaw.We will go into further detail for the yaw in Section 5.2 as it was not studied in Kay, Nieves-Chinchilla, et al. (2022).The color map is set so that white indicates little to no change, red is a positive change, and blue is a negative change in the final value of a parameter.
We find that the HSS can only cause a decrease in the transit time.The HSS either causes an acceleration or a decrease in the amount of drag so that the CME moves faster and reaches 1 au quicker.This causes a change of up to 9 hr for the fast CME and 5 hr for the average CME.These largest values occur when the CME is located in the plateau of the HSS for the entirety of its propagation.We saw this effect in Kay, Nieves-Chinchilla, et al. (2022), but now the CME can remain in the plateau for a longer time as the HSS also propagates out radially.
For the fast CME, we mostly see an increase in the axial aspect ratio (Figure 9c).This occurs in the plateau where the drag forces are minimal and the CME experiences less pancaking.For the average CME we see a decrease in pancaking toward the back of the plateau, but also see an increase toward the SIR.The drag forces increase when the CME emerges from the plateau to the SIR and the aspect ratio quickly begins decreasing.For the average case, the aspect ratio falls below the quiescent value for quite a few cases.The fast case experiences the same effects but the magnitude of the changes in this post-plateau "recovery" phase cannot fully counteract the effects from its time within the plateau for all but the closest initial HSS distance.In general, the internal CME properties mirror the evolution of the aspect ratio, there is less contraction and more expansion in the plateau, leading to lower/weaker plasma properties.
Once the CME reaches the SIR the drag increases, more pancaking and contraction occurs and the internal properties begin enhancing.The sheath also evolves in a systematic manner.The shock is weaker in the plateau since the background speed is closer to that of the CME.This causes less compression, so we see weaker values, such as shown for the magnetic field strength of the sheath driven by the fast CME in Figure 9e.At the SIR, the drag suddenly increases as the SW velocity decreases and we see a stronger shock and more compression.As the CME moves toward the upstream region the shock and corresponding sheath properties approach the values of the quiescent case.The average case cannot form a shock, or even a compression region in front of the CME, within the plateau as it propagates slower than the background SW in this region.The gray shaded region in Figure 9f corresponds to input parameters that have a CME unable of driving a shock at 1 au.Outside of this region, we find an enhancement region when the CME reemerges into the SIR or upstream region with faster speed than in the quiescent case, leading to an increase in the magnetic field strength, and a weaker region when the CME is within the HSS tail at 1 au.

Yaw
Since the quiescent case experiences no rotation in the longitudinally uniform background, the bottom row of Figure 5 shows the final yaw angle.For all CMEs, when the HSS front is initially at a distance close to the Sun then the lower-longitude flank of the CME is embedded in the HSS and the higher-longitude flank in quiescent SW.This causes more drag on the higher-longitude flank and a positive rotation in the yaw.Alternatively, if the HSS front is farther out then the opposite can happen, with only the higher-longitude flank embedded in the HSS and the lower-longitude flank experiencing more drag, resulting in a negative rotation.The exact distance to maximize a negative rotation varies with HSS size.In between these two extremes, we find a region where the CME is fully embedded within the HSS and little to no rotation occurs.A larger region of parameter space corresponds to this fully-embedded lack of rotation for the average CME because the CME is smaller in size than the fast CME.We also see a region where the HSS is sufficiently  far out that the CME is fully embedded in the downstream quiescent region, so again no rotation occurs.We see that the HSS can cause the fast CME to rotate by about ±7° for the fast case, and ±2° for the average case.The increased drag on the fast CME, due to its higher speed, causes there to be more rotation.
Figure 10 further explores a maximally rotating fast case.The left panel shows contours of the radial velocity in the equatorial plane, analogous to Figure 6.We show a single time step when the CME reaches 0.99 au.We also include the CME-driven sheath in Figure 10a using the sheath width and CME boundary and approximating a constant speed within the sheath.Again, this is a visualization of an analytic simulation, not a full fluid simulation, which is why the sheath abruptly ends at the flanks rather than smoothly transitioning back to ambient SW values.We include nine synthetic satellites, one at L1 and the rest separated by increments of 10°.All satellites are at 0.99 au.The colors match the in situ profiles shown in Figure 10b with the red satellite along the simulated Sun-Earth line, purple to pink representing increasingly positive satellite longitudes, and orange to yellow representing increasingly negative longitudes. Figure 10b has the same format as Figure 7 but shows physical units rather than being normalized.The short dashed portion of each profile corresponds to the ambient SW, the long dashed to the sheath, and the solid to the CME flux rope.
This case has a 5.1° rotation in yaw at the time of the first CME (not sheath) contact at L1.This can be seen by comparing to the satellite locations in Figure 10a.When the CME is just reaching the synthetic Sun-Earth L1 point, the satellite at −40° is already within the flux rope but the satellite at 40° is still within the sheath.The parameter space exploration showed that the introduction of a HSS could change the transit time by 9 hr for the fast CME, but this was only determined at the CME nose and did not account for any spatial variations.For this case, the CME reaches the satellite at −40° after 33.1 hr and the satellite at 40° at 38.4 hr.The corresponding sheath arrival times are 27.3 and 30.5 hr and the CME (not including shock) durations are 10.6 and 15.3 hr.Since the CME reaches the satellite at 40° later it has more time to develop a sheath and expand than at −40°.The forces are causing the CME to pancake, notice again the variation in velocities within the CME, but it is still overall expanding in size while the aspect ratio changes, leading to the increase in duration at the later contact.
The profiles in Figure 10b further illustrate this behavior.In our simulation, the temperature and number density are uniform throughout the CME.This means that all satellites see the same value if they are within the CME at a given time step/front distance.This is certainly an oversimplification of a real CME that likely shows much more variation along its toroidal axis.We still expect, for expanding CMEs, that earlier impacts should have higher density and temperature on average, but there will be more local variation about the general behavior.
While the density and temperature are uniform, our model does have some spatial variation in the flux rope magnetic field.The strength at the toroidal axis is uniform along the axis but the local magnetic field depends on the distance from the axis and the current aspect ratio of the cross section.Averaging over the full CME portion of the profile, we see the stronger magnetic field for earlier impacts.Within each profile, the field tends to increase as the satellite approaches the toroidal axis, then decrease as it moves away.
We see some variation in the direction of the magnetic field vector across different satellite longitudes.The vector is shown in radial-tangential-normal (RTN) coordinates.Even if there was no yaw rotation we would expect some longitudinal variations as the normal to the toroidal axis is, in general, only parallel to the radial vector at the CME nose.The normal component will not vary with longitude, but the equatorial projection of the vector will separate into both a radial and a tangential component, the exact nature of which depends on the CME shape, size, and distance.When the yaw is 0° the relative magnitude of the radial and tangential components will vary symmetrically for positive and negative longitude.
In this case, the toroidal axis lies within the equatorial plane so for a satellite impacting exactly at the nose, we would expect to first see only normal field with no radial or tangential component as the satellite begins sampling the fully poloidal field.Our satellite has an impact parameter of zero, if it did not go exactly through the center of the cross section then we would expect the poloidal field to have a radial component.As the satellite moves inward the poloidal field decreases, while the toroidal field, which points tangentially, increases.At the center of the cross section the field is entirely toroidal and tangential.As the satellite exits the CME the process is reversed, and if the CME is expanding all magnitudes will slowly decrease over time.
For this case, the positive yaw rotation cause the normal to the toroidal axis within the equatorial plane and toroidal direction to be more aligned with satellite radial and tangential directions as one moves toward more positive longitudes.The higher longitude profiles show very little radial component and, despite having impact occur later when B 0 is weaker, they have stronger toroidal field than the cases at lower longitudes.We see no variation in the normal field beyond that from the slight decrease in B 0 over time.
We see some variation in the speed with the magnitude increasing for impacts toward the flanks since this CME is pancaking.We see the largest values for the satellites at positive longitudes where the CME is not actually interacting locally with the HSS, but this is a result of oversimplifications in the model.The CME interacts with the HSS over time, continually causing more pancaking, which would lead to an increase in speed wherever the interaction is occurring.The model approximates the change in shape (and corresponding velocities) as symmetric about the nose because it is not yet sophisticated enough to handle that level of asymmetry in the flux rope shape.In a real scenario, we expect that there would be more distortion of the CME from its rigid (beyond aspect ratio) axial shape and more local variations so we would not see the enhanced speed at high longitudes.
The sheath portion of the profiles behaves similarly to density and temperature within the flux rope.The sheath evolves uniformly so its properties have the same magnitude at all locations at a given time/CME front position.
The magnetic field magnitude is uniform, but the orientation does vary according to the local orientation of the flux rope since we assume that the SW field drapes about the CME.We see that the normal component remains the same but the component in the equatorial plane rotates between the radial and tangential directions.

Observed Case
Our final example combining OSPREI and MEOW-HiSS is an observed case, previously simulated with OSPREI and a quiescent background in Ledvina et al. ( 2023)-hereafter L23.This CME erupted at 21:00 UT on 26 January 2022 from an active region on the far side of the Sun.As the event is back-sided, we see no evidence in Solar Dynamics Observatory (SDO) imagery, but Extreme UltraViolet Imager (EUVI), onboard Solar TErrestrial RElations Observatory (STEREO)-A, does observe post-eruptive arcades approaching the eastern limb.This CME is associated with the arrival of an event at Parker Solar Probe (PSP; Fox et al., 2016) on 29-30 January 2022.The PSP observations show a HSS stream following the CME, and we can still see a CH to the left of the source region when it rotates into the STEREO and SDO fields of view.We use the observations of the CH in SDO to measure the CH area (as described in Kay et al., 2023) and estimate the initial distance at the time of the CME eruption using the rate of solar rotation.
Before analyzing the in situ profiles, we present contours of the speed to visualize the CME-HSS interaction.
Figure 11 shows these contours in the meridional (left) and the equatorial (right) planes at three different time steps.The meridional cut is taken along the plane through the CME nose.The blue dot represents the location of PSP.Earth (not shown) is located roughly 150° clockwise from the longitude of PSP at the start of the simulation.
To visualize the HSS in the meridional plane we approximate it as having the same angular width and variation in the latitudinal direction as it has in the longitude in the equatorial plane.This assumption is only used for visualization purposes and is not a part of MEOW-HiSS generally.
The HSS starts just to the east of the CME, initially only interacting with the eastern flank of the CME (Figures 11a  and 11b).By the time the CME-driven sheath reaches PSP (Figures 11c and 11d), the HSS is catching up to this relatively slow CME so that the HSS front is about halfway between the Sun and the back of the CME along the radial line cutting through the CME nose.When PSP exits the CME (Figures 11e and 11f), the back portion of the CME is fully embedded in the HSS.
We now compare the observed and modeled in situ profiles.Figure 12 shows the PSP observations in black and include, from top to bottom, B, B R , B T , B N , v, T, and n.We identify a shock arrival at 0:00 UT on the 29th and the start of the CME (flux-rope-like component) at 9:00 on the same day.The exact boundary between the rear of the CME and front of the HSS is open to interpretation but we approximate it at 13:00 UT on 30 January 2022.One could argue for a boundary at 0:00 UT on 30 January based on v, n, and T, but the total B remains very smoothly varying beyond this time.We have also investigated the plasma beta (not shown), which essentially mimics the temperature profile, increasing above unity when the temperature begins to rise, and remaining high.The plasma beta is often a useful metric for determining flux rope boundaries but as a HSS pushes on a CME from behind, we expect heating will occur in its trailing portion.Furthermore, if the CME and HSS have oppositely directed magnetic field then reconnection can occur, mixing and heating the plasma as well as raising the beta.We suspect that the region between 00:00-13:00 UT is likely an interaction region between the CME and the SIR of the HSS.As we will see, the modeled flux rope and the following SIR do have opposite signs for B T so it is plausible that reconnection occurs and mixes the plasma.
Table 3 lists the OSPREI input parameters for both L23 and this work (labeled K23).In both cases, a full OSPREI simulation, including ForeCAT, was performed.The focus of L23 was how different magnetograms affect OSPREI results for four different events.They tuned each event to the in situ observations using a synchronic magnetogram, which updates a synoptic magnetogram using the most recent daily observations from the Earth's perspective.As this is a back-sided event, the synchronic magnetogram is not necessarily the best choice, so we instead use the synoptic magnetogram.We also use the most up-to-date version of OSPREI, incorporating all the changes in the ANTEATR component.This includes the sheath development, the 2-point drag, the HSS background, and the yaw rotation.We re-tune the OSPREI results by hand to find an approximate best fit with the new version incorporating all the changes.
Figure 12 includes the results from both L23 (yellow) and this work (K23, red), as well as several intermediate results to show how the different model improvements affect the in situ results.The line styles within each profile represent the same breakdown between ambient SW, sheath, and flux rope as in previous figures.The vertical dashed lines show the sheath, CME ejecta start, and CME ejecta end.We note that we are simply showing the model evolution, not trying to identify the most important factor in improving the fit to the in situ profile.The improved profile is a compound effect of all the changes and was found as such, using all the new features at once.When building to this case one step at a time the magnitude of each change depends on the order in which we introduce the changes, and the relative importance would likely vary between different CMEs.
In building to the new version, the first step (K23 0 , orange profile) is just a change in input parameters but no change in the OSPREI capabilities.The most significant change is in the internal temperature of the CME, here  we use a much lower value to decrease the excessive expansion, resulting in a shorter duration and a better match to the observations.The modeled in situ temperature is slightly lower than the observations for both the L23 and K23 0 cases, but the observations show many fluctuations that will never be captured with a uniform temperature model.The K23 input parameters also noticeably change the radial velocity profile in the corona (initial speed v 0 , the distance that linear acceleration starts, R 1 , the distance that linear acceleration ceases, R 2 , and the final coronal speed, vF), and slight changes in the initial position (latitude, longitude, and tilt), CME mass (M CME ), edge-on width (AW ⊥ ), internal magnetic field (B 0 ), adiabatic index (γ), and expansion factor (f Exp ).The net effect of these, when using the previous version of OSPREI, is a later arrival time, shorter duration, and higher internal magnetic field strength, density, and speed.Depending on the metric, this result is arguably not a better fit to the observations than L23, but we emphasize that this is just an interim step.This set of inputs is only optimal when all the changes are incorporated.
We then include the CME-driven sheath (K23 1 , light blue profile).We see minimal changes from K23 0 in the flux rope portion of the profile, mostly just a slight shift toward an earlier arrival time.However, the profile does now include a sheath, which fairly well reproduces the observed sheath duration, total magnetic field strength, velocity, and temperature.We see excessive compression in the density and the vector components are oversimplified as we are using a smooth profile to represent a highly turbulent region.Here is a clear example of where an interim step could produce a better fit with a different set of inputs, such as lowering B 0 in this case.Again, our focus is not in quantifying the importance of each step as we believe they are all important, but rather simply illustrating how OSPREI has evolved from the L23 case.Next, we update the drag to the new version (K23 2 , blue line), followed by including the HSS (K23 3 , maroon line) as well as the yaw rotation to reach the current OSPREI version (red line).These three lines are hard to distinguish between in Figure 12.In particular, this CME only rotates 2.5° so the maroon line is fully obscured by the red line.The only instance in which the maroon line is marginally noticeable, even when zooming in, is in the middle of the flux rope for B R as it falls slightly above the red line.We do see some difference between the blue and red cases, particularly near the rear of the CME.The change in drag moves the total magnetic field profile to a near exact match to the observations and the vector components do an adequate job of reproducing the average behavior of the observations.The fit is worse near the rear of the CME, where we expect there to be signatures of the interaction with the HSS.The updated drag causes a slightly earlier arrival and shorter duration than in K23 1 but the speed is still too low near the CME rear.
Including the HSS causes the CME duration to finally shorten to match that of the observations at least when including the "mixed" region of the observations.We also see an increase in the speed at the back of the CME, not only an increase relative to K23 2 but an actual increase within the profile between the middle and end.This is the signature of the HSS accelerating the CME.While we do capture some increase, the extent of it is still much weaker than that seen in the observations.Overall, the new changes in the model allow us to find an excellent match to the timing and properties of the sheath and front portion of the CME, but the model may still be a bit too simplistic to fully capture the dynamic interaction, including potential reconnection effects, at the rear of the CME.
We run an ensemble of OSPREI simulations varying the initial CME latitude, longitude, tilt, both angular widths, both aspect ratios, radial velocity, mass, internal magnetic field strength, and temperature, as well as the CH area and initial HSS front distance.This represents all of the input parameters with significant uncertainties in their initial values.As this is a highly multi-dimensional parameter space, it is possible that we were not able to find their optimal values in our initial manual exploration.We run an ensemble of 200 simulations with all of the chosen parameters varying simultaneously.We specify a range for each varied parameter, which sets the width of a Gaussian distribution centered about the ensemble seed, and from which the varied inputs are randomly drawn.
Our ensemble produces 200 different in situ profiles, all of which include portions of the ambient SW, sheath, and flux rope.To quantitatively define which profile is "best" we must develop some form of metric combining information about the timing and accuracy representing various plasma properties.We calculate the metric using the mean absolute fractional error in the total magnetic field, velocity, number density, and temperature between the observed sheath start and the end of the flux rope.We sum the fractional error for all four properties and add the absolute error (in days) for the arrival of the sheath front, CME front, and CME end.This metric is somewhat arbitrary and we have no justification for it being the best for all cases, but it is certainly a reasonable combination of the errors we wish to minimize.
Using this metric, we have identified the best fit. Figure 13 shows each of the ensemble members in gray, the seed in red (corresponding to the values listed in Table 3), and the best fit in purple.We note a minimal difference between the best fit and the seed case for this case.All panels, other than the temperature, show a slight shift in the flux rope profiles, but the best fit and seed profiles largely overlap.We do see a slight increase in the velocity in the best fit case within the mixed region at the back of what we have labeled as flux rope.Had we set the flux rope end to a different time, we would have likely found a different best fit.Further study needs to be performed to identify the optimal metric (and boundaries) for quantifying the quality of fit across a variety of observations.

Discussion
This work presents the coupling of two strategically simplified models in an attempt to reproduce an interaction that is very complicated in the real world.There are certainly a number of oversimplifications, but we see promising results in a variety of applications.We believe that the combination of OSPREI and MEOW-HiSS can be very powerful in the right circumstances, but some caution should be applied before blindly applying it to every case.
The biggest caveat is that OSPREI includes almost no local variations, the evolution all happens on a global scale.This is particularly an oversimplification for the yaw rotation.A CME is a coherent structure to some extent, but realistically information about what is happening at one flank does not immediately affect the behavior at the other flank and result in a solid body rotation.In the example of one flank embedded in an HSS and the other in quiescent SW, we would expect differing drag to cause the two flanks to propagate out differently so that one will move farther out relative to the other.Realistically this does not generate a net torque about the nose that rotates the CME as a rigid structure, it is a series of localized changes.Owens et al. (2017) show that a CME should be considered less and less of a coherent structure as it propagates outward while expanding and experiencing changes in the internal Alfvén speed.This suggests that the solid body rotation will become a worse approximation at farther distances.The solid body rotation within OSPREI should be thought of as a mathematical simplification of much more complicated behavior that can hopefully provide a useful approximation, much in the manner that a one-dimensional drag equation is used for the propagation of a CME, or MHD is used to represent the collective behavior of an untenable number of individual particles.We believe that OSPREI is largely providing a reasonable approximation of the overall behavior, but this needs to be confirmed with multipoint in situ observations.
With this in mind, we advocate intelligent application based on the relative location of a HSS-CME interaction and the observing satellite.If the interaction is happening close to the satellite's trajectory through the CME then we have relatively high confidence in the results.If the satellite encounter happens at one flank but the interaction occurs at the other we encourage caution, particularly in the in situ profile within the flux rope.As far as OSPREI and MEOW-HiSS can provide the general net evolution from a series of local changes then the arrival time may be fine, but we have already seen examples of the flux rope velocity seeming inappropriate for far satellite longitude with the fast CME in Figure 10.Using Figure 10, we estimate that the in situ profiles will be fine up to 10°-20° from the longitudinal boundary of the CME-HSS interaction, but this is a ballpark guess and needs further validation.
We also do not include any feedback between the rotating CME and the exterior SW.The yaw rotation should cause external compression and rarefaction, which will act to counteract and slow down the rotation.This means that our already small estimates of the rotation may likely be overestimates.This questions the necessity of actually including yaw rotation in future interplanetary CME simulations.We do not have reconstructions of CMEs in interplanetary space that account for any evolution in a yaw angle.Even if they did, this work suggests that the expected rotation would certainly be smaller than the uncertainty in the reconstruction.This work shows that while the CME-HSS interaction can certainly noticeably affect the interplanetary propagation of a CME, any large-scale rotations will likely only have minimal effects.

Conclusion
We have presented the coupling of two existing models, OSPREI and MEOW-HiSS.OSPREI simulates the Sun-to-Earth evolution of a CME, including a component modeling its interplanetary propagation, expansion, and deformation.Kay, Nieves-Chinchilla, et al. (2022) had used a static HSS background with OSPREI and found that the drag greatly decreases when a fast CME is embedded within the plateau and tail of a HSS.Once it reaches the SIR and upstream SW the drag suddenly increases and the CME properties begin approaching the values that would be expected in a quiescent background.Alternatively, a slow CME can be accelerated by a HSS.MEOW-HiSS models a time-dependent HSS background and can easily be coupled with OSPREI.We find that the same general results occur when a CME interacts with a time-dependent HSS but the motion of the HSS changes how long the CME spends within each region.
In addition to coupling OSPREI with MEOW-HiSS, we have made several improvements in the interplanetary component of OSPREI.These improvements take advantage of being able to use a higher dimensional SW background, as opposed to the previous option that only varied with radial distance.We have modified the drag force to better account for differences between the front and back of the CME, such as quiescent SW ahead and a HSS pushing from behind.We have also introduced a yaw-like rotation where the CME rotates about its nose based on a net torque from different drag forces on opposite flanks of the CME.
We present a series of studies with the coupled models.We first closely analyze the interplanetary evolution of both a fast and an average CME interacting with a HSS.This includes determining the propagation, expansion, and deformation of the CME, as well as the growth of a CME-driven sheath.We compare how these factors change between a quiescent background, a static HSS, and a time-dependent HSS.We then present an extreme case where we have optimized the inputs to have the HSS give the CME a strong push from behind.
We next look into the variation in the CME and sheath results for different HSSs, as defined by the initial distance of the time-dependent HSS front and the area of the CH that generates the HSS.We find similar results to Kay, Nieves-Chinchilla, et al. (2022) with respect to how CME properties vary from the quiescent case depending on which HSS region they are embedded in when they reach 1 au.The time dependence simply causes changes in the relation between the initial HSS properties and the corresponding final region for the CME at 1 au.
We also analyze the yaw rotation for different initial HSS properties.When one CME flank is embedded in the HSS and the other in quiescent SW we see a rotation in the yaw due to the flank on the HSS side propagating faster than the flank on the quiescent side.If the CME is fully embedded within the HSS then little to no rotation is observed as the forces are relatively balanced.We find a maximal rotation of about 7°, which can cause a difference of about 5 hr in the arrival time and the duration of opposite flanks.We conclude that large-scale yaw rotation can only introduce minimal changes in observed CME behavior.
We compare OSPREI with MEOW-HiSS results for a CME-HSS interaction that was observed by PSP on 29-30 January 2022.This case was previously simulated by Ledvina et al. (2023) with an older version of OSPREI and a quiescent background.We see that the improvements within OSPREI and the inclusion of MEOW-HiSS improve our ability to reproduce this event, finding excellent matches to the observed arrival time, CME duration, and total magnetic field profile.We also see good matches to the vector components of the magnetic field, temperature, and number density, but suspect some difference occurs near the back of the CME due to reconnection with the HSS, which is not captured by OSPREI.OSPREI reproduces the observed velocity very well in the front of the CME and even shows an enhancement toward the rear, but much weaker than seen in the observed case.
Overall the coupling of OSPREI and MEOW-HiSS shows promising initial results.Further validation against observed cases is required but we expect that this could be a powerful tool for forecasting complicated, interacting cases, which often tend to drive the most significant space weather events.

Figure 1 .
Figure 1.Example Mostly Empirical Operational Wind with a High Speed Stream results.Panel (a) shows the regions of the high speed stream using alternating blue and orange sections.Panel (b) shows contours of the velocity for the same case shown in panel (a).

Figure 2 .
Figure 2. (Top) Cartoon illustrating the change in the drag calculation.The gray arcs represent an edge-on view of a coronal mass ejection (CME).Previously, the drag force was calculated at a single point, illustrated with the single dark blue arrow at the single black dot in the top left CME.Now, we calculate the drag force at the front and back of the CME, indicated by the dark blue and light blue arrows.The average of these two values determines the net drag on the CME.The difference between them determines the amount of compression or expansion from the differential drag.(Bottom) Cartoon showing the yaw rotation.We determine the drag force at both edges (dark blue arrows).The difference between these determines a torque about the CME nose (black dot), which causes a rotation we refer to as a change in the yaw.

Figure 3 .
Figure 3. Profiles for the fast coronal mass ejection (CME) results with different solar wind (SW) backgrounds.Panel (a) shows an ambient background and (b) a static high speed stream (HSS) with area 8 × 10 10 km 2 and front at 1 au.Panels (c) and (d) show a time-dependent 8 × 10 10 km 2 background, aligned to have the CME either enter (c) or exit (d) the HSS at the same distance as in panel (b).Within each panel, the top row shows a graphical representation of the extent of the sheath (maroon) and the CME (gray) versus the distance of the CME front.The HSS is shown in a similar manner with the alternating light blue and orange showing the different regions within the Mostly Empirical Operational Wind with a High Speed Stream HSS (3 stream interaction region, 1 plateau, 1 tail).The second row shows the velocity of the sheath (maroon), CME front (black), and SW at the location of the CME front (blue dashed).The third row shows the edge-on angular width of the CME (black) and the axial aspect ratio (maroon).The fourth row shows the yaw angle of the CME (black) and the rotational speed (maroon).The fifth row shows the Alfvénic Mach number of the shock (black) and the number density within the sheath (maroon).Finally, the bottom panel shows the signed logarithm of the acceleration from various forces: the total radial acceleration of the CME (black), the bulk drag (purple), the pancaking from drag (maroon), the axial magnetic forces (dashed light blue), the magnetic pressure (dashed green), and the thermal pressure (dashed yellow).The fine dashed, horizontal line indicates the zero value.

Figure 4 .
Figure 4. Same format as Figure 3 but showing profiles for the average coronal mass ejection.

Figure 5 .
Figure 5. Results for the slow coronal mass ejection with a time-dependent high speed stream.The figure has the same format as an individual subfigure in Figure 3.

Figure 6 .
Figure 6.Contours of the speed in the equatorial plane.The top row shows the slow coronal mass ejection (CME) with the time-dependent high speed stream (HSS) and the bottom without any HSS.The boundary of the CME is outlined in white.For each row, we show three time steps roughly corresponding to the time of first contact at each of the three farthest satellites (pale blue dots).
(2022) (compare to the right panels within Figures 4-9 of

Figure 7 .
Figure7.Normalized in situ profiles comparing the results for the slow coronal mass ejection (CME) with and without the background high speed stream (HSS).From top to bottom, the rows show the normalized total magnetic field, B R , B T , B N , velocity, temperature, and number density.The pink, red, maroon, and brown profiles correspond to results using an HSS background, and respectively to increasing satellite distances (25 R ⊙ , 75 R ⊙ , 150 R ⊙ , and 213 R ⊙ ).The corresponding profiles for the ambient background are respectively the light blue, blue, navy, and dark blue lines.The thicker portion of each profile corresponds to the CME and the rest to solar wind (either ambient or HSS).For each satellite distance, both the ambient and HSS profiles have been normalized using the maximum value within the CME for the HSS case (e.g., both red and blue have been normalized by the maximum in the bold part of red).For the magnetic field components, we use the maximum in the total magnetic field strength for normalization.

Figure 8 .
Figure 8. Contours of various properties based on the coronal hole area and initial front distance of the high speed stream (HSS).(a) The region in which the front of the coronal mass ejection (CME) is located at the start of a simulation.(b) The longitude at which the front of the HSS is at 1 au at the start of the simulation.(c, d) The region in which the front of the CME is located once it reaches 1 au (fast CME in (c), average CME in (d)).

Figure 9 .
Figure 9.The effects of different coronal hole areas and initial front distances on the simulation results.The contours show the change in various values from their corresponding results in an ambient solar wind (SW) background.The left column shows results for the fast coronal mass ejection (CME) and the right the average CME.The rows show the changes in the transit time, cross-sectional aspect ratio, Mach number, and yaw.The gray region in (f) indicates where the CME moving slower than the background SW and no shock/compression exists in front of the CME at the time of arrival.

Figure 10 .
Figure 10.Results for the fast coronal mass ejection (CME) with the high speed stream (HSS) tuned to maximize the yaw rotation.The left panel shows contours of the total velocity in the equatorial plane, analogous to Figure 6.The different colored circles represent synthetic satellites used by ForeCAT In situ Data Observer.The right panel shows in situ profiles for all nine satellites with the color of each line matching the satellite position in panel (a).From top to bottom in panel (b), each row shows the total magnetic field, B R , B T , and B N field components, velocity, temperature, and number density.For each profile, the solid portion represents the CME, the long dashes the sheath region, and the short dashes the solar wind (ambient or HSS).

Figure 11 .
Figure 11.Contours of the total velocity in the meridional (left) and equatorial (right) planes at several steps during the simulation.The white line indicates the edge of the simulated coronal mass ejection (CME).The top row shows early in the evolution, the middle near the time of first contact with the CME, and the bottom near the time the satellite exits the back of the CME and enters the high speed stream.

Figure 12 .
Figure 12.Comparison of the simulations of the 26 January 2022 coronal mass ejection (CME) with Parker Solar Probe in situ observations (black lines).We shows incremental steps moving from the L23 results (yellow) to the final version developed for this work (red).The sequential steps are a change in input parameters (orange), the inclusion of the CME-driven sheath (light blue), the improvement in the drag calculation (blue), the inclusion of the background high speed stream (maroon), and finally the addition of yaw rotation (red).The order of the rows and the different line styles are the same as in Figure10.The vertical dashed lines represent the time of the shock arrival, the start of the CME ejecta, and the end of the CME ejecta.

Figure 13 .
Figure 13.Ensemble results using the red line from Figure 12 as an ensemble seed.The gray lines represent individual ensemble members and the best fit member is shown in purple.The rest of the Figure is in the same format at Figure 12.

Table 2
Variation in Final CME and Sheath Properties

Table 3
OSPREI Input Parameters for the Observed Case