The unexpected role of heliospheric boundaries in facilitating interstellar dust penetration at 1–5 AU

Aims. Interstellar dust (ISD) particles penetrate the heliosphere because of the relative motion of the local interstellar cloud and the Sun. The penetrated particles pass through the heliospheric interface, that is, the region in which solar wind and interstellar plasma interact. As a result, the ISD ﬂow is modiﬁed after the passage through this region under the inﬂuence of electromagnetic force. The main goal of this work is to show how the heliospheric interface a ﬀ ects the distribution of ISD particles near the Sun. Methods. We have developed a Monte Carlo model of the ISD distribution in the heliosphere. It ﬁrst takes the e ﬀ ects of the helio-spheric interface and the rotating heliospheric current sheet into account. The e ﬀ ects of the heliospheric interface were probed using a global heliospheric model. Results. The computation results show that the heliospheric interface strongly inﬂuences the distribution of relatively small (radius a = 150 − 250nm) astronomical silicates. The unexpected ﬁnding is that the heliospheric interface facilitates the penetration of a = 150nm particles at small heliocentric distances and, particularly, to the Ulysses orbit (1 − 5AU). We demonstrate that the deﬂec-tion of ISD particles in the outer heliosheath is the principal mechanism that causes the e ﬀ ects of the heliospheric interface on the distribution near the Sun. The computations with di ﬀ erent heliospheric models show that the distribution near the Sun is sensitive to the plasma parameters in the pristine local interstellar medium. Thus, we demonstrated that being measured near the Sun, the ISD may serve as a new independent diagnostics of the local interstellar medium and the heliospheric boundaries.


Introduction
The Sun moves through the surrounding interstellar medium with a bulk velocity of ∼26.4 km s −1 (Witte 2004).In the solar reference frame, the supersonic solar wind emitted from the Sun interacts with the supersonic flow of partially ionized local interstellar medium (LISM) plasma (Baranov & Malama 1993), and as a result, a complex structure forms that consists of two shock waves (the heliospheric termination shock, TS, and interstellar bow shock, BS) and tangential discontinuity (heliopause, HP).The interaction region between the two shock waves is called the heliospheric interface.However, the existence of a bow shock for the heliosphere has been questioned so far (see, e.g.Izmodenov et al. 2005;McComas et al. 2012).Within the heliospheric interface, the plasma parameters vary greatly, especially across the discontinuity surfaces.They create a specific magnetized environment that affects the motion of charged particles.
The same relative motion of the Sun in the LISM allows interstellar dust (ISD) particles to penetrate the heliosphere.The ISD particles are solid grains whose cores presumably consist of carbon and silicate-rich molecules (see, e.g., Draine 2003;Tielens 2008;Altobelli et al. 2016).The characteristic sizes of ISD particles are in the range of a few nanometers to a few micrometers (see, e.g., Weingartner & Draine 2001;Wang et al. 2015).They acquire positive electric charge through the combined effects of photoelectric and secondary electron emission (see, e.g., Kimura & Mann 1998;Godenko & Izmodenov 2023a), and their trajectories are therefore affected by electromagnetic fields.To compute the trajectories of ISD grains in the heliospheric interface, the distributions of the plasma parameters (velocity and magnetic field) should be obtained from models of the global heliosphere (see, e.g., Izmodenov & Alexashov 2023;Opher et al. 2020;Pogorelov et al. 2021).
For the first time, the existence of ISD particles in the heliosphere was undoubtedly discovered by the Ulysses spacecraft in 1993 (Grün et al. 1993).During the 16-year mission, Ulysses caught 6719 dust particles, of which approximately 600 units were identified to be of interstellar origin.Using these data, Strub et al. (2015) obtained the fluxes and directions of ISD particles along the Ulysses trajectory.Sterken et al. (2015) attempted to explain the experimental fluxes by numerical modeling.In the model, they performed Monte Carlo computations of ISD trajectories inside the TS, taking the time-dependent effects induced by the rotating heliospheric current sheet (HCS) into account.Using the model, they reproduced some parts of the Ulysses data, but not the whole dataset in one simulation run.Sterken et al. (2015) assumed that a possible reason for the obtained inconsistencies might be the effects of the heliospheric interface on the motion of ISD particles, which they did not consider in the model.
Two models of the ISD distribution in the heliosphere with the effects of the heliospheric interface are currently known.Slavin et al. (2012) computed the ISD trajectories and number density distributions using the heliospheric model from Pogorelov et al. (2008).However, they applied the simplified model of the heliospheric magnetic field (HMF) with only two stationary configurations (fixed focusing and defocusing), but not the time-dependent configuration.Alexashov et al. (2016) also presented a model with the effects of the heliospheric interface.They applied the heliospheric model from Izmodenov & Alexashov (2015) and also mainly considered two stationary configurations of the HMF.However, in Fig. 3 of their paper (bottom right panel), the authors demonstrated the results of time-dependent modeling, but did not analyze them because they mainly focused on the deflection of the stationary ISD flow at the heliospheric boundaries.To date, we therefore lack a complete study of the theoretical distribution of the ISD in the heliosphere where the heliospheric interface and time-dependent HMF are both taken into account.
The main goal of this Letter is to study the effects of the heliospheric interface on the distribution of ISD particles near the Sun.For this purpose, we develop a Monte Carlo model of ISD distribution in the heliosphere with the effects of the heliospheric interface and time-dependent HMF.Using the model, we show the modification in the distribution of number density and direction of ISD flow inside the heliospheric interface and the effects of the modified flow (compared with the uniform undisturbed flow) on the distribution near the Sun.The structure of the paper is as follows.In Sect.2, we briefly describe the developed model.In Sect.3, we present the modeling results.Section 4 concludes the paper.

