Magnetic field intermittency in the solar wind: PSP and SolO observations ranging from the Alfven region out to 1 AU

$PSP$ and $SolO$ data are utilized to investigate magnetic field intermittency in the solar wind (SW). Small-scale intermittency $(20-100d_{i})$ is observed to radially strengthen when methods relying on higher-order moments are considered ($SF_q$, $SDK$), but no clear trend is observed at larger scales. However, lower-order moment-based methods (e.g., PVI) are deemed more appropriate for examining the evolution of the bulk of Coherent Structures (CSs), $PVI \ge 3$. Using PVI, we observe a scale-dependent evolution in the fraction of the dataset occupied by CSs, $f_{PVI \ge 3}$. Specifically, regardless of the SW speed, a subtle increase is found in $f_{PVI\ge3}$ for $\ell =20 d_i$, in contrast to a more pronounced radial increase in CSs observed at larger scales. Intermittency is investigated in relation to plasma parameters. Though, slower SW speed intervals exhibit higher $f_{PVI \geq 6}$ and higher kurtosis maxima, no statistical differences are observed for $f_{PVI \geq 3}$. Highly Alfv\'enic intervals, display lower levels of intermittency. The anisotropy with respect to the angle between the magnetic field and SW flow, $\Theta_{VB}$ is investigated. Intermittency is weaker at $\Theta_{VB} \approx 0^{\circ}$ and is strengthened at larger angles. Considering the evolution at a constant alignment angle, a weakening of intermittency is observed with increasing advection time of the SW. Our results indicate that the strengthening of intermittency in the inner heliosphere is driven by the increase in comparatively highly intermittent perpendicular intervals sampled by the probes with increasing distance, an effect related directly to the evolution of the Parker spiral.


