MHD Modeling of the Background Solar Wind in the Inner Heliosphere From 0.1 to 5.5 AU: Comparison With In Situ Observations

The accurate prediction of solar wind conditions in the interplanetary space is crucial in the context of both scientific research and technical applications. In this study, we simulate the solar wind throughout the heliosphere from 0.1 to 5.5 astronomical units (AU) with our improved heliospheric magnetohydrodynamics (MHD) model during the time period from 2007 to 2017. The model uses synoptic magnetogram maps as input to derive the inner boundary conditions based on a series of empirical relations such as the Wang‐Sheeley‐Arge (WSA) relation. To test the performance of this model, we compare the simulation results with in situ measurements from multiple spacecraft including ACE/WIND, Solar TErrestrial Relations Observatory, Ulysses, Juno, and MErcury Surface, Space ENvironment, GEochemistry, and Ranging at different latitudes and heliocentric distances. There is an overall agreement between the model results and solar wind observations at different latitudes and heliocentric distances. Statistical analysis for Year 2007 reveals that our model can predict most of the corotation interaction regions, high‐speed streams, and magnetic sector boundaries at 1 AU. In addition, the bimodal structure of the solar wind for different latitudes is well reproduced by the model which is consistent with Ulysses data. This study demonstrates the capabilities of our heliosphere model in the prediction of the large‐scale structures of the solar wind in the inner heliosphere, and the model can be used to predict the ambient solar wind at locations of planets in the solar system such as Earth and Jupiter.


Introduction
Knowledge about the solar wind plasma properties and interplanetary magnetic field (IMF) plays a crucial role in the investigation of planetary magnetosphere and space-weather forecast. In the case of the Earth's magnetosphere, recurrent geomagnetic storms are strongly correlated with the high-speed streams (HSSs) and associated corotation interaction regions (CIRs) embedded in the background solar wind (e.g., Richardson et al., 2002;Tsurutani et al., 2006). In addition, the prediction of the ambient solar wind at locations of other planets in the solar system is also very useful and important in terms of both scientific issues and spacecraft operations. As an example, the role of the solar wind in the Jupiter's magnetosphere remains poorly understood (e.g., Bagenal et al., 2014) because of the absence of a solar wind monitor upstream of the planet. While the dynamics and configuration of Jovian magnetosphere are largely controlled by internal mechanisms associated with the internal plasma source and the planet's 10-hr period fast rotation Krupp et al., 2004;Vasyliunas, 1983), there are many magnetospheric processes with evidence of solar wind influence such as the opening and closing of magnetic flux in the outer magnetosphere (e.g., Cowley et al., 2003;Kivelson & Southwood, 2005;McComas & Bagenal, 2007;Wilson et al., 2018, and references therein). Thus, solar wind information around Jupiter is necessary for the study of Jovian magnetospheric dynamics.
However, in situ observations of the solar wind are merely available for a few points where spacecraft are located, most of which are around the Earth, for example, Advanced Composition Explorer (ACE) (Stone et al., 1998) and Wind (Gloeckler et al., 1995) spacecraft. Therefore, we have to rely on models to simulate the solar wind parameters in the heliosphere, especially in regions around planets' magnetospheres.
The treatment of the inner boundary conditions is of great importance to hybrid MHD models, because they heavily influence their prediction accuracy. The empirical formulas employed at the inner boundary vary significantly for the existing hybrid MHD models (Shen et al., 2018, and references therein). In general, the model is driven by the input of photospheric synoptic magnetic maps combined with a potential field source surface (PFSS) model (Altschuler & Newkirk Jr. 1969;Schatten et al., 1969), and the WSA model (Arge et al., 2003;Wang & Sheeley Jr. 1990) is used to specify the magnetic field and velocity at the inner boundary. The WSA-type model is the most widely used and successful empirical model that predict the state of solar wind based on the assumption of the solar wind expansion rate. It computes very fast but requires empirical calibration and can only work properly within its narrow range of validity (near the ecliptic plane around 1 AU) (Pinto & Rouillard, 2017). Other solar wind parameters including the density and temperature are then prescribed based on some assumptions, which address the empirical relationships between the density/temperature and the precalculated parameters (e.g., Detman et al., 2006;Odstrcil, 2003;Shen et al., 2018).
While a number of heliosphere models have been proposed for the solar wind at 1 AU, very few are developed for the prediction of the solar wind out of the ecliptic plane and at larger heliocentric distances. In this paper, we will present the results of the solar wind throughout the heliosphere from 0.1 to 5.5 AU at different latitudes with our improved heliospheric MHD model. Compared to previous hybrid MHD models such as Enlil (Odstrcil, 2003), here we implement a new treatment of the inner boundary conditions as similar to Shen et al. (2018). Specifically, the solar wind speed at the inner boundary is derived from the empirical WSA relation, the temperature is specified based on its empirical relation with the solar wind speed, and the magnetic field along with the density is obtained from the observations in the immediate past Carrington rotations (CRs) based on the fact that the overall heliospheric magnetic field and the solar wind energy flux vary weakly over a few solar cycles (Shen et al., 2018). Furthermore, the outer boundary in our model is extended to 5.5 AU, which is much larger than around 1.0 AU in most of the existing heliospheric MHD models. The paper is organized as follows. The model basics and the treatment of the inner boundary conditions are described in section 2. In section 3, we first compare our simulation results with spacecraft data for the solar wind near the ecliptic plane, out of the ecliptic plane, around Jupiter and inside of 1 AU, respectively. Then we take the Year 2007 as the time period to further evaluate the model performance. In section 4 we describe the aspects of the model improvement in future work. Finally, we conclude this paper with discussion and a summary in section 5.