Model
We assumed that ISD particles are the astronomical silicates MgFeSiO 4 (see, e.g., Greenberg & Li 1996).Carbon dust particles almost do not penetrate close to the Sun because the radiation pressure force dominates for these particles (see, e.g., Kimura & Mann 1999;Kimura et al. 2002).For simplicity, we also assumed that probed dust particles have a spherical shape.

Coordinate system
We first introduce the coordinate system we used in this study.This coordinate system is centered on the Sun, and the Z-axis is directed from the Sun toward the undisturbed interstellar flow (upwind direction).The interstellar velocity and magnetic field vectors form the so-called BV plane.The X-axis belongs to the BV plane and is perpendicular to the Z-axis, such that the projection of the interstellar magnetic field vector onto the X-axis is negative.The Y-axis completes the right-handed coordinate system.For reference, the direction of the solar north pole in this coordinate system is (0.6696, −0.7373, 0.089).

Forces
Three main forces determine the dynamics of ISD particles in the global heliosphere: the gravitational attraction to the Sun, the radiation pressure from solar photons, and the electromagnetic force.The drag force in the heliosphere is negligible compared with the radiation pressure (Gustafson 1994).The resulting motion equation is as follows: where G is the gravitational constant, M S is the solar mass, Q is the electric charge of dust particles, m is the particle mass, u is the particle velocity vector, r is the heliodistance, and u p and B are the surrounding plasma velocity and magnetic field vectors, correspondingly.The parameter β is the dimensionless parameter equal to the ratio of the magnitude of radiation pressure force to the magnitude of gravitational force.The value of β depends on the chemical and geometrical properties of the grain and on the characteristics of the solar radiation (see, e.g., Burns et al. 1979).
The charge Q = 4πε 0 Ua (ε 0 is the vacuum permittivity, and a is the grain radius) was computed using the values of the surface charge potential U from Godenko & Izmodenov (2023a).The local values of u p and B were taken from the global heliospheric model of Izmodenov & Alexashov (2015, 2020).The description of the heliospheric model is presented in Appendix A.