INTRODUCTION
Powered by the internal dynamics of the Sun, the Solar Wind (SW), a weakly collisional stream of charged particles, expands supersonically into the interplanetary medium carrying with it photospheric magnetic field lines, producing a magnetized sphere of hot plasma around the Sun, the heliosphere. (Parker 1958;Velli 1994). Magnetic field and velocity fluctuations resulting from dynamic processes (e.g., magnetic reconnection) in the solar corona or local dynamics (e.g., driven by stream interactions) that extend over a wide range in the frequency domain have been observed in the solar wind (Belcher & Davis Jr. 1971;Bruno & Carbone 2013). During the expansion, nonlinear interactions among the fluctuations result in a cascade of energy towards the smaller scales, and a character that resembles the hydrodynamic turbulence emerges (Coleman 1968;Roberts et al. 1992). Kolmogorov (1941) (hereafter K41), derived the relationship between scale dependent increment moments, or structure functions for longitudinal velocity increments, at spatial separation , δu q = |u(r + ) − u(r)| q , and the global energy dissipation rate S q = δu q ∼ q/3 ζq . (1) For a full derivation see (Rose & Sulem 1978;Frisch 1995). In Eq. 1, the scaling of field increments occurs with a unique exponent, ζ q = q/3, implying global scale invariance (self-similarity) of the fluctuations and a transfer energy rate that is independent of scale. The consequence of energy conservation in the inertial range can then also allow us to derive the relationship describing the distribution of energy among spatial scales in Fourier space as E(k) ∼ k −5/3 , corresponding to ζ 2 = 2/3. Similarly, in the solar wind, a wealth of information about turbulence can be obtained by studying the second statistical moment of the probability distribution function of the magnetic field fluctuations, or Power Spectral Density (PSD). Adopting Taylor's frozen-in hypothesis (Taylor 1938), where, V SW is the solar wind speed, the measured power spectral density in the frequency domain F (f sc ), is equivalent to the wavenumber power-spectrum E(κ) through Due to the nature of the physical processes taking place at different scales, the magnetic spectrum of the solar wind can be divided into several segments, each showing a power-law dependence over wavenumber, E(κ) ∝ κ −γ . At the largest scales, a spectral break separates the inertial from the injection scales, where the spectrum is characterized by a κ −1 dependence. In the fast solar wind, the break shifts to larger scales with increasing heliocentric distance (Bruno & Carbone 2013). At the intermediate scales, the spectrum steepens with the index, γ, taking values in the range 3/2 to 5/3 (Bavassano et al. 1982;Marsch & Tu 1990;Matthaeus & Goldstein 1982;Chen et al. 2020;Shi et al. 2021;Telloni et al. 2021). In the inertial range, fluctuations are described within the simplified nonlinear, incompressible magnetohydrodynamic (MHD) framework, and the energy that was injected into the system at the largest scales, probably of coronal origin, cascades via nonlinear, energy-conserving interactions among the oscillating modes to smaller scales (Matthaeus & Velli 2011). Once the cascade reaches ion scales, plasma dynamics is governed by kinetic processes, the spectrum steepens, and turbulent energy is converted to plasma heat through mechanisms such as ion cyclotron damping, kinetic Alfvén waves, kinetic scale current sheets, etc. (Dmitruk et al. 2004;TenBarge & Howes 2013;Karimabadi et al. 2013). The turbulent cascade is considered to be one of the main processes contributing to the non-adiabatic expansion, as well as the acceleration of the SW (see Matthaeus & Velli 2011, and references therein). Thus, for an accurate description of the dynamics of the heliospheric plasma, understanding the origin and evolution of MHD turbulence is crucial (Bruno & Carbone 2013). However, even though the analysis of conventional spectral properties can be informative, the second statistical moment of the probability distribution function of increments is only sufficient to fully characterize turbulence under the assumption of isotropic and scale-invariant fluctuations (Kolmogorov 1962). In practice, these conditions are in principle violated in space and astrophysical systems. As a matter of fact, departures from the linear scaling prediction in Eq. 1, have been observed for q greater than 3 (e.g., Burlaga (1991); Carbone et al. (2018); Chhiber et al. (2021a)), giving rise to the concept of intermittency in solar wind turbulence. These results indicate that to fully comprehend the statistics of turbulent fluctuations, the study of higher-order moments of the scale-dependent probability distribution function of increments is necessary (Frisch 1995).
To better understand intermittency, we may model the turbulent cascade as an effort of the system to approach thermal equilibrium . At shorter time scales, local turbulent relaxation may occur, giving rise to local correlations in MHD. The borders of such regions will typically not be relaxed but rather remain in a dynamic state, leading to local nonlinear interactions and processes such as magnetic reconnection or various types of instabilities. These boundaries correspond to coherent structures (CSs) , which in the case of the solar wind, can be either of coronal origin being passively advected by the SW (Borovsky 2021), or generated locally as an intrinsic feature of the ongoing nonlinear turbulent relaxation process (Matthaeus & Montgomery 1980;Veltri 1999;Greco et al. 2010;Matthaeus & Velli 2011;Matthaeus et al. 2015). Intermittency is associated with a fractally distributed population of small-scale CSs, superposed on a background of random fluctuations (Isliker et al. 2019;Chhiber et al. 2020;Sioulas et al. 2020). Even though they represent a minor fraction of the entire dataset (Osman et al. 2012;Sioulas et al. 2022b), CSs account for a disproportionate amount of magnetic energy dissipation, and have been shown to strongly influence the heating and acceleration of charged particles (Karimabadi et al. 2013;Osman et al. 2012;Tessein et al. 2013;Bandyopadhyay et al. 2020;Qudsi et al. 2020;Lemoine 2021;Sioulas et al. 2022a). The fundamental approach to studying intermittency involves the examination of the probability density functions (PDFs) of the dissipation rate. However, based on the Kolmogorov refined similarity hypothesis KRSH (Kolmogorov 1962), local averages of dissipation rate are related to the scale-dependent increments of the velocity field. Thus, intermittency is reflected on the PDF of field increments in the form of an increasing divergence with respect to a Gaussian distribution (i.e., PDFs display fatter tails) as increasingly smaller scales are involved (Castaing et al. 1990;Frisch 1995). This behavior, often referred to as multifractal, violates the concept of global scale invariance, a key assumption of the K41 theory, giving rise to the concept of local scale invariance, i.e., turbulence is characterized by a diverse set of fractals with varying scalings. At the same time, solar wind is an expanding medium. MHD fluctuations entering the super-Alfvénic wind in the trans-Alfvénic region, expected at ∼ 15 − 25R (DeForest et al. 2018), are modified in terms of structure and scale as the SW expands into the interplanetary medium driven by the turbulent cascade, as well as, shear at stream interfaces (Roberts et al. 1992) and other transients (Shi et al. 2022). It is, therefore, reasonable to expect that the statistical signatures of coherent structures evolve with heliocentric distance. Indeed, recent studies in the solar wind suggest a dynamic evolution of intermittency properties of MHD fluctuations indicating that the solar wind is an active turbulent medium involving both local and global dynamical processes that influence the higher-order statistics of fluctuations. Bruno et al. (2003) have utilized Helios data to examine the radial evolution of intermittency utilizing the flatness (i.e., SDK hereafter) of the magnetic field. Their analysis indicates a different behavior for slow and fast wind intermittency. More specifically, slow wind (V SW 500 km · s −1 ) was observed to display a higher degree of intermittency than the fast wind (V SW 600 km · s −1 ). Additionally, no radial dependence was observed for the slow wind, in contrast to an increase in intermittency with heliocentric distance for the fast solar wind. The distinct nature and radial evolution of intermittency were attributed to the different roles played by coherent non-propagating structures and by stochastic Alfvénic fluctuations for the two types of wind at different heliocentric distances. Turbulence in fast streams closer to the Sun is highly Alfvénic (i.e., the magnetic field and velocity fluctuations exhibit a high degree of correlation) and displays a self-similar (i.e., monofractal) character. However, during the expansion, due to nonlinear interaction amongst counter-propagating Alfvén waves, the fluctuations become decorrelated (Roberts et al. 1992;Chen et al. 2020;Shi et al. 2021) and the Alfvénic contribution, which tends to decrease intermittency because of its stochastic nature, is gradually depleted (Marsch & Liu 1993). On the contrary, advected structures tend to increase intermittency because of their coherent nature, while their relative contribution becomes more important with increasing heliocentric distance. As a result, the fractal nature of the magnetic field is modified, gradually approaching multifractal with increasing heliocentric distance. Slow wind does not show the same behavior since Alfvénic fluctuations have a less dominant role for this type of wind. The same line of reasoning was adopted by (Alberti et al. 2020;Telloni et al. 2021) to interpret the increasing deviation, with respect to the linear scaling expected from K41 theory, of the structure function scaling exponents (see Equation 1), as well as, Greco et al. (2012) who, utilizing the PVI method, observed an increase in the fractional volume occupied by coherent in the inner heliosphere. More recently, (Parashar et al. 2019;Cuesta et al. 2022a) have examined the relationship between SDK and R e , where R e is the Effective Reynolds Number, to show that regions with lower R e have on average lower kurtosis, at a fixed physical scale, suggestive of a less intermittent behavior. Even though Re is observed to decrease in the inner heliosphere, several effects overcome the relation of Re with intermittency, but at 1 au, a change of system dynamics begin to favor the effects from system size, resulting in progressively weaker intermittency at larger heliocentric distances, concurrent with a decreasing R e .
During its first ten encounters with the Sun, the Parker Solar Probe (P SP ) mission (Fox et al. 2016) has provided valuable measurements of solar wind particles and fields in the neighborhood of the solar wind sources. Aiming to approach the surface of the Sun by as close as 9.86 R , PSP offers unprecedented in-situ measurements in the vicinity of the Alfvén-zone, allowing us to study its influence on the evolution of spectral and intermittency properties of the field fluctuations (Kasper et al. 2021;Bandyopadhyay et al. 2022;Zhao et al. 2022). These observations will supplement simulations (Chhiber et al. 2022) and ultimately enable us to explore processes such as the heating of the solar corona and the acceleration of the solar wind in the vicinity of the Alfvén zone (Matthaeus & Velli 2011). In conjunction with the recently launched Solar Orbiter SolO (Müller et al. 2020), the synergy of the two missions offers a unique opportunity to explore the connection between the Sun and the heliosphere.
In this paper, we are interested in understanding the radial evolution of inertial range MHD turbulence and studying the basic features of scaling laws for solar wind fluctuations. We start our investigation by examining the radial evolution of intermittency without accounting for the anisotropy introduced with respect to the alignment angle, Θ V B . At a later stage, however, we show that accounting for anisotropy will complicate interpretation of the observations. For the purposes of our analysis, high-resolution magnetic field and particle data from PSP and SolO covering heliocentric distances 13 R 220 R are implemented. Our tools to study intermittency involve analytical methods such as the Partial Variance of increments (PVI), the scaling behavior of the high order moments of variations of the magnetic fields separated by a scale , or Structure Functions (SF s), and their respective scaling exponents, and finally the Scale Dependent Kurtosis (SDK) of the magnetic field.
The structure of this paper is as follows: Section 2 introduces the diagnostics of intermittency utilized in this study; Section 3 presents the selection of data (PSP, Solar Orbiter) and their processing; The results of this study are presented in Section 4: In Subsection 4.1, the radial evolution of magnetic field intermittency is investigated, and in Subsection 4.2 the dependence of intermittency on plasma parameters is examined; A summary of the results along with the conclusions is given in Section 5. When studying the radial evolution of intermittency in the solar wind, one fundamental issue arises: Due to the expansion of the solar wind, ion inertial length d i is expected to increase linearly as a function of heliocentric distance (Cuesta et al. 2022a). At the same time, PDFs of normalized magnetic field increments have been shown to display fatter tails when smaller spatial scales are considered (Bruno et al. 1999). This indicates that when the radial evolution of higher-order statistics is investigated, the use of a constant (i.e., not normalized) lag will result in calculating the intermittency diagnostics on progressively smaller spatial scales. It is thus becoming obvious that to reveal the underlying physical processes on a constant scale, it is important to normalize the spatial scales with physically relevant plasma parameters. For this reason, temporal scales are first converted to spatial scales by means of Taylor's Hypothesis (Taylor 1938 where V SW is the solar wind speed, and τ is the temporal lag. Taylor's hypothesis is founded on the idea that speeds of MHD wave modes (e.g., shear Alfvén modes, 1). Indeed, with the exception of a few sub-Alfvénic intervals during E 8 − E 11 of PSP, for the vast majority of the intervals, M A < 1 was satisfied. Note, however, that even when M A ∼ 1 the Taylor hypothesis may still be applied (Perez et al. 2021). Subsequently, the ion inertial length d i is utilized to normalize the spatial scales. The ion inertial length is defined as where c is the speed of light and ω pi is the proton plasma frequency, andñ p is the mean proton density for the respective interval. By means of integrating the conservation laws over a spherically expanding surface, a spatial scaling of n p ∝ R −2 is expected. Taking into account Equation 5, the ion inertial length is then expected to scale with distance as d i ∝ R. In Figure 1, the observed scaling of d i is presented, consistent with the theoretical expectation, as well as with previous Helios and Voyager observations (Cuesta et al. 2022a). The number of 5 hr-long intervals analyzed as a function of distance as well as the advection time of the solar wind is also presented in the same figure.