Model
Our heliospheric MHD model was first proposed by Florinski et al. (2013), where multiple populations of plasma and neutral particles, coupled via charge-exchange interactions, can be simulated simultaneously. The model was further developed by coupling the galactic cosmic ray propagation with magnetic turbulence transport and the MHD background evolution in the heliosphere by Guo and Florinski (2016). In this paper, we make an improvement to the model with a new treatment of the inner boundary conditions to simulate Space Weather 10.1029/2019SW002262 the solar wind large-scale structures in the inner heliosphere. In the following, we first introduce the basics of the model and then present the treatment of the inner boundary conditions.

Basic MHD Equations
Since our main focus is on the solar wind in the inner heliosphere, that is, inside of 5.5 AU, the effects of the neutral hydrogen atoms and pickup ions are neglected in the current model, which are thought to be important beyond 10 AU (e.g., Kim et al., 2016;Wang & Richardson, 2001). As a result, the equations in a conservative form are expressed as where U is the vector of conservative variables, F represents the corresponding flux, and Q is the source terms concerning gravity and inertia forces, as follows (in centimetre-gram-second units): Here, is density, u is velocity, B is magnetic field, g is the solar gravitational acceleration, I denotes the unit tensor, p = p th + B 2 ∕8 is the total pressure, p th is the thermal pressure of plasma, and the energy density e is defined as e = u 2 ∕2 + p th ∕( − 1) + B 2 ∕8 . The specific heat ratio is taken to be 1.46 everywhere in the simulation domain (e.g., Shiota et al., 2014;Shen et al., 2018). In addition, when solving the equations above in a reference frame rotating with the Sun, additional terms concerning Coriolis and centrifugal forces should be introduced in the momentum and energy equations. In these terms, = Ω 2 d + 2u × Ω is the inertia force characterized by the Sun's angular velocity vector Ω and a radial position vector d (orthogonal to the rotation axis). Although the Sun's observed photosphere rotates differentially, we assume that the inner boundary rotates rigidly with a period of 24.47 days for simplicity (Snodgrass & Ulrich, 1990;Wiengarten et al., 2014).
Two alternative methods are implemented to reduce the numerical error caused by the divergence of B in this model. One is the source term cleaning method (Powell et al., 1999), and the other employs a generalized Lagrange multiplier for a mixed hyperbolic-parabolic correction (Dedner et al., 2002). Here we use the generalized Lagrange multiplier method for a more conservative magnetic field than the other one. The global MHD model solves the governing equations with a conservative finite-volume method with a second-order unsplit total-variation-diminishing-like scheme combined with second-order Runge-Kutta time integration. To achieve second-order spatial accuracy, the interface values are derived from the limited piecewise linear reconstruction, and the fluxes are calculated based on a one-dimensional (1-D) Riemann problem at each cell interface, where Harten-Lax-van Leer contact (Li, 2005) and Harten-Lax-van Leer discontinuities (Miyoshi & Kusano, 2005) approximate Riemann solvers are utilized. A detailed description of the numerical algorithms can be found in Florinski et al. (2013).

Simulation Grid
Our three-dimensional (3-D) grid system is composed of a two-dimensional (2-D) geodesic unstructured grid on a sphere and a concentric nonuniform radial grid. The 2-D geodesic grid is a Voronoi tessellation of a sphere produced from a dual triangular (Delaunay) tessellation, which is generated by a recursive subdivision of an icosahedron (Florinski et al., 2013). Such grid system can naturally avoid the singularity in the spherical coordinate system and can help to increase the computational stability as well as the integration time step.
Here, we use spherical coordinates (r, , ) defined in a reference corotating with the Sun: The origin is located at the center of the Sun, r is the heliocentric distance, is the polar angle (with the north pole corresponding to = 0), and corresponds to the Carrington longitude. The simulation domain covers 0.1 AU ≤ r ≤ 5.5 AU, 0 • ≤ ≤ 180 • , and 0 • ≤ ≤ 360 • . In the simulations presented here, the total grid number of the nonuniform radial grid is 1,004, with grid sizes increasing from 0.67 R s (R s , radius of the Sun) at the inner boundary to 1.72 R s at ∼5.5 AU. For the 2-D geodesic grid, Level 7 hexagonal geodesic grids with a total of 40,962 faces are utilized, which is equivalent to a grid resolution of Δ = 1.25 • , Δ = 1.25 • , corresponding to a time resolution of 2.1 hr.