Simulation method
To compute the ISD distributions, we applied a Monte Carlo approach.The computational domain consists of cubic cells with a space resolution of 1 AU × 1 AU × 1 AU and a time resolution of 1 yr.We started the motion of ISD particles from the plane perpendicular to the upwind direction (Z = Z 0 ).In this plane, we assumed that the spatial distribution of ISD particles is uniform and that all dust grains have the same velocities equal to the velocity of the undisturbed LISM.The sampling procedure for each dust grain and the details of the Monte Carlo computations are presented in Appendix B.
To study the effects of the heliospheric interface, we considered two cases of the location of the initial plane: (1) Z 0 = 1000 AU (with the heliospheric interface), and (2) Z 0 = 50 AU (without the heliospheric interface).We further denote these cases as model 1 (Z 0 = 1000 AU) and model 2 (Z 0 = 50 AU) for convenience.Figure 1 demonstrates the schematic picture in the BV plane.Model 2 is analogical to the model of Sterken et al. (2015).

Results of the modeling
We first demonstrate the effects of the heliospheric interface on the distribution of ISD particles along the trajectory of Ulysses. Figure 2 shows the counts of ISD particles per day obtained by models 1 and 2. The computation results show substantial differences for the smallest (a = 150 nm) grains.In particular, the count values for model 2 are zero for almost the entire mission, while for model 1, there are significant peaks in the middle (∼1997-1998) and close to the end (∼2002-2006) of the mission.Thus, the heliospheric interface facilitates the penetration of particles that are small enough at the heliodistances of the Ulysses orbit (1-5 AU).This somewhat contradicts the concepts that the heliospheric boundaries filter these particles out (Linde & Gombosi 2000;Slavin et al. 2012).
For the medium-sized grains (a = 250 nm), the corresponding distributions presented in Fig. 2 are similar, but differ noticeably in the middle of the mission (∼1997−1999).For the largest grains (a = 500 nm), the results are almost the same, that is, the heliospheric interface does not influence the trajectories of these particles.We note that the modeled results for medium-sized grains reproduce the Ulysses experimental data well for the second half of the mission (∼1999−2006).In the first half, however, the discrepancies between the modeled and experimental results are substantial.However, the Ulysses experimental counts we present refer to dust grains of all sizes.For a correct comparison with the Ulysses experimental data, the modeling results averaged over the ISD size distribution should therefore be used, but this is beyond the scope of the paper.
The effects of the heliospheric interface are more pronounced for smaller particles because the magnitude of the electromag- To study the differences between the results obtained by models 1 and 2, we considered the distributions of ISD particles in the plane Z = 0, which almost coincides with the plane of the Ulysses trajectory.Figure 3 presents the distributions of the number density for a = 150 nm ISD particles in the plane Z = 0 in 2004−2005 obtained by models 1 and 2. The period L4, page 3 of 9 2004−2005 was chosen because during this period, the differences between the two considered models are greatest.
The distribution in the plane Z = 0 is not homogeneous in either case.The areas of the ISD particle accumulations appear due to specific configurations of the HMF.The magnetic force changes sign at the HCS, and as a result, caustics of the ISD particle trajectories form.The ISD number density increases in the vicinity of these caustics.This effect was discussed by Mishchenko et al. (2020) for the stationary case and by Godenko & Izmodenov (2024) for the nonstationary case.These areas of the ISD accumulation appear in the plane Z = 0 for the two considered models, but the locations, shapes, and sizes of these areas differ significantly.In 2004−2005, the Ulysses trajectory crossed the region of increased number density for model 1, but does not do so for model 2. The fact that the Ulysses experimental data also demonstrate an increase of ISD counts during this period could serve as evidence of the heliospheric interface effects.
Inside 50 AU, the two considered models are identical, so that the differences in the ISD distribution may only appear because the ISD flow is modified beyond 50 AU.In the plane Z = 50 AU, outer heliosphere effects may be pronounced in (1) the inhomogeneous distribution of the ISD number density, and (2) deflection of the ISD flow (see Fig. C.1).Additional numerical experiments (see Appendix C) show that the deflection of the ISD flow is the main reason for the significant differences (between the results obtained from models 1 and 2) in the ISD distribution near the Sun that were discussed above.Our numerical calculations show that the velocity deflection mainly appears in the outer heliosheath (between the HP and the BS).In the inner heliosheath, the fast periodic variations in the HMF polarity lead to a fast change in the direction of the electromagnetic force, so that the net effect on the ISD velocity is weak (see, e.g., Fig. 7 from Godenko & Izmodenov 2024).
To explore the ISD velocity deflection in the outer heliosheath, we considered its Y component.Figure 4 demonstrates the evolution of the magnetic field and plasma and dust velocity components along the trajectory of a typical a = 150 nm ISD particle from the plane Z = 1000 AU up to the intersection with the heliopause.The Y-component of electromagnetic force in the LISM is approximately proportional to the following term: because the charge potential U does not change substantially in this region (see, e.g., Kimura & Mann 1998;Godenko & Izmodenov 2023a).The term (v d,x − v p,x ) • B z in this expression can be neglected because the B x component is greater than B z (in magnitude), and the Z-velocity components (dust and plasma) are greater than the X-velocity components (see Fig. 4).In the undisturbed LISM, the ISD velocity is approximately equal to the plasma velocity, so that F mag,y ≈ 0. In the outer heliosheath, where plasma decelerates (|v p,z | < |v d,z |), F mag,y becomes positive, and therefore v d,y increases.The B x component is the projection of the magnetic field vector on the axis perpendicular to the upwind direction (in the BV plane), that is, it depends in particular on the angle α BV between the interstellar magnetic field and the plasma velocity vectors.Therefore, the distribution near the Sun could also be sensitive to the parameter α BV .
To verify this hypothesis, we performed computations using the global heliospheric models of Izmodenov & Alexashov (2020) with different values of α BV .Figure 5 presents the computed counts per day on Ulysses for a = 200 nm particles for  three different values of α BV : 20 • , 40 • and 60 • .In some segments of the Ulysses trajectory, the differences between the corresponding results are significant.In particular, in the middle of the mission (1997−1998), the peak value of counts for α BV = 60 • is four months earlier than for α BV = 20 • and 40 • .In the second half of the mission, there are some points (≈2001 and 2006) where the quantitative differences of 20−30% between the results obtained by the heliospheric models with different values of α BV appear.The described differences are related to the different levels of deflection in the outer heliosheath.
Because the distributions of the ISD particles near the Sun are sensitive to the heliospheric model, the ISD data might therefore potentially be used as a new distant diagnostics of interstellar parameters.The upcoming missions that will study the ISD experimentally, such as Interstellar Mapping and Acceleration Probe (McComas et al. 2018) or Interstellar Probe (Brandt et al. 2023), will open opportunities for using this approach.