Partial Variance of Increments
The boundaries of CSs are associated with spatial variations or reversals of the local magnetic field. In recent years, a variety of methods, suitable for the detection of sharp gradients in a turbulent field, have been proposed (Bruno et al. 1999;Hada et al. 2003;Khabarova et al. 2021;Pecora, F. et al. 2021). A convenient statistical tool to perform this study, is the Partial Variance of Increments P V I (Greco et al. 2008). The PVI method has been used in the past in a variety of space plasma environments to determine the portion of the data corresponding to the underlying CSs (Tessein et al. 2013;Chasapis et al. 2015;Bandyopadhyay et al. 2020;Chhiber et al. 2020;Vasko et al. 2022;Lotekar et al. 2022). Assuming the validity of Taylor's hypothesis (Taylor 1938), the P V I index at time t, for lag = −V SW τ , is given by (Greco et al. 2008): where, |δB(t, τ )| = |B(t + τ ) − B(t)| is the magnitude of the magnetic field vector increments, and ... stands for an average over a suitably large window that is a multiple of the estimated correlation time for the magnetic field. Greco et al. (2018) have shown that as the PVI index increases, the identified events are more likely to be associated with Non-Gaussian structures that lay on the "heavy tails" observed in the PDF of scale-dependent increments, suggesting that coherent structures correspond to events of index P V I ≥ 2.5. By further increasing the threshold value θ, one can then identify the most intense magnetic field discontinuities like current sheets and reconnection sites ).

Structure Functions and Scaling Exponents
Another method for assessing intermittency in a time series is based on estimating a sequence of q th order moments of the magnetic field increments. The physical significance of this method is founded on the susceptibility of higher-order moments to concentrations of energy dissipation related to coherent structures and, from KRSH to extreme values of the magnetic field increments. We can estimate the component increments of the magnetic field at time t as where, i = X, Y, Z. Temporal lags may be recast to spatial lags by means of Taylor's approximation, and the q th order structure-function for the magnetic field can be estimated through where, |δB(t, )| = ( i δB 2 i ) 1/2 , is the magnitude of the vector magnetic field increments, and ... T stands for averaging over an interval of duration T. For increment scale , the moments are expected to display a power-law dependence, S q B ( ) ∝ ζ(q) . Thus, after the moments are calculated, power-law fits may be applied on the curves over different ranges of spatial scales to obtain the scaling exponents ζ q of the structure functions. Based upon the assumptions that the statistical properties of the turbulent fields are locally homogeneous and isotropic, i.e., the energy dissipation rate within the inertial range is constant on average, both the K41 theory (Kolmogorov 1941) in hydrodynamics, as well as, the Iroshnikov-Kraichnan (Iroshnikov 1963) model for MHD turbulence predict a linear scaling of ζ(q) with order q in the fully developed regime, where ζ q = q/3, and ζ q = q/4 respectively. However, as denoted by Landau (Oboukhov 1962;Kolmogorov 1962), irregularity of energy dissipation is expected to alter the scaling exponents of field increments with . More specifically, in the case where statistically depends on scale due to the mechanism that transfers energy from larger to smaller eddies, should be replaced by and Equation 1 should be recast to Expressing q/3 via a scaling relation with we obtain and thus where, ζ(q) = q/3+τ q /3 is generally a nonlinear function of q. Thus, when the scaling exponents ζ q of the structure functions are considered, the different local subsets of fractal dimension within a turbulent field are reflected on the departures from the linear scaling. This departure implies the violation of global scale invariance which in turn indicates a process characterized by multifractal statistics and intermittency. In recent years, it has been observed that for turbulence that is either not fully developed or is in a system of finite size, one does not have direct access to the scaling exponents ζ(q). In such cases, the extended selfsimilarity (ESS) hypothesis (Benzi et al. 1993) posits that it is still possible to investigate the functional form by plotting structure functions of different order against each other. This implies that the scaling of structure functions of order q may relate better to the behavior of another structure-function of order p than to spatial lag . As a result, ESS establishes the following scaling for the structure functions

Scale Dependent Kurtosis
An intermittency-affected generic time series exhibits alternate intervals of very high activity followed by extended periods of quiescence. Thus, intermittency in a signal is manifested in the form of a decrease in the fraction of volume occupied by structures at scale with decreasing scale. Thus, due to it's relationship to the scale dependent filling fraction F ( ) for structures through the Scale Dependent Kurtosis (SDK), defined as can be utilized to characterize the intermittency of a statistically homogeneous signal (Frisch 1995). An increase in K( ) with the involvement of smaller and smaller scales is indicative of a signal that exhibits activity over only a fraction of space, with the fraction decreasing with the scale under consideration. For a scalar that emerges from an additive random process subject to a central limit theorem (i.e., follows a Gaussian distribution), the SDK is independent of scale and attains a constant value K( ) = 3, indicating the selfsimilar character of the fluctuations. On the contrary, the PDFs of intermittent fluctuations progressively deviate from a Gaussian distribution (i.e., distributions display fat tails) at smaller scales (Frisch 1995). As a result, intermittency in a generic timeseries is manifested in the form of a monotonically increasing SDK with the involvement of smaller spatial/temporal scales. Additionally, when comparing two different time series, the one for which SDK grows more rapidly will be considered as more intermittent. Note that in the case where SDK fluctuates around a value different from 3, fluctuations are still characterized as selfsimilar but not Gaussian (i.e., formally referred to as Super-Gaussian). In this case, the fluctuations are not considered intermittent. This can be better understood when considering the way PDFs of increments are modified with scale. For instance, one may consider increments of a field φ that follow a given scaling By introducing a change of scale, → κ , where κ > 0, we get the following transformation According to this relationship, increments estimated at different scales, are characterized by the same statistical properties (Frisch 1995) This means that if κ is unique, the PDFs of the normalized increments (e.g., rescaled by their standard deviations), δφ (x) = (φ(x + ) − φ(x))/ (φ(x + ) − φ(x)) 2 1/2 collapses to a single PDF highlighting the self-similar (fractal) nature of the fluctuations. On the other hand, intermittency implies multifractality and, as a consequence, an entire range of values for κ. It is thus reasonable to expect that over the scales for which the PDF of increments collapse on to each other, the SDK will fluctuate around a constant value. analyzed. Magnetic field data were obtained at a cadence of 4.58 samples per second ( ∼ 0.218s resolution). However, it was found that the higher cadence is mostly offered close to the perihelia of PSP, while for periods where PSP is further away from the Sun, the cadence is reduced to 0.42s. For plasma moments, the cadence strongly depends on the interval studied, ranging from 0.218 − 0.874s during the encounters and to ∼ 27.9s at larger heliocentric distances for SPC, while for Spani, the median cadence is ∼ 3.5 for the encounters and ∼ 28s further away from the sun.
In the second step, magnetic field and particle data from the SolO mission from 2021-01-01 to 2021-12-01 are also employed. Magnetic field measurements from the Magnetometer (MAG) instrument (Horbury et al. 2020), downloaded from the ESA Solar Orbiter archive, have been utilized. Particle moment measurements for our study are provided by the Proton and Alpha Particle Sensor (SWA-PAS) onboard the Solar Wind Analyser (SWA) suite of instruments (Owen et al. 2020).