Inner Boundary Conditions
To simulate the solar wind in the interplanetary space, the inner boundary of our MHD model is set at 0.1 AU (21.5 R s ), which is assumed to be beyond the Alfvén critical surfaces (e.g., Jian et al., 2011;Odstrcil, 2003). Note that previous observations and model predictions have shown that the slow solar wind streams can sometimes be subcritical and far from their asymptotic state even at 30 R s (e.g., Chhiber et al., 2019), which is larger than the radial distance of the inner boundary in the simulation. As a result, the solar wind is super-Alfvénic everywhere in the simulation domain and the perturbations cannot travel toward the Sun. Therefore, all the physical quantities at the inner boundary should be imposed before the simulation. Since the inner boundary is located far from the solar surface, the physical processes in the corona with low plasma are not simulated, so that the time step determined by the Courant-Friedrichs-Levy condition is computationally acceptable.
Due to the lack of observational data at the inner boundary (0.1 AU), the synoptic maps of photospheric magnetic field (1.0 AU) are used as the input to calculate the inner boundary conditions based on a series of empirical relations such as the WSA relation. Specifically, the employed empirical methods to specify the eight MHD parameters (N, V x , V , V z , P, B x , B , and B z ) at the inner boundary are as follows.
First, the radial magnetic field at the inner boundary is calculated based on the PFSS model with the input of synoptic maps from Global Oscillation Network Group (GONG) project (http://gong.nso.edu/). The CR averaged synoptic maps from GONG project are provided continuously since September 2006 (ftp:// gong2.nso.edu/mnt/oQR/mqs/). In the PFSS model, a current-free corona is assumed with a source surface (typically set at 2.5 R s ) introduced, beyond which the magnetic field is supposed to be radially orientated, so that the Laplace's equation can be solved to reconstruct the coronal magnetic field. Here we only keep the polarity of the magnetic field from the PFSS model and use the observational data at 1 AU instead to specify the radial magnetic field B r (e.g., Shen et al., 2018;Wiengarten et al., 2014), which is written as where mean(B r 1AU ) is the average radial magnetic field observed at 1 AU from the OMNI database during the past three CRs, and R b = 0.1 AU is the location of the inner boundary.
Second, the solar wind radial velocity at the inner boundary is obtained from the empirical WSA relation. The WSA model is based on the empirical relationships of the solar wind speed with the magnetic flux expansion factor ( s ) near the Sun (Wang & Sheeley Jr. 1990) as well as the minimum angular separation ( b ) between an open field foot point and its nearest coronal hole boundary (Arge et al., 2003). Thus, the solar wind radial velocity V r at 5R s can be written as follows: where V s is the minimum possible solar wind speed, while V sets the maximum solar wind speed, a 1 − a 4 are four additional free parameters. Here we set V = 675 km/s, a 1 = 0.22, a 3 = 1.0, and a 4 = 1.0 for all simulations in this paper, while V s and a 2 are two free parameters varying for different CRs as similar to Shen et al. (2018). The V s has a range from 250 to 300, and a 2 has a range from 2.0 to 4.0. The two coronal magnetic field parameters s and b can be derived from the PFSS model with input of synoptic maps mentioned above. Finally, the solar wind radial velocity at the inner boundary (0.1 AU) can be calculated from the formula (4) with a subtraction of 50 km/s to account for the acceleration in the heliosphere (McGregor et al., 2011).
The solar wind speed at the inner boundary is very important for determining other parameters like plasma density and temperature. Based on the solar wind speed calculated above, we can obtain the solar wind number density using an empirical model (e.g., Le Chat et al., 2012;Shen et al., 2018): The empirical model is based on the fact that the solar wind energy flux F e = V r varies weakly over the whole solar cycle and is largely independent of the solar wind speed and latitude (Le Chat et al., 2012). Where N 0 and V 0 are number density and velocity at 1 AU obtained from the OMNI

10.1029/2019SW002262
database, G is the gravitational constant, M s is the solar mass, and R s is the solar radius. Set V 0 = 750 km/s, and N 0 can be readily obtained from the average energy flux at 1 AU during the past three CRs.
The plasma temperature T can be derived from the T-V relationship at 1 AU [Le Chat et al., 2012, and references therein], T(K) ∼ 1∕2V r 2 (km/s). The temperature at the inner boundary is then deduced by normalizing the temperature at 1 AU to the inner boundary at 0.1 AU by the power law T ∼ 1 r 2( −1) of the case = 1.46.
Finally, the plasma thermal pressure P = 2Nk B T is determined. Note that here we assume that the solar wind proton temperature is equal to the electron temperature for simplicity.
Since we model the solar wind in the corotation frame, the other two components of the solar wind velocity and magnetic field are determined by Until now, all the eight MHD parameters are determined for the inner boundary, which are assumed stationary in the corotation frame during each CR. Thus, a quasi-steady state can be reached for each CR. The initial conditions of the heliosphere are given as follows: Although the initial conditions do not affect the steady state of the simulation, an appropriate setting of the initial conditions can help to reduce the time costs for a steady state. As a result, it takes about 21 hr of runtime on a 500-core cluster to reach a quasi-steady state.

Simulation Results
In this section, we present simulation results of our heliospheric MHD model and make comparisons with multispacecraft observations including ACE/WIND, Solar TErrestrial Relations Observatory (STEREO)-A and STEREO-B, Ulysses, Juno, and MErcury Surface, Space ENvironment, GEochemistry, and Ranging (MESSENGER). In detail, the in situ measurements of the solar wind parameters at Earth are obtained from OMNI website (https://omniweb.gsfc.nasa.gov), which are mainly derived from Wind and ACE spacecraft located at the L 1 Lagrangian point. Another valuable in situ measurements at 1 AU are from the STEREO spacecraft mission. The two nearly identical spacecraft, one of which was ahead of the Earth in its orbit (STEREO-A) and the other was trailing behind (STEREO-B), were designed to provide stereoscopic measurements of the Sun and its coronal mass ejections (Driesman et al., 2008). In addition to the in-ecliptic observations, Ulysses orbiting the Sun in a highly tilted trajectory to the solar equator (Wenzel et al., 1989) is also employed in order to study the Sun at various latitudes. All of the in situ data mentioned above with a time resolution of 1 hr can be obtained from the CDAWeb (https://cdaweb.gsfc.nasa.gov). For the validation of the solar wind around Jupiter, we use the data derived from Jovian Auroral Distributions Experiment ion sensor (JADE-I) on board the Juno mission provided by Wilson et al. (2018), and the magnetic field data are from the magnetometer investigation . In addition, the magnetic field data from MESSENGER spacecraft are used to validate the model performance inside of 1 AU (Anderson et al., 2007). Figure 1 shows the availability status of each spacecraft during the period, ranging from 25 August 2006 to 1 January 2018, when the GONG synoptic maps are available. In the following sections, we first compare our simulation results with spacecraft data for CR2060, CR2056, CR2177-2178, and days of year (DOYs) 60-245 in 2007, where multispacecraft observations are employed to study the solar wind near the ecliptic plane, out of the ecliptic plane, around Jupiter and inside of 1 AU, respectively. It should be noted that here we treat the solar wind around Jupiter and inside of 1 AU as separate cases because of their larger/smaller distance from the Sun compared with other spacecraft, though Juno and MESSENGER during this time were near the ecliptic plane as well. Then we take the Year 2007 as an example to further evaluate the model performance.

Comparison With In Situ Observations During Separate CRs
To validate our global MHD model at different latitudes with different heliocentric distances, four typical cases are discussed here for the evaluation of the model quality.
As an example, Figure 2 displays the 3-D view of the global configuration of the inner heliosphere for CR2056 in heliographic inertial (HGI) coordinates: (a) the distribution of the radial velocity on the equatorial plane and (b) the large-scale heliospheric current sheet inferred from B r = 0. The HGI coordinates are Sun centered and inertially fixed with respect to an X axis along the intersection line between the ecliptic plane and the solar equatorial plane. The Z axis is along the solar rotation axis (+Z for northward) and the Y axis completes the right-handed set. As shown in Figure 2a, the initially three thin HSSs begin to expand along the Parker spirals with a larger radial extent with increasing heliocentric distance. In addition, the heliospheric current sheet displayed in Figure 2b, which is often compared to a ballerina's skirt, is the sector boundary of the IMF, through which the polarity of the magnetic field is reversed.

Solar Wind Near the Ecliptic Plane
Since the ACE/WIND and STEREO-AB spacecraft are always near the ecliptic plane, synchronized multipoint observations for the solar wind at 1 AU are available for most of the time (see Figure 1). The in situ observations for even larger heliocentric distances near the ecliptic plane can also be obtained from Ulysses during the periods of its equator crossings. Here we select the case of CR2060 to investigate the solar wind near the ecliptic plane, when all the spacecraft were near the ecliptic plane.  in Figure 4, where the solar wind profiles shown are (from up to bottom) radial velocity (V r ), number density (N), plasma temperature (T), total magnetic field (B), and radial magnetic field (B r ). Also shown at the top left corner in each panel is the correlation coefficient (cc) of the corresponding solar wind parameter. Note that for Ulysses data, there are two estimates of the proton temperature, that is, T-large and T-small. T-large is calculated from the integration of the ion distribution in three-dimensional velocity space over all energy channels and angle bins that are statistically above noise, while T-small is estimated by summing over angle the observations at a fixed energy. In general, T-large and T-small will bracket the true temperature. During CR2060, all the spacecraft are near the solar equator (deviate 22.0 • in latitude at most for Ulysses).
As seen in Figure 4, the simulated solar wind velocity agrees quite well with in situ measurements for all spacecraft, where the ccs for OMNI, Ulysses, STEREO-A, and STEREO-B are 0.64, 0.75, 0.73, and 0.69, respectively. A very good match can be found along the orbit of Ulysses as shown in Figure 4b, where the observed HSS centered around DOY 242 is captured by our MHD model in terms of both the magnitude and stream width. For the spacecraft in the ecliptic plane including ACE/WIND, STEREO-A, and STEREO-B, the first HSS (e.g., STEREO-A, centered around DOY 241) is also well reproduced, while the latter two HSSs (e.g., STEREO-A, centered around DOYs 247 and 251) are somewhat underestimated. Interestingly, only one prominent HSS is detected for Ulysses in both the spacecraft data and model results during this CR, which is because of its off-equatorial orbit. During this CR, Ulysses left the equatorial plane and moved to higher latitudes (see Figure 3b), so that the latter two HSSs observed by other spacecraft were missed. Other solar wind properties including N, T, and B also show a good agreement in the trend with the observations, especially in the compression regions associated with CIRs where high-speed solar wind overtakes slow speed solar wind. In addition, the magnitudes of these quantities are of the same order as those in the observations. Some observed small-scale fluctuations, however, are not well reproduced in the MHD model, which we attribute to the limited grid resolution of the model. As an indicator for crossings of the sector boundary, the radial magnetic field agrees quite well with the observations especially for the spacecraft near the equator. The ccs of Br for OMNI, STEREO-A, and STEREO-B are, respectively, 0.59, 0.66, and 0.66, while the cc of Ulysses is only 0.36 owing to its polar-orbiting trajectory.

Solar Wind Out of the Ecliptic Plane
As we mentioned earlier, Ulysses, as a polar-orbiting spacecraft, can observe the solar wind at high latitudes.
Here we select the case of CR2056, when Ulysses was at high latitudes. During this CR, the two STEREO spacecraft and the ACE were all near the equatorial plane, while the Ulysses was at a higher latitude ranging from −62 • to −50 • south in latitude with a distance of about 1.69-1.85 AU from the Sun (see Figure 5). Figure 6 shows the results of CR2056 with the same format as Figure 4. Through the course of CR2056, three prominent HSSs were observed by OMNI, STEREO-A, and STEREO-B, all of which are well captured by the model. Since the longitudinal separation of STEREO-A (STEREO-B) from ACE was less than 6.2 • (3.1 • ), the arrival time of the HSSs for different spacecraft is very similar. Again, the comparison results show that the solar wind velocity V r gives the best predictions with the ccs larger than 0.81 for near-equatorial spacecraft, that is, ACE, STEREO-A, and STEREO-B. A comparison with these near-equatorial spacecraft data for other simulated solar wind parameters can be also found in Figure 6. During this CR, Ulysses was at higher latitudes (from −62 • to −50 • south) with a heliocentric distance of about 1.69-1.85 AU, which is farther from the Sun than that of the Earth. As a result, the ambient solar wind observed by Ulysses, predominantly coming from the southern coronal holes, was as expected relatively fast, uniform, and rarified through the entire CR. Furthermore, the observed radial magnetic field Br displays a positive polarity for most of the time, which is the main feature of the IMF at high latitudes near the southern pole during Solar Cycle 24. As seen from Figure 6b, all the features for solar wind at high latitudes mentioned above are well reproduced by our model, except for an additional velocity dip centered around DOY 128, which was not observed by Ulysses. The magnitudes of the simulated parameters are of the same order as the spacecraft data for all quantities. The cc of each quantity is quite low because the dominating small-scale fluctuations are missed by the model.

Solar Wind Around Jupiter
One of the main objectives of Juno mission is to explore Jupiter's polar magnetosphere and intense aurora by taking advantage of Juno's close-in polar orbits . During Juno's approach to Jupiter, the JADE-I instrument designed to measure Jovian aurora and magnetospheric ions was turned on measuring the solar wind ions for ∼40 days prior to its arrival at Jupiter (Wilson et al., 2018). This study provides valuable solar wind observations around Jupiter for both researchers of magnetospheric physics and solar wind modelers. The JADE-I began to work on 15 May 2016, that is, DOY 136, and was finally switched off on 25 June 2016, that is, DOY 177 in preparation for Jupiter orbit insertion. This period is covered by CR2177-2178, when the solar wind in situ observation at 1 AU from OMNI and STEREO-B are also available (see Figure 1). Figure 7, during DOYs 136-177 in 2016, Juno was near the solar equator at a latitude of −5.6 • with a heliocentric distance of about 5.42 AU. The ACE/WIND and STEREO-A were on the opposite sides of the Sun. The Sun's rotation angle during the time of the solar wind travels from the solar surface to the observer can be approximated as follows:

As shown in
where Ω is the angular velocity of the Sun, r is the heliocentric distance of the observer, that is, the spacecraft, r s is the location of the solar source surface, and V r is the solar wind radial velocity. Here we assume a constant solar wind radial velocity. According to this formula, the transit time for the solar wind from the solar surface to Juno can be readily calculated by assuming a constant radially propagating velocity of 400 km/s. For Juno at around 5.42 AU, the transit time of the solar wind is about 23.5 days, that is, ≈ 310.8 • , which is near the period of a CR of 27.3 days. A comparison of the simulated solar wind parameters with in situ measurements for ACE/WIND ( Figure 8a) and STEREO-A (Figure 9b) is displayed with the same format as Figure 4. The corresponding CR numbers are labeled with colored horizonal bars at the top of each panel. It should be noted that the CR2177 is from DOY 131 to 158, while the CR2178 is from DOY 158 to 185, so that the time range shown in Figure 8 consists of two partial CRs. As shown in Figure 8a, the simulated solar wind velocity for OMNI agrees quite well with the in situ data except for some periods as denoted with shading in Figure 8a. The main observed HSSs are generally captured by the model with the right magnitudes and trends. Additionally, the radial magnetic field B r indicating the sector boundary crossings is well reproduced with a cc of 0.65. Other solar wind properties show a general agreement with the OMNI data as well.
Similarly, a good agreement with the in situ data for STEREO-A can be found except for the period between DOYs 143 and 154 in Figure 8b. As we can see in Figure 1, the solar activity during CR2177-2178 was much higher than that of the Year 2007, which means that the temporal variations of the Sun may be stronger and more solar wind transients were present. Therefore, the CR-averaged inner boundary employed here is not appropriate any more. We also make a comparison of our simulation results (red) with Juno observations (blue) from DOY 136 to DOY 177 in Figure 9.   From top to bottom, the radial velocity V r , number density N, temperature T, solar wind ram pressure Pram, total magnetic field B, and the total pressure P t . Note that the axes of N, T, and P ram are in logarithm scales. The blue dashed vertical lines indicate the shocks/waves observed by Juno, while the red vertical lines represent the shocks from simulation results. The gray shading denotes the domain with large discrepancies between the simulation results and the spacecraft data.
last two shocks between the model and Juno data at around DOYs 169 and 174. Since the solar wind transients such as an ICME and small-scale perturbations are not currently included in our model, we do not expect a perfect consistency between models and observations.

Solar Wind Inside of 1 AU
To validate the model performance inside of 1 AU, we make a comparison of our model results with in situ measurements obtained from the MESSENGER spacecraft during its cruise phase approaching Mercury. Limited by the data availability, here we choose a period between DOY 60 and DOY 245 in 2007, when the solar activity was quite low. During this period, the spacecraft was almost on the equatorial plane with a latitude ranging from −5 • to 4 • , and its heliocentric distance was about 0.31-0.90 AU.
Since the plasma moments data like density and temperature is not available during this period, we only present the magnetic field data comparison here. Note that no data are available between DOYs 183 and 201. As shown in Figure 10, the magnitudes of the modeled magnetic field are of the same order as the spacecraft data. The radial magnetic field, indicating the crossings of the sector boundary, is well reproduced in both the trends and the magnitudes by the model. However, there are also some discrepancies between the model results and observations in some periods such as DOYs 85-102 and DOYs 228-237, particularly for the total magnetic field shown in Figure 10b. These discrepancies may be caused by the limitation of the model such as the grid resolution or the temporal effects, which are not included in the inner boundary conditions.

Comparison With In Situ Observations in 2007
As shown in Figure Figure 11 presents the results of the solar wind parameters for simulation results (red) and OMNI data (blue) during 2007. As shown in Figure 11, all the modeled quantities show good agreement with in situ data in terms of both the magnitudes and the trends during the whole year of 2007. The ccs of the V r , N, T, B, and B r are 0.63, 0.35, 0.44, 0.34, and 0.59, respectively. Similar to the results of the three CRs discussed in section 3.1, the velocity and the radial magnetic field give better predictions than other solar wind parameters. Comparisons of the simulated solar wind parameters with the observations from STEREO-A and STEREO-B are also displayed in Figures 12 and 13, respectively, where the captions are the same as Figure 11. Similar to the results of OMNI, the model results agree quite well with the observations in terms of both magnitudes and trends, particularly for the solar wind velocity and radial magnetic field.
To give a detailed comparison on the CIRs, we inspect all the CIRs for both the MHD model and in situ OMNI data in 2007. Similar to the method used in Jian et al. (2015), here we use the solar wind speed to identify every major CIR during the year. The method is based on Owens et al. (2005) and MacNeice (2009), and a detailed description can be found in Jian et al. (2015). Here we only review the major steps of the method applied in this paper.
(1) Mark all the points whose velocity is faster than one day earlier by at least 50 km/s but ignore the isolated points.
(2) Group each bunch of points as a distinct high-speed enhancements (HSE) and mark the start time t 0 and the end time t 1 of each HSE.
(3) Find the minimum velocity V min during the period from t 0 −2 days to t 0 , and the maximum velocity V max during the time from t 0 to t 1 +1 day for each HSE. (4) Locate the last time reaching V min (t 3 ) and the first time reaching V max (t 4 ), marking them as the start and end of a CIR. (5) Combine CIRs separated by less than 0.75 days and update t 0 ∼ t 4 by Steps 3 and 4. (6) Delete the CIRs with the duration shorter than 0.5 days, with V min faster than 500 km/s or V max slower than 400 km/s or speed increase less than 100 km/s. (7) Delete CIRs that cross two CRs. As displayed in Figure 14, the solar wind velocity from OMNI data (blue) and model results (red) with markers of CIRs and stream interfaces (SIs) are shown. In panels (a) and (b), the regions colored with green curves indicate the CIRs, and the blue and red dashed vertical lines show the SIs identified from OMNI  Table 1. The negative Δt = −6.5 hr means that our model tends to predict an earlier arrival of the CIRs. A detailed comparison of the IMF polarity is shown in Figure 15, where the Br displayed in panel (a) from the simulation and in situ data are both smoothed by using a 2-day running average. Panel (b) illustrates the IMF polarity identified based on the algorithm of Jian et al. (2015). Here the positive polarity means the IMF is outward from the Sun, and the spiral angle is defined in the ecliptic plane with respect to the radial direction, where 0 • corresponds to the radial outward direction. Note that the amplitudes of red lines are scaled to be smaller than the blue ones for a better comparison. As shown in this figure, the modeled Br agrees quite well with the OMNI data in both magnitudes and trends. From the IMF polarities shown in panel (b), 43 sector boundary crossings are identified by OMNI data, where 37 of them are captured by our model. Therefore, our model successfully catches 86% of the polarity reversals and misses 14% of them. A summary of the typical statistical results for sector boundary crossings are listed in Table 1. The simulation results also show good agreement with the 10.1029/2019SW002262  Figure 16, where the heliocentric distances of Ulysses along with the latitudes are displayed in panel (a), and the panels (b)-(f) are in the same format as Figure 11. As seen from this figure, the model is in a good agreement with the Ulysses data during the whole year of 2007 and reproduces a bimodal solar wind structure, which is consistent with the observations at different latitudes. As expected, the solar wind at high latitudes is dominated by the fast, uniform, and tenuous wind, while for lower latitudes the solar wind is much more dynamic with higher density and lower speed. Furthermore, all of the plasma properties and the magnetic field of the model agree well with the observations in both the trends and magnitudes except for some regions around midlatitudes, with a correlation coefficient of cc = 0.88 for solar wind velocity, cc = 0.65 for density, cc = 0.46 for temperature, cc = 0.61 for total magnetic field, and cc = 0.76 for radial magnetic field. During the period between DOYs 110 and 160, when the spacecraft was at midlatitudes around 40-60 • , the amplitudes of the simulated solar wind fluctuations are much larger than those of the observation, which is directly caused by the large gradient of solar wind parameters at the transition region between the fast and slow solar winds at the inner boundary (e.g., Guo & Florinski, 2014). In addition to the results for 2007, we also display the distribution  Figure 17. The 13-month smoothed sunspot number is shown with gray curves scaled by a reversed axis on the right-hand side. As seen from this figure, the main trends of the ccs of the solar wind parameters are correlated with the solar activity; that is, the lower the solar activity is, the higher ccs are obtained. There are many reasons responsible for this variation. Specifically, during years with higher solar activity, more solar wind transients like ICMEs occur, which will affect the model performance since these effects are not currently included in our model. Besides, the photospheric synoptic maps and the WSA model, which we used for the derivation of the inner boundary conditions, will get worse when approaching solar maximum. This shows that the solar wind prediction during solar maximum is more complicated in our simulation cases. Therefore, a more comprehensive modeling of the solar wind transients should be considered for periods approaching solar maximum in our future simulations. For all 11 years, the best model results are obtained for the parameters V r and B r , which is consistent with the results obtained for 2007.

Aspects of the Model to Be Improved in Future Simulations
In the simulations presented here, we use the inner boundary conditions derived from a series of empirical relationships with the input of GONG synoptic maps. As demonstrated in this paper, such boundary condi- tions can reproduce the solar wind parameters in reasonable agreements with the spacecraft data for most of the time. However, we also found some occasional discrepancies between the simulations and observations, particularly for the solar wind at larger heliocentric distances like 5.4 AU (i.e., Jupiter's orbit). Besides, the model performance at higher latitudes becomes less effective in the simulation cases of this work. Thus, the model concerning the solar wind at higher latitudes needs to be improved in future studies. Since the GONG synoptic maps used here to drive the simulation are integrated from the observation at Earth over each CR, an implicit assumption is that the photospheric magnetic field does not vary significantly during this period. Thus, the employed inner boundary conditions, which are assumed stationary over the whole CR, are not appropriate for periods approaching solar maxima when temporal variations dominate. Furthermore, for the solar wind at larger distances, the assumption of the static inner boundary becomes less reasonable sometimes because the start time of a CR is in the reference of the Earth, while the source of the solar wind backtracked to the solar surface along the Parker spirals may be from the inner boundary of the previous CR (see Figure 7).  Currently, our model is more suited for past events analysis, because inner boundary conditions of the model are calculated based on the CR averaged synoptic maps. Therefore, the simulated solar wind parameters upstream of the planets such as Jupiter can be used to investigate the solar wind-magnetosphere interaction by coupling with the magnetospheric MHD models. At present, every MHD model for Jovian magnetosphere utilized the idealized and simplified solar wind conditions for lack of a solar wind monitor around Jupiter (Wang et al., 2018;Zhang et al., 2018, and references therein). Thus, the coupling of the heliospheric model presented here with our Jovian magnetosphere model (Wang et al., 2018) will be considered in our future simulations.

Summary and Discussion
In this paper, we simulated the background solar wind in the inner heliosphere ranging from 0.1 to 5.5 AU with our improved heliospheric MHD model. The model employs a geodesic unstructured grid system on a sphere, which can avoid the singularity problem occurred in normal polar spherical grids that can make the calculation time steps unacceptably small and can also help increasing the computational stability. To simulate a more realistic solar wind in the interplanetary space, a treatment of the boundary conditions with input of the synoptic maps from GONG project is implemented for the inner boundary located at 0.1 AU. The model results are compared with in situ measurements from multiple spacecraft including MESSEN-GER, ACE/WIND, STEREO-A, STEREO-B, Ulysses, and Juno, covering nearly all latitudes and various heliocentric distances.
To investigate the solar wind near the ecliptic plane, out of the ecliptic plane, at the distance around Jupiter, and at the distance inside of 1 AU, four cases, that is, CR2060, CR2056, CR2177-2178, and DOYs 60-245 in 2007 are chosen when multispacecraft data are available for synchronized observation. For the solar wind near the ecliptic plane around 1 AU (near the Earth) and 1.40 AU (near the Mars's orbit at 1.42 AU), the modeled solar wind parameters agree well with the in situ observations from all the four spacecraft: ACE/WIND, STEREO-A, STEREO-B, and Ulysses, in terms of both the magnitudes and trends. As for the solar wind at higher latitudes ranging from −62 • to −50 • with a distance of about 1.69-1.85 AU from the Sun, the model can reproduce the right magnitudes and the general trend of all the parameters, though many of the small-scale perturbations are not well captured in our model. Limited by the solar wind measurements at larger distances for outer planets like Jupiter, we only chose the period when Juno data are available before its arrival to Jupiter. According to the comparison results, our model reproduces the right orders of magnitude for all plasma quantities. The general trends of the solar wind parameters roughly follow the in situ measurements after DOY 150, though there are also large discrepancies in the trends between the simulation results and in situ data during periods of DOYs 141-155 and DOYs 160-165. The discrepancies may arise from the existence of large-scale fluctuations and the solar wind transient events like the ICME in DOY 141-144. The model can also capture shocks associated with CIRs. As for the solar wind inside of 1 AU, the simulated magnetic field including the magnitude and the radial component follows the general trend of the observations. Similar to other cases mentioned in this study, the radial magnetic field component has better prediction accuracy than that of the field magnitude.
A statistical study of the solar wind simulation was also conducted for 2007 to 2017. Because the data availability for Ulysses during 2007 is better than any other years, we selected the Year 2007 as an example to compare the model results with observations for all the spacecraft mentioned above except Juno. The results demonstrate that our model can capture many of the characteristic solar wind structures including the CIRs, HSSs, and sector boundaries at 1 AU, and an overall agreement can be found between the modeled solar wind parameters and in situ observations at different latitudes and heliocentric distances. For the 11-year simulation, the model gives the best prediction in 2007 and 2008, during when the solar activity was quite low. The main trends of the ccs of the solar wind parameters are roughly consistent with the solar activity. This shows that solar wind transient events may greatly influence the performance of the model prediction, and a more comprehensive modeling of the solar wind transients should be considered for periods approaching solar maxima in future simulations.
Another heliospheric model for the prediction of the solar wind near Venus and Mars was proposed by Shiota et al. (2014), where the numerical results showed reasonable agreement in the trends with the observations. However, only the velocity and magnetic field polarity were displayed for comparison in their work. Since all the simulated solar wind parameters near 1.4 AU near the ecliptic plane agrees quite well with in situ data from Ulysses for CR2060 as discussed in section 3.1.1, our model is likely to be capable of predicting the solar wind near Mars as well. Other studies also compared the model results with Ulysses data for the validation of the model (Jian et al., 2011;Wiengarten et al., 2014). Jian et al. (2011) inspected the two CIRs with Enlil model from the Community Coordinated Modeling Center from CR2016 to CR2018 when ACE and Ulysses were in latitudinal alignment at a distance of 1.0 and 5.4 AU, respectively. They found their model could capture the two CIRs at both distances with some time shifts, though the temperatures were underestimated for both the CIRs and ambient solar wind. Although we cannot test our model in the same period as Jian et al. (2011) limited by the availability of GONG synoptic maps (see Figure 1), our model results at around Jupiter (∼5.4 AU) show the right magnitudes of all solar wind parameters compared with Juno observations. With a focus on CIRs, Wiengarten et al. (2014) also investigated CR2060 to validate their model by comparison with STEREO-A and STEREO-B and Ulysses. The agreement for the results of STEREO-A and STEREO-B is comparable to our model presented here, while for Ulysses, there exist considerable discrepancies between the simulation results and in situ data in Wiengarten et al. (2014).
Compared with the preexisting models, our model extends the previous work on the prediction of the solar wind to higher latitudes and larger distances with good performance. It has the same high prediction accuracy near Earth. The consistency of the simulation results with in situ measurements suggests that our model can be used not only for the prediction of the steady solar wind at 1 AU, for example, the arrival time and magnitudes of CIRs, but also for the solar wind prediction at larger distances extending to the orbit of Jupiter.