Conclusions
We have studied the effects of the heliospheric interface on the distribution of ISD particles near the Sun.For this purpose, we developed a Monte Carlo model of the ISD distribution based on the global heliospheric model of Izmodenov & Alexashov (2015, 2020).The conclusions are summarized below: -The heliospheric interface strongly modifies the distribution of small enough (a = 150−250 nm) astronomical silicates near the Sun but almost does not affect the distribution of larger particles (a = 500 nm).The reason is that the magnitude of the electromagnetic force, which is dominant in the heliospheric interface, is approximately inversely proportional to the squared radius of the grains.-The unexpected result is that the heliospheric interface facilitates the penetration of a = 150 nm particles into the inner heliosphere, which somewhat contradicts previous concepts that the heliospheric boundaries filter small dust particles out (Linde & Gombosi 2000;Slavin et al. 2012).The reason for this inconsistency is the fact that our ISD distribution model takes the time-dependent rippled shape of the HCS induced by the 22-year solar cycle and by the 25-day solar rotation into account, in contrast to previous modelers, who only considered the stationary focusing or defocusing configuration of the HMF.Penetrating the heliosphere, in the inner heliosheath, ISD particles move through regions with alternatively changing polarities of the HMF, so that the net effect of the electromagnetic force on the dust particle trajectories is weak.This allows dust grains to move almost freely from the heliopause up to the entrance into the inner heliosphere.Therefore, they are not filtered out at the heliospheric boundaries.-The effects of the heliospheric interface are mainly pronounced through the deflection of the ISD flow in the outer heliosheath under the influence of the interstellar magnetic field.For a = 150 nm ISD particles, this deflection reaches 15−20% of v LISM in the direction perpendicular to the BV plane.We showed that as a result of this deflection, the ISD particles approach to the Sun more closely, and in particular, cross the Ulysses trajectory.When the deflection in the outer heliosheath is not taken into account, small ISD particles almost entirely miss the Ulysses orbit.-The distribution of ISD particles near the Sun is sensitive to the interstellar plasma parameters.This fact provides a new independent approach for distant diagnostics of the LISM and heliospheric boundaries.As follows from Eq. ( 1), to compute the trajectories of dust particles in the heliosphere, the distributions of the surrounding plasma velocity u p and magnetic field B vectors should be provided from a model of the global heliosphere.We applied the global heliospheric model of Izmodenov & Alexashov (2015, 2020).The model is based on a self-consistent solution of the ideal MHD equations for plasma particles and the kinetic Boltzmann equation for neutrals.A summary of technical details and boundary conditions used for the heliospheric model can be found in Appendix A of Izmodenov & Alexashov (2020).For computation results presented in our study, the interstellar magnetic field parameters (boundary conditions) were: B = 3.75 µG, α BV = 60 • , where α BV is the angle between the interstellar plasma velocity and magnetic field vectors.The only exception is Fig. 5, where the heliospheric models with other boundary conditions were probed.