Data processing
In order to account for gaps in the magnetic field timeseries, the mean value of the cadence between successive measurements δτ has been estimated for each interval. Subsequently, intervals have been divided into three classes (a) δτ ≤ 250ms (b) δτ ≤ 500ms (c) δτ ≥ 500ms. For the first two classes of intervals magnetic field data have been linearly resampled to a cadence of dt = 250ms and dt = 500ms respectively, while the remaining intervals have been discarded. This decision was based on the observation that when resampling to dt = 450ms, the minimum spatial scale of 20d i could not be achieved for a minor fraction of the intervals at distances R ≤ 20R . To confirm that the different cadence does not affect the results of our study, the analysis was repeated by resampling all magnetic field data to dt = 450ms, with qualitatively similar results.
To obtain the plasma parameters from PSP, either SPAN or SPC data are utilized depending on the quality and cadence of the data for the interval. Subsequently, a Hampel filter was applied to all particle moments timeseries to eliminate spurious spikes and outliers exceeding three standard deviations from a moving average window spanning 200 points (Davies & Gather 1993). Finally, in order to maintain a sufficient statistical sample within any given interval, intervals that were found to have more than 5% of the values missing in the magnetic field or 10% in the particle timeseries have also been discarded. For the purposes of the SDK and SF's analysis, magnetic field and particle data have been divided into intervals of duration d = 5 hr. On the other hand, due to the nature of the PVI analysis (i.e., our measure of intermittency is provided in the form of a timeseries) after PVI was estimated, data have been divided into intervals of duration d = 15 min. The smaller duration intervals were chosen to mitigate the effects of mixing different types of solar wind, as well as to allow us to study intervals for which the alignment angle, Θ V B , estimated as where, || · ||, is a vector magnitude, · denotes ensemble average, and B, V SW are the magnetic field and In this section, we investigate the evolution of the magnetic field kurtosis as a function of lag and advection time of the solar wind. Following the methods described in Section 2.4, the SDK for the magnetic field magnitude was estimated for 3100 individual intervals of duration d = 5 hr. In Figure 2a the evolution of SDK as a function of advection time τ adv of the solar wind is presented. To emphasize the trend, the average of 100 intervals that fall within the same τ adv bin is estimated, and the mean value of τ adv for each bin is reflected on the line color. At the largest scales, the kurtosis is near Gaussian, and no clear trend is observed. At smaller lags, ∼ 5 · 10 3 d i , the lines intersect and beyond this point separate. In particular, for 5 · 10 3 d i , an in-crease in kurtosis is observed with increasing τ adv . For a more quantitative comparison, the evolution of the maximum value of the kurtosis (thereby K max ) as a function of τ adv is presented on the right-hand side inset of Figure  2a. The blue dots indicate K max estimated for individual intervals, and the red line is the binned mean of the same quantity. Uncertainty bars indicate the standard error of the sample, σ i / √ n, where σ i is the standard deviation and n is the number of the samples inside the bin (Gurland & Tripathi 1971). Additionally, the spatial scale at which K max is observed, i.e. (K max ), expressed in units of d i , is presented in the left inset figure. Even though there is a considerable scatter in the data, the maximum shifts towards smaller lags with increasing τ adv . Moreover, as τ adv increases, the peaks of SDK are progressively shifted to larger and larger values. As noted in the introduction, intermittency manifests itself as a growing kurtosis with decreasing spatial scale. Consequently, the increase of K max indicates a radial strengthening in intermittency within the inner heliosphere. Figure 2 offers a different perspective on the evolution of SDK as a function of τ adv for different spatial scales. The data points were binned according to , and τ adv , and the median value inside each bin was calculated, which is reflected in the colors and written in the plot. The bracketed numbers in the plots are the number of data points inside each bin. Note that bins that include less than 10 data points have been discarded. From this figure, we can understand that at small lags, there is a clear upward trend as a function of τ adv , whereas, for larger spatial scales, such trend becomes progressively less obvious.

Structure Functions (SFs).
As outlined in Section 2.3, a convenient way to describe the scaling variation of PDFs of the field fluctuations is by observing the deviation of the structurefunction scaling exponents from the linear dependence on the order. Note that for a more direct comparison of the scaling exponents with the K41 linear scaling, we adapt the ESS hypothesis (see Section 2.3). We thus proceed by dividing ζ(q) for the different lag ranges by ζ(3) for the respective range. Two d = 5 hr-long intervals, sampled by PSP, and SolO, were randomly selected, and the structure functions S q B ( ) up to sixth order were calculated. Power-law fits have been applied on each q th order structure-function in the ranges between (20−10 2 d i ), (10 2 −10 3 d i ) and (10 3 −10 4 d i ) respectively. The resulting power-law exponents ζ(q) are presented in Figure 3. The same process was then repeated for each of the 3100 intervals, and ζ(q)/ζ(3) as a function of q for the three different spatial scales is portrayed in Figure  4. It is important to note that the statistical accuracy of higher-order moments is affected by sample size. As a general rule, the highest order that can be computed reliably is q max = log(N ) − 1, where N is the number of samples (Dudok de Wit et al. 2013). In this case, since the majority of the intervals contain N ∼ 72000 samples, we can estimate q max ≈ 4, meaning that higher-order statistics should be interpreted with caution. Consequently, a shaded gray area has been added to all the figures to indicate scaling exponents recovered for moments that were determined with questionable accuracy.
In Figure 4a, the scaling exponents obtained by applying a best-fit linear gradient in log-log space over a window that spans between (20 − 10 2 d i ) are shown. Over this range of spatial scales, which, roughly speaking, corresponds to the transition region (Bowen et al. 2020), the scaling exponents obtain the highest observed values, indicating the presence of relatively stronger gradients in the magnetic field. Additionally, a roughly linear, i.e., monofractal, but super-Gaussian scaling at lower τ adv , which after normalization with ζ(3), closely resembles the K41 predicted curve, is obtained. However, as τ adv increases, the lines exhibit a more concave behavior at large q. This result is in qualitative agreement with the results obtained through the SDK method and suggests the strengthening of transition range intermittency as a function of the advection time. In the inertial range, (10 2 − 10 3 d i ) the normalized scaling exponents exhibit a concave scaling that strongly deviates from the K41 curve and does not show a dependence on the advection time of the solar wind. As a matter of fact, the concave scaling is observed at all times, indicating that the inertial range fluctuations exhibit a multifractal character even in the vicinity of the Sun. A similar behavior (i.e., multifractal scaling that does not evolve radially) is observed at yet larger, nonetheless still inertial, scales between (10 3 − 10 4 d i ).
To provide a quantitative context for the radial evolution of intermittency with respect to the advection time of the solar wind, the quantity D( , τ adv ) = ( ( z(q, ) z(q=3, ) − z k41 (q)) 2 ) 1/2 was estimated. This quantity is essentially the distance between two curves. The first curve corresponds to the scaling exponents z(q, ) estimated for 5 hr-long intervals by applying a moving power-law fit to each of the structure functions up to order q = 6, normalized by z(q = 3, ); the second curve corresponds to the K41 prediction for the scaling exponents z (q) k41 ). More specifically, the scaling exponents ζ(q, ) have been obtained by applying a moving power-law fit on each of the structure functions up to sixth order over a range of scales in the spatial domain ∈ [x i , 3x i ], where x i is the starting point of the powerlaw fit; Each scaling exponent, ζ(q, ), has been normalized by z(q = 3, ) estimated over the same range of the respective interval. As a result, for a given interval and over a certain range of spatial scales, six normalized scaling exponents ζ(q, )/z(q = 3, ) have been obtained. Subsequently, the square of the deviation from  To further quantify the divergence from the linear scaling predicted by the K41 theory of isotropic turbulence, z (q) k41 = q/3, and as a measure of intermittency at a given scale, the quantity D( , τ adv ) = ( ( z(q, ) z(q=3, ) − z k41 (q)) 2 ) 1/2 is presented. The data points were binned according to , and τ adv , and the median value inside each bin was calculated, which is reflected in the colors and written in the plot. The bracketed numbers in the plots are the number of data points inside each bin. Note that bins, including less than 10 data points, were discarded. the K41 prediction, z (q) k41 ), was estimated. Finally, the square root of the sum for q = 1, ..., 6 was calculated resulting in D( , τ adv ). By sliding the moving fit window over the spatial domain, several estimates of D( , τ adv ) have been obtained for different spatial scales over the same interval. The process was then repeated for all the available 5 hr-long intervals. Data points were subsequently binned according to and τ adv , and the median value of D( , τ adv ) inside each bin was calculated, which is reflected in the colors. The results of this analysis are illustrated in Figure 4. The picture that emerges fits well with the evolution of SDK (see Figure 2) since the general trend is towards stronger intermittency at larger τ adv for 100d i and no dependence on the advection time for 100d i . Obviously, the statistical trend is heavily influenced by the higher-order moments. The normalized scaling exponents, which relate to the moments of order q ≤ 3, remain almost linear (see Figure  4), hence contributing a negligible amount to the difference between the two curves. It is also important to note that the same process was repeated by only taking into account scaling exponents up to 4 th order, and the results were found to be qualitatively similar.