Appendix A: Heliospheric model
Figures A.1 and A.2 demonstrate the distributions of plasma parameters (velocity and magnetic field components) from the model with B = 3.75 µG, α BV = 60 • along the upwind direction and in the planes Y = 0 and X = 0, correspondingly.One can see that plasma parameters vary substantially in the global heliosphere, especially across the discontinuity surfaces.In particular, near the heliopause, the magnetic field increases in magnitude due to compression of the plasma flow in this region.This pile-up of the magnetic field determines the proposed effects of the heliospheric interface on the motion of charged particles.
It is important to note that a self-consistent solution in the model of Izmodenov & Alexashov (2015, 2020) is obtained in the frame of a "monopole approach" which assumes a single polarity of the HMF in the entire heliosphere.This approach is chosen because the terms responsible for the influence of the magnetic field in the ideal MHD equations do not depend on the orientation of the magnetic field.Therefore, the ideal MHD solution does not depend on the magnetic field polarity.Also, in the frame of the ideal MHD, the HCS is assumed as a discontinuity, where the magnetic field changes its orientation.The structure and width of the HCS are not considered in the heliospheric model.This assumption has been confirmed so far by Voyager magnetic field measurements in the solar wind and in the inner heliosheath (see, e.g., Burlaga et al. 2022) that demonstrate a quick change of the magnetic field orientation rather than a thick (comparing to the size of the heliosheath) transition region.
However, the motion of interstellar dust strongly depends on the polarity of the HMF.In this Letter, we calculate the time evolution of the HCS in the entire heliosphere using the kinematic approach by Czechowski et al. (2010).In the frame of this approach, it is assumed that at the source surface (in our study, the sphere of three solar radii), the HCS is a circle that separates the regions of opposite polarities.Observations show that the inclination of this circle to the solar equatorial plane changes at an approximately constant rate throughout the 22year solar cycle (Hoeksema 1995).From the source surface, the HCS expands with the solar wind velocity and acquires a rippled shape due to the 25-day solar rotation.To calculate the polarity of the HMF at any given point and moment of time, we trace the trajectory of the solar wind particle from this point and moment of time back to the source surface using the plasma bulk velocity from the global heliospheric model (for details, see Appendix A of Alexashov et al. 2016).Since in the ideal MHD the HMF is frozen into the solar wind plasma, polarity does not change along the solar wind trajectory.At the source surface, the polarity value is determined by the location of the point with respect to the HCS, which is the circle inclined to the solar equatorial plane and turned around the solar rotation axis.The inclination to the solar equatorial plane is obtained from the data of Wilcox Solar Observatory1 .
Figure A.3 presents the HCS shape inside the TS in the plane Y = 0 for different moments of time throughout the solar cycle.The HCS in panel A corresponds to focusing configuration of the HMF for positively charged particles because, in this case, the electromagnetic force is directed toward the solar equatorial plane at any point.On the contrary, panel C demonstrates defocusing configuration where the electromagnetic force repels positively charged particles from the solar equatorial plane.In panels B and D, the transition phases of HCS are presented, during which particles move almost along straight lines because the multiple variations of the HMF polarity cancel the effects of electromagnetic force out (see, e.g., Godenko & Izmodenov 2024).Namely, the demonstrated variations of the rotating HCS make the distribution of ISD particles in the heliosphere essentially time-dependent.After the termination shock, plasma decelerates, and as a result, the unipolar regions of the HMF start to alternatively change even more frequently.Figure A.4 demonstrates the unipolar regions separated by the HCS near the upwind direction in the inner heliosheath (between the TS and HP).When approaching the heliopause, the averaged distance between the consequent regions of the same polarity decreases, and therefore, the HCS loses its rippled structure, and, in particular, the dissipation of the HMF could appear.However, in this Letter, we studied the transport of ISD particles through the HMF, where the proposed dissipation near the heliopause is not considered.
Thus, the magnetic field used for the modeling of ISD trajectories inside the heliosphere is a product of the corresponding component from the global heliospheric model of Izmodenov & Alexashov (2015, 2020) and the kinematically computed polarity function at the given point and moment of time.
where v x,i j , v y,i j , v z,i j are the velocity components of the i-th simulated particle when it is inside the j-th computational cell.
To compute counts per day along the Ulysses trajectory (Fig. 2 and 5), we used the obtained values of number density and bulk velocity at the corresponding computational cells containing the Ulysses orbit.The value n d LIS M = 2 km −3 was applied, which provides the modeled counts per day comparable with the Ulysses experimental data, but, actually, the value of n d LIS M is uncertain so far.The analysis of Ulysses data provides a way to estimate correctly the value of n d LIS M , but it is beyond the scope of the paper.To study the influence of the ISD velocity deflection, we again consider the motion of ISD particles from the plane Z = 50 AU, but this time take the deflection of the ISD velocity vector into account using the results from model 1 (panels B-D of Fig. C.1).The ISD number density distribution in the plane Z = 50 AU is assumed to be unmodified by the interface (that is, uniform).Panel D of Fig. C.2 shows the computed number density in the plane Z = 0 in this case of modified velocity components.Panels E and F demonstrate the differences between the results from panel D and the results obtained by models 2 and 1, respectively.One can see that the distribution presented in panel D is visually similar to the distribution from model 1.Although there are some regions of nonzero differences, the results of computations show that the deflection of the ISD velocity vector reproduces the main effects of the heliospheric interface on the distribution of ISD grains near the Sun.starts to significantly deflect.For instance, the magnitude of the Y-velocity component reaches almost 20% of v LIS M .Then, after crossing the heliopause (in the inner heliosheath), velocity components start to oscillate around the values reached after the passage through the outer heliosheath due to alternatively changed unipolar regions of the HMF.However, the amplitude of oscillations is relatively small (±5% of v LIS M ), so namely the deflection of particles in the outer heliosheath mainly produces the resulted deflection of the flow in the plane Z = 50 AU.