Partial Variance of Increments (PVI).
In this section, we examine the radial dependence of intermittency by considering the evolution of the fractional volume with respect to the overall fluctuations occupied by coherent structures identified by means of the PVI method. As suggested in (Greco et al. 2008), we consider a coherent structure any event for which the corresponding PVI index attains a value of P V I > 2.5. To estimate the PVI timeseries we follow the method outlined in Section 4.1.3 and employ non-overlapping 10 hr-long intervals sampled by PSP and SolO throughout the inner heliosphere. To ensure estimating the PVI on a constant plasma scale, we adapt a lag normalized by the ion inertial length estimated at the respective interval (see Section. 4.1.3). After the PVI timeseries is estimated, data are further divided into 15-minute intervals, and the mean plasma advection time τ adv for each interval is estimated. As a result of this process, a total number of ∼ 30, 000, intervals of duration d = 15 minutes are collected throughout the inner heliosphere. For a given interval, each data point is assigned a PVI value, and the fraction of the data points exceeding a  certain PVI threshold f P V I≥θ , where, θ = 1, 3, 6, is estimated. In Figure 6a,b,c, we use 25 bins, to present the average value of f P V I≥θ per bin plotted against τ adv , for PVI estimated with a lag of = 20, 500, 1000 d i , respectively. Uncertainty bars indicate the standard error of the sample. It is readily seen that as we move to smaller spatial lags (i.e., going from right to left panel), the fractional volume occupied by the extreme events is increasing. This result is consistent with the elevated probability density of extreme events when the PDFs of the magnetic increments are considered and indicates the presence of large gradients in the magnetic field at the smaller scales. In Figure 6a, for lag = 20d i and P V I greater than 3 (f P V I≥3 , orange line), no clear trend is observed with τ adv . This result introduces a paradox to our analysis as it seems to contradict the evolution of SDK at the smaller scales, shown in Figure 2a. However, when the evolution of f P V I≥6 with lag = 20d i (red line, panel a) is considered, a clear upward trend as a function of τ adv is observed. Thus, a natural hypothesis to explain the apparent contradiction is that the evolution of scale-dependent kurtosis (SDK) is dominated by the presence of the extreme events (P V I 6) that usually lay on the tails of the PDFs of normalized magnetic increments. The relationship between f P V I≥θ and K max is further examined in Appendix 6. From this analysis, we can understand that the fluctuations characterized by a smaller PVI index, albeit accounting for the bulk of the distribution of data points, have a relatively minor contribution to the final SDK values; On the other hand, even though the high PVI value events occupy only a small fraction of the dataset, they can significantly impact the behavior of SDK, due to the susceptibility of the fourth-order moment, found in the numerator of Equation 14, to extreme increments.
The trend of f P V I≥6 is also consistent with the evolution of SDK at larger scales, shown in Figure 6b,c, as in both cases, the line practically remains flat as a function of τ adv . A different trend is observed for f P V I≥3 , and = 500, 1000d i . In both cases, an abrupt increase in f P V I≥3 (yellow line) is observed up to ∼ 20 − 25 hours. This feature is of particular importance, as it could be related to the crossing of PSP through the Alfvén region and is further discussed in Section 5. Beyond this point, a slight upward trend is observed at both lags for subsequent times. Another interesting feature in Figures  6a,b,c is the evolution of fluctuations with P V I < 1, as in all three cases a monotonic increase of f P V I≤1 with τ adv is observed. Since both f P V I≤1 and f P V I≤3 are following an upward trend, we can understand that the fluctuations in the 1 P V I 3 are gradually getting depleted as the solar wind expands. In particular, for , the depletion process is gradual with a decrease of ≈ 2.5 %, observed between 5 and 130 hours. On the other hand, at the largest scales following an abrupt decrease of ≈ 5 % up to ≈ 25H, the f P V I<1 practically remains constant over the ranges examined.