Fig. 1 .
Fig. 1.Schematic picture in the BV plane.Models 1 and 2 correspond to the models with Z 0 = 1000 AU and Z 0 = 50 AU, respectively.The black lines match the discontinuity surfaces from the global heliospheric model.The green arrows show the direction of the plasma velocity and the magnetic field vectors in the undisturbed LISM.

Fig. 3 .
Fig. 3. Number density distributions for a = 150 nm ISD particles in the plane Z = 0 in 2004−2005.Panels A and B correspond to models 1 and 2, respectively.Panel C shows the difference between the results from panels A and B. The black line matches the Ulysses orbit, and the white segment of the orbit corresponds to 2004−2005.

Fig. 4 .
Fig. 4. Plasma (solid lines) and dust (dashed lines) velocity (top panel) and magnetic field (bottom panel) components along the trajectory of a typical a = 150 nm dust particle.The left boundary corresponds to the intersection with the heliopause (Z ≈ 112 AU).The vertical dotted black lines match the weak interstellar bow shock.

Fig. A. 1 .
Fig. A.1.Plasma parameters (velocity and magnetic field components) distributions along the upwind direction from the heliospheric model.The vertical dotted lines match the location of discontinuity surfaces along the upwind direction.Although the interstellar plasma flow is supersonic (M ≈ 2), the existence of interstellar bow shock has been questioned so far because of the influence of the interstellar magnetic field (see, e.g.Izmodenov et al. 2005;McComas et al. 2012).However, the heliospheric model used provides the weak shock wave in the LISM.
Fig. A.2.The B x component distributions in the plane Y = 0 (left panel) and in the plane X = 0 (right panel) from the heliospheric model.The white lines match the location of discontinuity surfaces (TS and HP) in the corresponding planes.