ΘV B
Two important factors need to be considered when studying the radial evolution of intermittency in the solar wind: (a) MHD turbulence in the solar wind has a well-known tendency to develop and sustain several manifestations of anisotropy, e.g., wavevector anisotropy, variance anisotropy, etc. (Oughton et al. 2015). One of these is the anisotropy in magnetic field intermittency introduced by the presence of the background solar wind flow. For parallel intervals (i.e., Θ V B ≈ 0 • , or equivalently Θ V B ≈ 180 • ), the statistical signature of the magnetic field fluctuations is that of a non-Gaussian globally scale-invariant process, in contrast to multi-exponent statistics observed when the local magnetic field is perpendicular to the flow direction (Horbury et al. 2008;Osman et al. 2012). (b) Because of the conservation of magnetic flux (Parker Spiral), the radial component of the magnetic field decreases faster than the transverse component.
Consequently, as radial distance increases, so does the number of perpendicular intervals. On the other hand, as shown in the inset of Figure 7, the fraction of parallel/anti-parallel intervals is monotonically decreasing. Therefore, for a complete understanding of the radial evolution of intermittency in the solar wind, an analysis that takes into account both τ adv and Θ V B is required. In this section, we examine the radial evolution of anisotropic intermittency by means of the PVI and SDK methods. As a first step, the PVI method and the 15-min intervals described in Section 2.2 are adapted. Having confirmed that the anisotropy is symmetric with respect to Θ V B = 90 • , we proceeded by not distinguishing between parallel and anti-parallel directions. As a result, intervals with an estimated Θ init V B ≥ 90 • have been recast to Θ V B = 180 • − Θ init V B . We, therefore, require that the alignment angles lie within a range between 0 • and 90 • degrees. In Figures 8a,b,c, the fraction of PVI events at a given PVI threshold, f P V I≥θ , is plotted against the Θ V B angle for PVI estimated with lag = 20, 500, 1000 d i , respectively. For clarity, data have been binned into 45 linearly spaced bins in the Θ V B domain, and each dot indicates the mean value of f P V I≥θ within the bin. Error bars are also shown, indicating the standard error of the mean. For = 20d i (Figure  8a), the fraction of random fluctuations with P V I < 1 is rapidly decreasing as intervals with greater Θ V B angles are considered. The opposite trend is observed in the fraction of magnetic increments with P V I ≥ 1. As a matter of fact, the anisotropy grows stronger as higher PVI thresholds are considered. For instance, an increase in f P V I≥θ of at least one and two orders of magnitude is recovered between the lowest and highest Θ V B angles for PVI thresholds θ = 3, 6, respectively. Note, however, that regardless of the PVI threshold value, the increasing trend is halted at Θ V B ≈ 50 • . Beyond this point, no statistically significant differences in f P V I≥θ are observed at the largest Θ V B angles. A similar degree of anisotropy, as a function of Θ V B , is recovered for larger lags ( = 500, 1000d i ) shown in Figure 8b,c, respectively. Note, however, a deviation from the trend for f P V I≥6 indicating that the degree of anisotropy is lessened at progressively larger spatial scales.
In Figure 9, the evolution of anisotropic intermittency is examined as a function of τ adv . For this reason, the data points were binned according to Θ V B and τ adv , and  the mean value inside each bin was calculated, which is reflected in the colors and presented in the plot. The bracketed numbers in the plots are the number of data points inside each bin. Note that bins, including less than 10 data points, were discarded. In Figure 9a, the dependence of f P V I≥3 as a function of Θ V B and τ adv for lag = 20d i is illustrated. Two major but contradicting conclusions can be drawn from this figure: (1) When the evolution of intervals that belong to the same Θ V B bin is considered, a monotonic decrease in f P V I≥3 is observed for all Θ V B rows. (2) The average f P V I≥3 , with regard to Θ V B , shows a negligible, slightly positive trend as a function of τ adv . This seemingly inconsistent result can be addressed by considering the radial evolution of the Parker spiral. Closer to the Sun, the solar wind speed and background magnetic field tend to be aligned, i.e., the intervals tend to concentrate around Θ V B ∼ 0 • (180 • ). As we move further away from the Sun, this angle shifts towards the perpendicular direction, i.e., Θ V B ∼ 90 • . However, as shown in Figure 9, perpendicular intervals are typically associated with higher f P V I≥3 values. Thus, despite the gradual decrease in the fraction of coherent structures with P V I ≥ 3 as a function of τ adv for a constant Θ V B angle, on average, the fraction of the entire dataset shows signs of a very subtle increase. When considering the evolution of coherent structures of P V I, ≥ 6, Figure 9b, a slightly different evolution may be noticed. In particular, not a clear trend is observed for intervals of constant Θ V B . Additionally, as pointed out in Figure  8a the degree of anisotropy with regards to the Θ V B angle is strengthened when higher PVI thresholds are considered. Therefore, by applying the same logic as outlined before for P V I ≥ 3, we can explain the apparent increase in the fraction of coherent structures as a function of τ adv .
We move on to examine the evolution of K max as a function of τ adv and Θ V B . In order to mitigate the effects of mixing different types of solar wind, the duration of the intervals used has been reduced to d = 30 min- utes. It is important to note that even though the radial trend of kurtosis is not affected (i.e., the maximum of the kurtosis is observed to increase with increasing τ adv regardless of interval size), the curves are shifted vertically to larger values when larger averaging windows are considered. This may be attributed to the fact that by increasing the interval size, more and more extreme events are taken into consideration during the averaging process. Since these events have been shown to strongly affect the SDK (see Section 2.2), an increase of SDK for larger averaging windows is to be expected. The results of this analysis are illustrated in Figure 10. As expected, K max follows a qualitatively similar trend with f P V I≥6 . More specifically, intervals for which the magnetic field and solar wind speed tend to be aligned exhibit lower K max values. Moreover, no clear trend may be observed with τ adv when examining intervals with a similar Θ V B .
As a result of our analysis, it is apparent that understanding the physical mechanisms driving the evolution of intermittent properties in the magnetic field of the solar wind requires making a distinction between the effects of mixing strongly and less intermittently perpendicular and parallel intervals, respectively, as opposed to the evolution of turbulence during the expansion due to the local plasma dynamics.

Dependence of intermittency in plasma
parameters.

Solar Wind Speed
In this section, the relationship between solar wind speed V SW and magnetic field intermittency is investi- Figure 12. Fraction of datapoints with P V I value exceeding a given threshold, f P V I≥θ , where θ = 1, 3, 6 as a function of solar wind speed, VSW .
gated. Similarly to Section 4.1.4, to mitigate the effects of mixing different types of solar wind, the duration of the intervals used has been reduced to d = 30 minutes. Also, note that the intervals that are associated with fast solar wind V SW 600 km s −1 comprise only a minor fraction of our dataset. Moreover, the majority of these intervals were observed during the latest perihelia of PSP and thus in proximity to the Sun. Taking these arguments into account, we can understand that the study of the radial evolution of the fast wind is not feasible with our current dataset. Nevertheless, the study of intermittency properties as a function of V SW is still possible since a considerable number of intervals with V SW in the range 200 km s −1 V SW 600km s −1 have been sampled by both PSP and SolO throughout the inner heliosphere. We begin our analysis by considering the relationship between K max with τ adv and V SW , estimated for respective intervals. The results of this analysis are presented in Figure 11. In accordance with (Bruno et al. 2003;Weygand et al. 2006), we find that the kurtosis for the magnetic fluctuations in the slow solar wind exhibit higher peaks when compared to those of the fast solar wind. As a matter of fact, K max almost monotonically decreases as a function of V SW when intervals sampled at the same τ adv column are considered. Additionally, regardless of the solar wind speed, K max , increases as a function of τ adv . This result comes in disagreement with (Bruno et al. 2003), as it indicates the radial strengthening of intermittency regardless of the solar wind speed. Additionally, the relationship between the fractional volume occupied by coherent structures, f P V I≥θ , identified by using the PVI method (see Section 2.2), and V SW is examined. The results of this analysis are illustrated in Figure 12, for PVI threshold θ = 1, 3, 6. It is readily seen that fast solar wind is characterized by an elevated number density of magnetic increments with PVI index greater than unity, θ ≥ 1 (cyan line). Strictly speaking, and following the definition of (Greco et al. 2008), only events of P V I 2.5 correspond to coherent structures and consequently strengthen the intermittent character of the magnetic fluctuations. However, this result is of interest as it was recently shown (Sioulas et al. 2022b) that the number density of structures with PVI index greater than unity f P V I≥1 is very strongly correlated with the temperature of protons T p in the solar wind. At the same time, one of the clearest correlations between plasma parameters in the solar wind is the one between proton temperature with solar wind speed (Perrone et al. 2019). It is thus quite probable that events of PVI index greater than unity not only contribute to magnetic energy dissipation and the heating of the ambient plasma environment but at the same time are partly responsible for the acceleration of the solar wind. Moving on and considering the events of θ ≥ 3 (shown in yellow), a picture that contradicts our conclusions from the SDK analysis outlined earlier emerges. In particular, within the error bars, no statistically significant differences in the fractional volume occupied by coherent structures can be observed between fast and slow solar wind streams. As a matter of fact, one could even argue that a slight increase of f P V I≥3 can be observed with increasing solar wind speed. A different picture emerges when we consider the highest PVI threshold P V I ≥ 6. More specifically, f P V I≥6 is progressively reduced when faster solar wind streams are considered. This result provides a natural explanation for the lower SDK peaks observed at faster solar wind streams since, as already discussed in Section 2.2, the number density of PVI greater than 6 events, f P V I≥6 is tightly correlated with K max . We move on to investigate the evolution of f P V I≥θ as a function of V SW and τ adv . The results are presented in Figure13 for P V I ≥ 3, and P V I ≥ 6 respectively. For P V I ≥ 3, though on average slightly higher values of f f V I≥3 , may be observed at greater τ adv , there is, strictly speaking, not a clear horizontal trend. On the contrary, for P V I ≥ 6, an increasing trend is observed for most of the rows (i.e., streams of similar solar wind speed) in qualitative agreement with the increasing K max reported in Figure 12.

Normalized cross helicity.
In this section, the correlation between the normalized cross-helicity and intermittency, as indicated by the scale-dependent kurtosis (SDK) of the magnetic field magnitude, is examined. The normalized cross-helicity σ c is defined as: where z ± = v ± v a , v a = δB/ µ 0 m p n p , δB = B − B , and is the ensemble average. Note that for this analysis, the length of the interval has been reduced to d = 30min, to ensure that σ c does not vary significantly within the interval. For each interval, the median of σ c has been estimated, and intervals with standard deviation of σ c greater than 0.2 have been discarded. In Figure 14a, the dependence of SDK as a function of |σ c | is illustrated. Note that each line corresponds to the average of 100 intervals that fall within the same |σ c | bin. As shown in the right inset figure, the maximum value of kurtosis is decreasing with increasing |σ c |, indicating that Alfvénicity is negatively correlated with intermittency. However, it has been shown that the Alfvénic character of the field fluctuations in the solar wind strongly decreases with radial distance (Chen et al. 2020;Shi et al. 2021). Therefore, to distinguish between the effects of radial evolution and decrease in σ c , we show in Figure 14b, the dependence of K max as a function of |σ c | and τ adv . It is readily noticed, that on average, highly Alfvénic intervals, exhibit lower K max values. Nevertheless, it is evident that for a any |σ c | row in Figure 14, an increasing trend is observed for K max at larger τ adv values. The increase may be explained by the fact that there is still a mixing of parallel and perpendicular intervals outlined in Section 4.1.4.