Fig
Fig. A.3.The HCS shape inside the TS in the plane Y = 0 for different moments of time throughout the solar cycle: 1986 (panel A), 1992 (panel B), 1998 (panel C), 2004 (panel D).The red and blue colored areas match the regions of negative and positive polarities, respectively.The boundary line between the red and blue areas is the HCS.

Fig
Fig. A.4.The unipolar regions of the HMF separated by the HCS in the inner heliosheath near the upwind direction.The white lines match the location of discontinuity surfaces from the heliospheric model.The domination of red colored regions in the inner heliosheath is an optical illusion.Actually, the total sizes of blue and red colored domains are almost the same.

Appendix C :
The effects of the ISD flow modified at 50 AU on the dust distribution near the SunIn this Appendix, we demonstrate how modifications of ISD distribution in the plane Z = 50 AU influence the ISD distribution near the Sun (in the plane Z = 0).Figure C.1 presents the distributions of number density and velocity components in the plane Z = 50 AU in 1995-1996 obtained from model 1 (Z 0 = 1000 AU).The 1995-1996 period in the plane Z = 50 AU was chosen to be in accordance with the 2004-2005 period in the plane Z = 0 because ISD particles require approximately 9 years to cover the distance of 50 AU along the Z-axis direction.One can see that the number density and velocity components are both significantly modified in the heliospheric interface.For model 2 (Z 0 = 50 AU), the distributions of number density and velocity components in the plane Z = 50 AU are uniform:n d /n d LIS M = 1, u/v LIS M = (0, −1, 0).We first study the influence of the modified dust number density in the plane Z = 50 AU on the distribution in the plane Z = 0.For this purpose, we consider the motion of particles from the plane Z = 50 AU but take the inhomogeneous distribution of number density into account using the results from model 1 (panel A of Fig. C.1).The velocities of all ISD particles in the plane Z = 50 AU are assumed to be unmodified by the interface (that is, equal to the velocity of undisturbed LISM).Panel A of Fig. C.2 shows the computed number density in the plane Z = 0 in this case of modified number density.Panels B and C demonstrate the differences between the results from panel A and the results obtained by models 2 and 1, respectively.One can see that the results presented in panel A almost do not differ visually from the results obtained by model 2. The reason is that the particles that crossed the plane Z = 0 inside the area shown in Fig. C.2 come from the area in Fig. C.1 where differences in density are not strong.These areas in the plane Z = 50 AU are shown by black rectangles in Figure C.1.

Fig. C. 1 .
Fig. C.1.Number density (panel A) and velocity components (panels B-D) for a = 150 nm ISD particles in the plane Z = 50 AU in 1995−1996 obtained from model 1.The black rectangles correspond to the area from where particles cross the area [−10 AU, 10 AU]× [−10 AU, 10 AU] in the plane Z = 0.

Fig. C. 3 .
Fig. C.3.Velocity components along the trajectory of a typical a = 150 nm particle from the plane Z = 1000 AU up to the intersection with the plane Z = 50 AU.The vertical dotted red lines match the intersections with the heliopause and the weak interstellar bow shock.