SUMMARY AND CONCLUSIONS
This study has tried to address the following question: How do the statistical signatures of turbulence and intermittency evolve as the solar wind expands in the inner Heliosphere?
As already discussed in Section 1, intermittency lies at the heart of MHD turbulence in the solar wind. Thus, an improved understanding of its radial evolution can offer insights into some of the major open problems in the field of space physics, including the origins of the fluctua-tions and coherent structures observed in the solar wind; the influence of local and global dynamics in the evolution of the higher-order statistics; and ultimately into fundamental questions, such as the generation, acceleration and adiabatic expansion of solar and stellar winds. For this purpose, we have analyzed high-resolution magnetic field and particle data from the first 11 orbits of the Parker Solar Probe mission, as well as Solar Orbiter observations, ranging from the vicinity of the Alfvén region (R ≈ 13.7 R ) out to 1 au (R ≈ 215 R ). Our study has been made possible by a variety of statistical tools, such as the Scale Dependent Kurtosis, the normalized scaling exponents of the Structure functions, and the PVI method that enable us to exploit the property of PDFs of intermittency affected magnetic fluctuations to be increasingly flared out at progressively smaller scales.
The main findings of our study can be summarized as follows: (1) When methods utilizing higher-order moments are considered (e.g., SDK, SF q ), a strengthening of smallscale intermittency is observed with increasing advection time of the solar wind. Closer to the Sun, fluctuations of spatial scale ≈ 20 − 10 2 d i , exhibit a monofractallike but Super-Gaussian scaling that gradually evolves into multifractal as the solar wind expands into the interplanetary medium. Deeper in the inertial range, a multifractal scaling is observed that does not exhibit clear signs of radial evolution.
(2) The PVI method provides a different perspective on the evolution of intermittency. For lag = 20d i , the fraction of the dataset occupied by coherent structures, f P V I≥3 , displays a very subtle upward radial trend, whereas a more obvious increase is observed for f P V I≥6 . At larger spatial scales = 5 · 10 2 , 10 3 d i , the opposite trend is observed as f P V I≥6 is, within the error bars, independent of radial distance, while an increasing trend is observed for f P V I≥3 . It is important to note, however, that even though the trend remains positive at larger τ adv , the biggest gain is observed for τ adv 35 hrs.
In an effort to explain the disparity between SDK and PVI on the radial evolution of intermittency, the relationship between f P V I≥θ and the maximum values of SDK, K max was examined. We have shown that the fractional volume of events with P V I ≥ 6 is strongly correlated with K max . In light of this result, we can understand that methods relying on the estimate of higher-order moments as a measure of intermittency will mostly be affected by the extreme events that lie at the very tails of the distribution of increments. Such events are usually characterized by PVI values of the order of P V I 6 and constitute only a minor fraction 0.2% of the fluctuations observed in the solar wind. As a result, higher peaks in SDK may still be observed at larger τ adv even though f P V I≥3 stays constant as long as f P V I≥6 radially increases. However, to fully characterize the radial evolution of intermittency, one has to also take into account the evolution of coherent structures with P V I ≥ 3, as these structures have been shown to dissipate a considerable amount of magnetic energy in the solar wind (Osman et al. 2012). In other words, a comprehensive analysis of the radial evolution of intermittency in the solar wind requires the use of lower-order moment-based methods, such as PVI.
(3) CSS can both decay and reform due to local plasma dynamics during the expansion, with in situ generation being more efficient at the larger scales ( Figure 6). Nevertheless, the existence of passively advected coherent structures of Solar origin cannot be ruled out.
An observation that warrants a brief discussion is the abrupt increase incoherent structures of P V I ≥ 3 at the largest spatial scales in the vicinity of the Alfvén region. Recently Tenerani et al. (2021) have analyzed PSP, Helios, and Ullyses data to show that the evolution of the occurrence rate of Switchbacks in the solar wind is scaledependent as the fraction of longer-duration switchbacks increases with radial distance, whereas it decreases for shorter switchbacks. The PVI method is agnostic to the nature of the discontinuities, meaning that coherent structures may be identified by PVI as long as there are strong gradients in the magnetic field. As a result, several types of coherent structures such as current sheets, vortices, reconnection exhausts, and switchbacks may be identified with the PVI method. In this sense, one contributing factor to the increasing fraction of the dataset occupied by coherent structures might be the increasing trend of longer-duration switchbacks associated with increasing solar wind advection times. At the same time, several mechanisms, including stream-stream dynamic interactions, parametric decay instability of large amplitude Alfvén waves, (e.g., Biskamp & Müller 2000;Malara et al. 2001;Wan et al. 2009) might coexist simultaneously, resulting in the generation of several types of CSS.
(4) In agreement with earlier studies, we identify a strong anisotropy in intermittency with respect to the angle between the background magnetic field and the solar wind flow. Intermittency is weaker at Θ V B ≈ 0 • and is progressively strengthened at larger angles. More specifically, peaks in the SDK (K max ) are shifted upward, and an increase is observed in the fraction of the dataset occupied by coherent structures (P V I ≥ 3) when intervals with increasingly larger Θ V B angles are considered. The anisotropy is more pronounced at higher PVI thresholds but becomes weaker at progressively larger spatial scales ; (5) Even though at the smallest scales ( = 20d i ), the fraction of the dataset occupied by coherent structures (P V I ≥ 3) radially decreases for intervals with fixed Θ V B , on average (i.e., averaging over all Θ V B bins at a given τ adv column in Figure 9) the fraction of measured coherent structures increases in the inner Heliosphere. This is because closer to the Sun, the solar wind flow is statistically (anti)parallel to the magnetic field (i.e., Θ V B ≈ 0 • (180 • ). However, due to the radial evolution of the Parker Spiral, the fraction of observed parallel intervals gradually decreases. Note that the changing fraction of parallel vs. perpendicular intervals is due to the limitations of the single-point measurements by PSP and the use of the Taylor hypothesis and not a reflection of the radial evolution of turbulence. Taking the Θ V B anisotropy into account (see point (3) ), we can understand that the mixing of highly intermittent perpendicular and relatively less intermittent parallel intervals blur the averaged behavior of the radial evolution of intermittency. Due to the fact that the anisotropy is stronger at a higher P V I threshold, the averaged fraction of events with P V I ≥ 6 shows a more prominent positive radial trend with τ adv . This increase is also reflected in K max , as already discussed in (2).
(6) Solar wind of lower speed exhibits higher SDK peaks and is characterized by a higher fraction of f P V I≥6 events. However, no statistically significant differences are observed in f P V I≥3 as a function of solar wind speed. A strengthening of intermittency with respect to the advection time τ adv is observed regardless of solar wind speed. (see Figures 11, 12, 13) (7) A negative correlation is observed between the absolute value of the normalized cross-helicity and intermittency of the magnetic field. That is, Alfvénic intervals statistically display lower levels of intermittency as indicated by the maximum value of the SDK. (see Figure 14) As already discussed in point 4, the mixing of parallel and perpendicular intervals in the inner Heliosphere will result in a subtle radial increase in intermittency. However, perpendicular intervals are expected to progressively dominate with increasing heliocentric distance due to the conservation of magnetic flux (i.e., Parker Spiral). It is, therefore natural to expect that when the fraction of parallel intervals becomes statistically insignificant, the decreasing trend in the fraction of the dataset occupied by coherent structures will become apparent. Several observational studies have indicated that intermittency is expected to become progressively weaker with increasing heliocentric distance beyond 1AU (Parashar et al. 2019;Cuesta et al. 2022a). Based on our analysis, we propose that the decreasing trend in intermittency beyond 1AU can be attributed to the fact that the mixing effect is diminished due to the dominance of perpendicular intervals.
Finally, our results indicate that when it comes to analyzing the radial evolution of turbulence and intermittency, monitoring the changes in sampling direction is crucial. As the interplanetary magnetic field follows the Parker spiral, its angle with the spacecraft sampling direction will also vary as the distance from the Sun increases, which will then have an effect on measured turbulence characteristics. This effect needs to be disentangled from observations before the nature of the radial evolution of turbulence can be revealed. Previous studies using PSP data to analyze the radial evolution of intermittency have not taken this effect into consideration. However, our analysis indicates that obtaining such information is essential for understanding the more complex dynamics of the solar wind in the inner Heliosphere and can facilitate improvements to simulations of the solar wind (see also, Zhao et al. 2020;Chhiber et al. 2021b;Cuesta et al. 2022b). Our results will further the understanding of how CSS are generated and transported in the solar wind and will guide the development of future solar wind turbulence models.
We would like to thank Professor Tim Horbury and the the magnetometer instrument team of the Solar Orbiter mission. This research was funded in part by the FIELDS experiment on the Parker Solar Probe spacecraft, designed and developed under NASA contract NNN06AA01C; the NASA Parker Solar Probe Observatory Scientist grant NNX15AF34G and the HERMES DRIVE NASA Science Center grant No. 80NSSC20K0604. Solar Orbiter is a space mission of international collaboration between ESA and NASA, operated by ESA. Solar Orbiter Solar Wind Analyser (SWA) data are derived from scientific sensors which have been designed and created, and are operated under funding provided in numerous contracts from the UK Space Agency (UKSA), the UK Science and Technology Facilities Council (STFC), the Agenzia Spaziale Italiana (ASI), the Centre National d'Etudes Spatiales (CNES, France), the Centre National de la Recherche Scientifique (CNRS, France), the Czech contribution to the ESA PRODEX programme and NASA. Solar Orbiter SWA work at UCL/MSSL is currently funded under STFC grants ST/T001356/1 and ST/S000240/1. In this Appendix, we investigate the relationship between the maximum of SDK, K max , and the fraction of PVI events that exceed a threshold θ, f P V I≥θ for θ = 3, 6. For this reason, a new set of runs was initiated. The maximum value of kurtosis as well as f P V I≥3 , and f P V I≥6 , was estimated for intervals of duration d = 5hr. In Figure 15 the maximum value of Kurtosis K max is plotted against f P V I≥3 (blue line), and f P V I≥6 (red line). As already shown in Figure 6, f P V I≥θ for θ = 3, 6 obtain values that are separated by at least one order of magnitude. Consequently, for a more direct comparison, f P V I≥θ values were normalized by their corresponding maximum observed value. Though, in both cases, the underlying scatter plot does not exhibit a clear linear relationship (especially for f P V I≥3 ), a linear fit was applied to provide a rough estimate for the relationship between K max and f P V I≥θ . From this analysis, we can see that the slope for f P V I≥6 is significantly steeper than the one corresponding to f P V I≥3 . Similar ratios for the slopes were obtained when the quantities were normalized by the mean or their standard deviation. Figure 15. The maximum of SDK, Kmax, is plotted as a function of the fraction of PVI events that exceed a threshold θ, f P V I≥θ for θ = 3, 6, shown in the top and bottom pannel respectively. For a more direct comparison of the two curves f P V I≥θ was normalized to the maximum value f max P V I≥θ estimated for the respective threshold.