{'fraction_non_alphanumeric': 0.05385379173389599, 'fraction_numerical': 0.022502172024326673, 'mean_word_length': 4.219681264576315, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 17, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 16, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Diffraction of light by plasma in the solar system

We study the propagation of electromagnetic (EM) waves in the solar system and develop a Mie theory that accounts for the refractive properties of the free electron plasma in the extended solar corona. We use a generic model for the electron number density distribution and apply the eikonal approximation to find a solution in terms of Debye potentials, which is then used to determine the EM field both within the inner solar system and at large heliocentric distances. As expected, the solution for the EM wave propagating through the solar system is characterized by a plasma-induced phase shift and related change in the light ray's direction of propagation. Our approach quantitatively accounts for these effects, providing a wave-optical treatment for diffraction in the solar plasma. As such, it may be used in practical applications involving big apertures, large interferometric baselines or otherwise widely distributed high-precision astronomical instruments.

We study the propagation of electromagnetic (EM) waves in the solar system and develop a Mie theory that accounts for the refractive properties of the free electron plasma in the extended solar corona. We use a generic model for the electron number density distribution and apply the eikonal approximation to find a solution in terms of Debye potentials, which is then used to determine the EM field both within the inner solar system as well as at large heliocentric distances. As expected, the solution for the EM wave propagating through the solar system is characterized by a plasmainduced phase shift and related change in the light ray's direction of propagation. Our approach quantitatively accounts for these effects providing a wave-optical treatment of scattering in the solar plasma. As such, it may be used in a number of practical applications involving big apertures, large interferometric baselines or otherwise widely distributed high-precision astronomical instruments.

I. INTRODUCTION
The propagation of light in a refractive medium is highly dependent on the frequency of the electromagnetic (EM) wave and the properties of the medium [1,2]. This is why astronomical observations conducted in the solar system are influenced by the solar plasma. With the advent of solar system exploration by space probes, the influence of plasma on interplanetary radio communication links was studied extensively [3][4][5] and it is well characterized [6][7][8][9]. Similar efforts were conducted to account for the solar plasma in astronomical observations conducted at optical and IR wavelengths as well as γ-rays and X-rays (see [10] and references therein).
Many astrophysical phenomena require precision observations, in which the solar plasma contributes significantly. The effects of the solar plasma are especially prominent at radio wavelengths. It is necessary to account for these effects in diverse applications such as very long baseline interferometry and spacecraft navigation. These effects are also relevant to experiments aiming to achieve very high magnification via gravitational lensing. One example is the solar gravitational lens (SGL), where bending of light by the gravitational field of the Sun is used to achieve extreme light amplification and angular resolution [11], which, in the foreseeable future, may become a means to obtain direct megapixel scale imaging and spectroscopy of Earth-like planets orbiting nearby stars [11,12]. Therefore, an appropriate description of light propagation in the refractive medium of the solar plasma is an important problem. To solve this problem, we need to develop an approach to model the optical effects of solar plasma, starting from first principles and describing its influence on astronomical measurements.
Propagation of light with wavelength λ in the vicinity of a large sphere with radius R ≫ λ is typically described using the geometric optics approximation. With this approach, we can study trajectories of individual light rays and describe the plasma-induced phase shift and related frequency change [4,5,[7][8][9], but this is not sufficient to account for the plasma's influence on light amplification and resolution. For this, a wave-theoretical approach is needed. In practical applications involving modern astronomical instruments, including, for instance, large aperture telescopes and long baseline optical interferometers that conduct high-precision observations on the refractive background of solar plasma and interplanetary medium, a wave-theoretical treatment of light propagation is critical.
Recently [13], we presented a wave-optical description of the shadow cast by a large opaque sphere. Specifically, we considered the scattering of a high frequency monochromatic EM wave by the large sphere and developed a Mie theory that accounts for the presence of this obscuration. We were able to determine the EM field in the shadow in the wave zone behind the sphere. Based on that treatment, we demonstrated that there is no EM field deep in the shadow. In the present study, we rely on the tools and methods developed in that paper to describe light propagation in the vicinity of the Sun, also taking into account the interplanetary medium. Our main current concern is the effect on the EM field due to the dispersive nature of the solar plasma as the field propagates through the solar system. This paper is organized as follows. In Section II we develop an analytical approach to account for the solar corona that we model as a free electron nonmagnetic plasma. For this, we use the most generic model for an electron plasma and develop a solution for Maxwell's equations characterizing the propagation of a monochromatic EM wave in such a refractive medium. In Section III we derive a solution of the EM field equations in terms of Debye potentials. We also investigate those potentials at various distances from the Sun. In Section IV we introduce the eikonal approximation to deal with the long-range component of the scattering potential. We account for the short-range potential as we obtain a full solution for the EM field. In Section V we complete the solution by setting up appropriate boundary conditions that characterize the EM field and describe the field outside the region dominated by solar plasma. In Section VI we extend the discussion to all regions of the solar system, characterized by distance from the Sun. Finally, in Section VII we discuss the results and their practical applications. To make the main results of the paper more accessible, we placed some technical material in the appendices. In Appendix A, we present Maxwell's equations in terms of Debye potentials. Appendix B discusses the Wentzel-Kramers-Brillouin (WKB) approximation.

II. THE EXTENDED SOLAR CORONA
To describe the propagation of an EM wave in the solar system, we first need to introduce a model for the solar plasma and the interplanetary medium that would cover the heliocentric ranges of interest.
Thermonuclear reactions occurring inside the Sun result in the emission of large amounts of energy [14]. Much of this energy is released in the form of EM radiation, mainly as visible light. However, the Sun also emits a stream of charged particles, known as the solar wind. The solar wind is bimodal: a mixture of several types of wind, fast (∼ 3 × 10 6 km/h) at high latitudes of the Sun, above 40 • , and slow (∼ 1 × 10 6 km/h) at latitudes between 22 • S and 21 • N. The solar wind is ionized: electrons and protons are separate, yielding a gaseous medium in which electrons are free, with no restoring force due to nearby atomic nuclei. This plasma extends to the outer solar system.
For an EM wave of angular frequency ω propagating through plasma, the dielectric permittivity of the plasma, in general, is defined as [2]: where e is the electron's charge, m e is its mass and n e = n e (t, r) is the electron number density. The quantity ω p is known as the electron plasma (or Langmuir) frequency. It is reasonable to assume that the solar plasma is nonmagnetic, i.e., its magnetic permeability is µ = 1. Therefore, in order to evaluate the plasma contribution to Maxwell's equations, we need to know the electron number density along the path. In general, of course, the electron plasma density shows temporal variability. Thus, we start by decomposing the electron number density n e into a static part n e (r) plus a temporal fluctuation δn e : n e (t, r) = n e (r) + δn e (t, r). ( Much of our knowledge about the solar plasma comes from the tracking of spacecraft in the inner solar system [4-7, 9, 14-16]. Distant spacecraft provide information about the extent of solar plasma as we approach interstellar space [17][18][19]. The heliopause at ∼ 130 AU is the last frontier of the heliosphere, the region of space dominated by the solar wind. At this distance, the momentum density of the solar wind is no longer sufficient to repulse the rarefied hydrogen and helium that is found in interstellar space. The region just inside the heliopause is called the heliosheath: the turbulent region where the solar wind is slowed and compressed by interstellar pressure. The inner boundary of the heliosheath, the termination shock, represents the region where the solar wind first collides with the interstellar medium. In what follows, to describe the plasma distribution throughout the entire solar system (Fig. 1), we assume that the electron number density in the solar corona and the solar wind is static, spherically symmetric 1 and may be parameterized in the following generic form: where α i and β i are empirically determined values and R ⊙ is the solar radius. Note that i ≥ 2, which is needed to replicate the 1/r 2 behavior of the solar wind at large distances from the Sun, where r ≫ R ⊙ . The value R ⋆ represents the heliocentric distance to the termination shock 2 , which we take to be at R ⋆ ≃ 100 AU, roughly corresponding to the distance to the inner boundary of the heliosphere. The symbol n 0 represents the electron number density in the interstellar medium. The presence of this term in the model is for completeness only as it does not scatter light due to its assumed uniform and homogenous behavior. Specifically, for distances r > R ⋆ , trajectories of light rays do not change direction as their phase is uniformly delayed due to the uniform background given by n 0 in (3). Therefore, without loss of generality, we can take n 0 = 0 or reinstate its nonzero value if needed. Note that n e (r) is not continuous at R ⋆ , reflecting the abrupt change in plasma density at the termination shock as reported by the Voyager 1 spacecraft [22]. As a light ray crosses the termination shock at R ⋆ and proceeds into the inner solar system, it is now affected by the second term in (3). For the range of heliocentric distances R ⊙ ≤ r ≤ R ⋆ light rays are refracted by the solar plasma with their phase being delayed and their trajectories bent. Finally, as light reaches the surface of the Sun at r = R ⊙ , it is absorbed by the Sun resulting in a geometric shadow behind the Sun [13]. Expressions (1)- (3) represent what we call the extended solar corona model introduced within the entire solar system and extending beyond the termination shock. The solar plasma modeled by (3) has a variable, negative index of refraction. As such, it has the effect of deflecting the wavefronts of light passing by the Sun.
The steady-state behavior of the solar plasma is known reasonably well. There are several plasma models found in the literature that we can utilize (see discussion in [23]). One specific example of the model (3) with particular values for α i and β i is given by [5] (see also Fig. 2): where r ≥ R ⊙ . The coefficients in this model are determined empirically by processing the tracking data using radio communication links to interplanetary spacecraft. In fact, the model given by Eq. (4) was used to process Cassini tracking data [24,25]. While other models exist [9], they are generally compatible with the Cassini model (4). Although results obtained in this paper will be valid for a generic model (3), we will use the empirical model (4) for estimation purposes. Evaluating (4) for rays of light passing near the Sun with the smallest impact parameter, b = R ⊙ , we see that the electron number density, at most, would be of the order of n e (r) 10 10 cm −3 , which implies a frequency of ν p = ω p /2π = n e e 2 /πm e 1 GHz. For optical frequencies (ν = c/λ ∼ 300 THz) and at the smallest impact parameter, we see that (1) contributes at most to the order of (ω p /ω) 2 ×10 −11 , even though for radio frequencies (ν ∼ 10 GHz) this ratio is much higher: (ω p /ω) 2 ×10 −2 . As the subject of our interest is visible or near infrared light, in the present discussion we only consider terms that are linear with respect to the plasma contribution and neglect higher order terms.
The plasma frequency ω 2 p in Eq. (1), in the case of the spherically symmetric plasma distribution model (3), in the range of heliocentric distances, R ⊙ ≤ r ≤ R ⋆ , has the form Distance from the Sun (r/R ⊙ )

FIG. 2:
The electron number density model (4) (thick blue line) given by [4,5,15]. The leftmost part of the curve, at short heliocentric distances dominated by terms with higher powers of (R⊙/r) corresponds to the visible solar corona [16]. The thin dotted line shows the contribution of the inverse square term, which dominates beyond a few solar radii. The lightly shaded region on the left represents the solar interior. The approximate location of the termination shock is also marked.
This generic spherically symmetric model for the plasma frequency in the extended solar corona allows us to study the influence of solar plasma on the propagation of EM waves throughout the solar system in the range of heliocentric distances given by R ⊙ ≤ r ≤ R ⋆ .

III. THE EM FIELD AND DEBYE POTENTIALS
We begin our derivation of the EM field equations in the presence of plasma with presenting Maxwell's source-free field equations in their well-known form: These equations describe the EM field in a medium with a variable index of refraction [26], which is the case of the solar plasma with a variable electron number density. Equations (6)-(7) capture the contribution of the solar plasma to the propagation of light in the vicinity of the Sun. Following closely the derivation presented in [11], we now consider a solution to these equations. Assuming, as usual [26], the time dependence exp(−iωt) and taking k = ω/c, the time-independent parts of the electric and magnetic vectors must satisfy Maxwell's equations (6)-(7) in their time-independent form: In the case of a static, spherically symmetric plasma distribution, solving (8)-(9) is most straightforward. Following the derivation outlined in [26] the essence of which is captured in Appendix A, we obtain the complete solution of the field equations (8)- (9) in terms of the electric and magnetic Debye potentials [26], e Π and m Π: where the potentials e Π and m Π satisfy the following wave equations: We may simplify equations (16) and show that, in the case of the weakly interacting, spherically symmetric free electron plasma of the extended solar corona, they may be written in an identical form. First of all, using (1) together with (3) for ǫ, while setting µ = 1, we can rewrite the left equation in (16) as the equation that describes scattering in the presence of the plasma: A similar equation (but without the last, fourth term inside the curly braces, as µ = 1) may be obtained for m Π from the second equation (16). However, we now show that the last term inside the curly braces in (17) may also be omitted. For this, we note that ω = kc and also observe that ω 2 p is expressed in terms of various inverse powers of r, as shown in (5). We therefore introduce the static, spherically symmetric plasma potential V p (r), that has several terms that decay either as r −2 or faster: The two terms in the curly braces of (18) represent the repulsive potentials due to plasma that, based on the model given by Eq. (5), vanish as r −2 or faster. The second plasma term in this expression is dominated by a factor of (kR ⊙ ) −2 , which, given the large value of the solar radius, makes its contribution negligible, especially at optical wavelengths (λ ∼ 1 µm), for which (kR ⊙ ) −2 ∼ 5.23 × 10 −32 . Therefore, the term ∝ (ω 2 p ) ′′ in (18) may be neglected. 3 Although the remaining terms are also small, they may contribute to the phase shifts of the scattered wave and, therefore, they may affect the diffraction of light by the Sun. Thus, they will be considered. Therefore, the plasma potential, V p (r), in (17) has the following from As a result, and taking into account that the magnetic permeability µ is constant, both equations (16) take an identical form: where the quantity Π represents either the electric Debye potential, e Π/ √ ǫ, or its magnetic counterpart, m Π/ √ µ, namely Π(r) = e Π/ √ ǫ; m Π/ √ µ and plasma potential V p (r) is given by (19).
Equation (20), together with the potential given in Eq. (19), may now be used to determine the solution for the Debye potentials. Together with (10)- (15), these Debye potentials determine all the components of the EM field.
Note that (20) resembles the time-independent Schrödinger equation describing scattering problems in quantum mechanics [27][28][29]. Interestingly, various forms of the power law potential (19) appear in many problems of modern atomic physics, related to the scattering of light on a cloud of cold atoms [20,21]. Our method is developed from a generic case of spherically symmetric potentials, many of which are found in the literature describing atomic collisions [30][31][32][33]. The tools developed in our present paper may also be applicable to these problems in atomic physics.

IV. SOLUTION FOR THE EM FIELD
As we discussed above, to find the solution to the Maxwell equations (6)-(7), we first have to solve (20) for the Debye potential and then use the result in (10)- (15) to obtain each component of the EM field. Typically [26], in spherical polar coordinates, the solution to Eq. (20) is obtained by separating variables: with coefficients that are determined by boundary conditions. Direct substitution into (17) reveals that the functions R, Θ and Φ must satisfy the following ordinary differential equations: The solution to (24) is given as usual [26]: where m is an integer, and a m and b m are integration constants. Equation (23) is well known for spherical harmonics. Single-valued solutions to this equation exist when α = l(l + 1) with (l > |m|, integer). With this condition, the solution to (23) becomes Θ lm (θ) = P We now focus on the equation for the radial function (22), where, because of (23), we have α = ℓ(ℓ + 1). As a result, (22) takes the form To determine the solution to (27), we first separate the terms in the plasma potential (19) by isolating the 1/r 2 term from the rest of the terms in the plasma potential (calling it the short-range potential V sr ) and present (19) as with where µ 2 is 4 the strength of the 1/r 2 term in the plasma model at r = R ⊙ . Using the values from the phenomenological model (4), we can evaluate this term: µ 2 ≃ 5.89 × 10 15 . Also, we note that the range of V sr is very short. In fact, as we see from Fig. 2, this potential provides a negligible contribution after r ≃ 8R ⊙ . Nevertheless, as it propagates through the solar system, light acquires the largest phase shift as it travels through the range of validity of this potential. Thus, it is important to keep V sr in the model. The separation of the terms in the plasma potential (28)-(29) allows us to appropriately present the radial equation (27) as where the new index L is determined from Requiring that when µ → 0, the new index L must behave as L → ℓ, we find the solution to (31): When µ/ℓ ≪ 1, this solution behaves as For a typical region where the plasma potential (3) is present, the value of ℓ may be estimated using its relation to the classical impact parameter, namely ℓ = kb ≥ kR ⊙ = 4.37 × 10 15 . Therefore, the quantity µ/ℓ ≤ 1.75 × 10 −8 is indeed small, justifying the approximation (33).

A. Eikonal solution for Debye potential
Equations similar to (20) are typical for many problems of nuclear scattering. However, no exact solution is known for an arbitrary potential V p . This motivated the development of various approximation tools [20,21]. One such approximation is known for the case of short-range potentials that decay faster than 1/r 2 . In this case, one may introduce a small parameter and present the total solution to (20) as a series expansion with respect to this parameter. This method is called the Born approximation. The method uses the radial Green's function solution to (27) (obtained while setting V p = 0) to determine each successive term in the expansion. The final solution determines the cumulative phase shift for the EM wave as it traverses the area where the short-range scattering potential is present. Since this is a Born-type approximation for the phase shifts relative to the plasma-free wave, the relevant approximation is referred to as the distorted-wave Born approximation. It determines the additional phase shift due to the short-range potential V sr [34,35]. However, it is known that for potentials that behave as 1/r 2 , this approximation leads to divergent results, as such potentials do not decay fast enough with distance. This is precisely our situation, where the plasma potential contains the 1/r 2 terms. Thus, these methods are not very useful for our purposes.
In our approach, we address the scattering of high frequency EM waves on the plasma-induced potential V p that is bounded by the heliocentric distance to the heliopause, R ⋆ from (3). As this distance is finite, the divergence is avoided. In this case and for the case of high energy scattering, we implement the so-called eikonal (or high-energy) approximation [36][37][38][39][40][41]. In this approximation, the short-range plasma potential contributes a phase shift to the EM wave which can be directly calculated.

Solution with short-range potential Vsr absent
Eq. (20) can be solved numerically, but with a great deal of effort, especially at large energies. An exact closed form solution for Eq. (20) for the general case V sr = 0 does not exist. However, a number of approximation methods to solve equations of this type were developed for scattering problems in quantum mechanics. At large incident energies, for a wavefront moving in the forward direction, a very useful approximation becomes available. This is the eikonal approximation [36][37][38][39][40][41]. The eikonal approximation is valid when the following two criteria are satisfied [41]: kb ≫ 1 and V sr (r)/k 2 ≪ 1. In our case, both of these conditions are fully satisfied, indeed, the first condition yields kb = 4.37 × 10 15 (λ/1 µm)(b/R ⊙ ) ≫ 1 and also, taking the short-range plasma potential V sr from (29), we evaluate the second condition as V sr (r)/k 2 ≤ V sr (R ⊙ )/k 2 ≈ 4.07 × 10 −13 (λ/1 µm) 2 ≪ 1.
To develop a solution to (20) using the eikonal approximation, we first note that when the short-range potential V sr is absent, (30) takes the form The solution to this equation is well known and is given in terms of Riccati-Bessel functions [13,26]: where the subscript (2) simply stands for the solution to (34) that includes the inverse-square term, 1/r 2 . With the solution for R (2) L is known, we combine results for Φ(φ), Θ(θ), given by (25) and (26), to obtain the corresponding Debye potential, Π (2) (r), in the form where L = L(ℓ) is given by (32) and µ ℓ , a m , b m are arbitrary and as yet unknown constants to be determined later. This solution is well-known and can be studied with available analytical tools (e.g., [26]).
Examining (20), we see that Π (2) (r) is a solution to the following wave equation: which is the equation for the "free" Debye potential, Π (2) (r), which is "unperturbed" by the short-range plasma potential, V sr .

Eikonal wavefunction
We may now proceed with solving (20), given the relevant form of V sr , (29), first representing this equation as To apply the eikonal approximation to solve this equation, we consider a trial solution of (38) in the form where Π (2) (r) is the "free" Debye potential (36). In other words, in the eikonal approximation the Debye potential Π (2) (r), becomes "distorted" in the presence of the potential V sr given in Eq. (29), by φ, a slowly varying function of r, such that When substituted in (38), the trial solution (39) yields As Π (2) (r) is the solution of the homogeneous equation (37), the first term in (41) is zero. Then, neglecting the third term because of (40), we have As we discussed above, the plasma contribution is rather small and it is sufficient to keep the terms that are of the first order in ω 2 p /ω 2 . Thus, to formally solve (42) we may present the solution for Π (2) (r) from (36) in the series of the small parameter µ/ℓ, which enters there via index L as shown by (32). Then, to solve (42) it is sufficient to take only the zeroth order term (i.e., with µ = 0) in Π (2) (r). It is easier, however, to obtain such a solution directly from (37) by setting µ = 0, which yields the well-known free space solution, e i(k·r) . As a result, to the accuracy needed to solve (42), we have Therefore, with (43) in hand, we may present (42) as We may now compute the eikonal phase. For this, we need to introduce a derivative along the propagation path. As we discussed in [11], we represent the unperturbed trajectory of a photon as where n is the unit vector on the incident direction of the photon's propagation path, such that k = kn, n 2 = 1, and r 0 represents the starting point (see Fig. 1). Following [11,42,43], we define b = [[n × r 0 ] × n] to be the impact parameter of the unperturbed trajectory of the light ray. The vector b is directed from the origin of the coordinate system toward the point of the closest approach of the unperturbed path of light ray to that origin. We will use the z = z(t) coordinate of the ray as the parameter along the path: which may take both negative and positive signs. These quantities allow us to rewrite (45) as As was shown in [11], the differential operator on the left side of (44) is the derivative along the photon's propagation path, namely (k·∇) = kd/dz, where z being the parameter taken along the path, which from (47) is given as r = (b, z). As a result, for (44) we have the solutions of which are That is, we have the following two particular eikonal solutions of (38) for Π(r): where we introduced the eikonal phase Given V sr (r) from (29), we reduced the problem to evaluating a single integral to determine the Debye potentiual Π(r) from (39), which is a great simplification of the problem. Given the fact that b is constant and by taking the short-range plasma potential V sr (r) from (29), we evaluate (51) as where we introduced the function Q βi (z), which, with z = (k · r) = √ r 2 − b 2 , is given as being the hypergeometric function [44]. For r = b, the function (53) is well-defined, taking the value of Q βi (0) = 0, for each β i . For large values of r, for any given value of β i , the function Q βi (z) rapidly approaches a limit: where B[x, y] is Euler's beta function. For the values of β i used in the model (4) for the electron number density in the solar corona, these values are: Note that the quantities Q βi for β i > 2 are always small, 0 ≤ |Q βi | < 1, and as functions of r, they reach their asymptotic values quite rapidly after r ≃ 3.2b. Next, we place the source at a large distance from the Sun: |z 0 | ≫ R ⋆ (see Fig. 1). Then, from definition (53) and the asymptotic behavior given by (54), we have Q βi (z 0 ) = −Q ⋆ βi . As a result, we express the total eikonal phase shift acquired by the wave along its path through the solar system (52) as Expression (56) is the total phase shift induced by the short-range plasma potential along the entire path of the EM wave as it propagates through the solar system. One may see that, as the light propagates from the source to the point of closest approach to the Sun, it acquires the first part of the phase shift, i.e., the term proportional to Q ⋆ βi in (56). As it continues to propagate, the second term in (56) kicks-in, providing an additional contribution.
Substituting the solution that we obtained for the total eikonal phase shift ξ b (r) of (56) in (50) results in the desired solution for the Debye potential Π(r). Effectively, this solution demonstrates that the phase of the EM wave is modified by the short-range plasma potential, as expected from the eikonal approximation. Although (50) is the solution to (38), it still has arbitrary constants µ ℓ , a m , b m present in (36). These constants must be chosen to satisfy a particular boundary value problem that we set out to solve: Obtaining the solution for the EM field as it propagates through the solar system with the refractive medium given by (3).

B. Solution for the radial function RL(r)
We now proceed, following [11], to solve (20) with the help of (21). A particular solution Π is obtained by multiplying together the functions given by (25), (26) and that for R L from (30). We then obtain a general solution to (20). Thus, if the solution for R L is known, we may combine results for Φ(φ), Θ(θ), given by (25) and (26), to obtain the Debye potential in the form where L = L(ℓ) is given by (32) and µ ℓ , a m , b m are arbitrary and as yet unknown constants. Thus, our immediate task is to solve (30).
In the plasma-free case, the entire plasma potential V p is absent, thus L = ℓ. The solution in this case is known and, for the case of EM waves diffracted (i.e., obscured) by a large sphere, was given in [13]. In this case, in order to determine the coefficients µ ℓ in (57), we choose R ℓ (r) to be the regular Bessel function ψ ℓ (kr), and require the resulting EM field to match the incident plane EM wave. As a result, in the vacuum, the solutions for the electric and magnetic potentials of the incident wave, e Π 0 and m Π 0 , may be given in terms of a single potential Π 0 (r, θ) (see [11] for details): Considering the plasma, we notice that, for large r, the potential V sr (r) in (30) can be neglected and this equation reduces to the vacuum discussed in [13] with the solution given by (58). The solution of (30) that is regular at the origin can thus be written asymptotically as a linear combination of the regular and irregular Riccati-Bessel functions ψ L (kr) and χ L (kr), respectively [20,21,45,46], which are solutions of (30) in the absence of the potential V sr (r). Asymptotically these functions behave as Hence, we look for a solution satisfying the boundary conditions where n is a normalization factor. The quantity δ ℓ introduced in these equations is the phase shift due to the presence of the short-range potential V sr (r). We note that δ ℓ vanishes when the short-range potential is absent and, thus, it contains all the information necessary to describe the scattering by V sr (r). We generalize expressions (61)-(62) in the case when the plasma potential has an additional 1/r 2 term, leading to the substitution ℓ → L.
We can satisfy the conditions (61)-(62) by choosing the function R L (r) in a form of a linear combination of the two solutions (50). Furthermore, we will rely on the two solutions in the form of the incoming and outgoing waves [47], given by the functions ζ where ζ [13]) with their asymptotic behavior given by In addition, ξ b (r) is the eikonal phase shift that is accumulated by the EM wave starting from the point of closest approach, r = b. The expression for the quantity is obtained directly from (51) by setting z 0 = 0 (or, equivalently, from (56) by dropping the Q ⋆ βi -term), which results in Clearly, (63) satisfies the radial equation (30) with condition (40). One may verify that R L in the form of (63) satisfies the conditions (61)-(62). Indeed, as the plasma potential exists only for R ⊙ ≤ r ≤ R ⋆ (which is evident from (1) and (3)), the eikonal phase ξ b is zero for r < R ⊙ . Therefore, as r → 0, the index L → ℓ and the radial function (63) becomes R L (r) → ψ ℓ (kr). However, we know that ψ ℓ (kr) obeys the condition (61). Next, we consider another limit, when r → ∞. Using the asymptotic behavior of ζ (±) L from (64), we see that, as r → ∞, the radial function obeys the asymptotic condition (62), taking the form where the phase shift δ ℓ is given by the eikonal phase ξ b introduced by (51). As a result, we established that the radial function (63) represents a desirable solution to (30) inside the termination shock, 0 ≤ r ≤ R ⋆ .
We may further simplify the result (63), putting it in the following equivalent form: which explicitly shows the phase shift induced by the short-range plasma potential.
To match the potentials (58) with those of the incident and scattered waves, the latter must be expressed in a similar series form, but with arbitrary coefficients. Only the function ψ L (kr) may be used in the expression for the potential inside the sphere since χ L (kr) diverges at the origin. On the other hand, the scattered wave must vanish at infinity, and the functions ζ (+) L (kr) will impart precisely this property. These functions are suitable as representations of a scattered wave. For large values of the argument (kr), the result behaves as e ikr and the Debye potential will be Π ∝ e ikr /r for large r. Thus, at large distances from the sphere (i.e., beyond the termination shock boundary) the scattered wave is spherical, with its center at the origin r = 0. Accordingly, it will be used in the expression for the scattered wave, i.e., in the trial solution for the Debye potentials of the scattered wave for r > R ⋆ .
Collecting results for the functions Φ(φ), Θ(θ), given by (25), (26), correspondingly, and R L (r) = ζ (+) L (kr)e iξ b (r) from (50), we have the Debye potential for the scattered wave: where a ℓ , a ′ m , b ′ m are arbitrary and yet unknown constants and the relation between L and ℓ is given by (32). Representing the potential inside the sphere via ψ L (kr) is appropriate. Hence the trial solution to (20) for the electric and magnetic Debye potentials inside the termination shock (r ≤ R ⋆ ) relies on the radial function R L (r) given by (66) and has the form where b ℓ , a m , b m are arbitrary and yet unknown constants. We recall that Π relates to the electric, e Π/ √ ǫ, or magnetic, µ , the expression that was introduced after (20). The boundary conditions (A30) set at the termination shock boundary, r = R ⋆ , for the electron plasma distribution (1) and (3) with n 0 = 0 and thus with ǫ(R ⋆ ) = µ(R ⋆ ) = 1 written more fully are From the symmetry of the geometry of the problem [26] and by applying the boundary conditions (69)-(72), we would have to set the constants a m and b m for the electric Debye potentials as a 1 = 1, b 1 = 0, a m = b m = 0 for m ≥ 2; and for the magnetic Debye potentials as a 1 = 0, b 1 = 1, a m = b m = 0 for m ≥ 2, with identical values for a ′ m and b ′ m . As a result, the solutions for the electric and magnetic potentials of the scattered wave, e Π s and m Π s , may be given in terms of a single potential Π s (r, θ) (see [11] for details), which is In a relevant scattering scenario, the EM wave and the Sun are well separated initially, so the Debye potential for the incoming wave can be expected to have the same form as in the pure plasma-free case that is given by (58). Therefore, the Debye potential for the inner region has the form: with the potential Π in given as We thus expressed all the potentials in the series form given in (57), and any unknown constants can now be determined easily. If we substitute the expressions (58), (73) and (74)-(75) into the boundary conditions (69)-(72), we obtain the following linear relationships between the coefficients a ℓ and b ℓ from (73) and (75), correspondingly: where R L (r) is from (66) and ′ = d/dr. From the definition of the eikonal phase (51) we see that For the electron plasma distribution (1) and (3), the value ξ ′ b (R ⋆ ) is extremely small and may be neglected. Rescaling, for convenience, a ℓ and b ℓ by introducing α ℓ and β ℓ as from (76)-(77), we have: Equations (80)-(81) may now be solved to determine the two sets of coefficients α ℓ and β ℓ : Taking into account the asymptotic behavior of all the functions involved, namely (64) for ζ (+) L and (59)-(60) for ψ L and χ L , we have the following solution for the coefficients α ℓ and β ℓ : with ξ ⋆ b = ξ b (R ⋆ ) and δ * ℓ being the phase shift induced by the plasma to the phase of the EM wave propagating through the solar system as measured at the termination shock, δ * ℓ = δ ℓ (z ⋆ ). As expected, when the plasma is absent, L = ℓ and ξ b = 0, the total plasma phase shift vanishes, resulting in δ ℓ = 0. However, in the case of scattering by the plasma, ξ * b = ξ b (R ⋆ ) = 0 and δ ℓ is important. Also, for large heliocentric distances along the incident direction, for which r ≫ b, and certainly for the region outside the termination shock, r > R ⋆ , the eikonal phase shift ξ ⋆ b = ξ b (R ⋆ ), given by (65), together with (54), will be which, for any given b, is a constant value. In the case when µ/ℓ ≪ 1 and (33) is valid, expression (84) for the plasma-induced delay, to order O(µ 4 ), takes the form We can evaluate the contribution of the plasma to the phase of the EM wave as the wave traverses the solar system. In the case of the electron number density model (4) and from (84), the plasma phase shift δ * ℓ in (86) is given as with η 2 , η 6 and η 16 having the form where, to derive the expression for η 2 , we used µ 2 from (29) and approximated (84) for the case of µ/ℓ ≪ 1 by using (33) with Q ⋆ βi in the incident direction is given by (55). Note that this approach results in the additional factor of π/2 (which came from the first term in (84)) that is characteristic to the eikonal approximation (see discussion in [20,21]). To derive η 6 and η 16 , we used (85). The empirical model for the free electron number density in the solar corona (4) results in the following value for the constants η 2 , η 6 and η 16 in (88): Starting from b ≃ 3.65R ⊙ , the contribution from the η 2 term becomes the largest among the terms in (87) and rapidly becomes dominant over the remaining terms for larger impact parameters. With this analysis and by using the value for a ℓ from (79), together with α ℓ and δ ⋆ ℓ from (84), we determine that the solution for the scattered potential (73), for r > R ⋆ , takes the form Similarly, substituting the value for b ℓ from (79), together with β ℓ from (84), we determine the solution for the inner Debye potential (84) in the form As the solar plasma effect is rather weak, we may use the asymptotic expressions for ψ L , χ L and ζ (±) L for r ≥ R ⊙ . Therefore, the radial function R L (r) from (63) (or, equivalently, from (66)), in the region of heliocentric distances R ⊙ ≤ r ≤ R ⋆ , may be given as where δ ℓ (r) has the form given by (84) where the eikonal phase ξ b is given in its original form (65) which depends on the heliocentric distance, namely Similarly to (86), in the case when µ/ℓ ≪ 1 and (33) is valid, expression (93), to O(µ 4 ), takes the form Thus, in the eikonal approximation, distance dependence in the plasma delay comes from the terms in the short-range plasma potential V sr for which i > 2. The term with i = 2 (i.e., the first term in (86) and (94)) provides no distance dependence. The physical interpretation of this observation follows, in part, from (4). It can be seen that near the solar surface, b R ⊙ , the potential is dominated by terms containing higher powers of (R ⊙ /r). The inverse square term contributes an approximately uniform background potential that, at these small heliocentric ranges, is several orders of magnitude smaller compared to the other terms (see Fig. 2). This inverse square term becomes dominant only at greater distances from the Sun, where the scattering potential is already several orders of magnitude smaller than it is near the solar surface. As a result, outside the Sun, we may present (91) in the following equivalent form: Note that this solution is valid, in principle, even inside the opaque Sun. Indeed, because of the plasma model (1) and (3), the phase shift vanishes, δ ℓ = 0, and (95) reduces to the plasma-free solution (58). The first term in (95) is the Debye potential of an EM wave propagating in a vacuum but modified by the plasma in the solar system. The second term represents the effect of the solar plasma on the propagation of the EM waves inside the termination shock, 0 < r ≤ R ⋆ . Notice that, as the distance increases, this term approaches the form of the Debye potential Π s for the scattered EM field given by (90). Proper accounting for such a dependence makes it possible to compare high-precision observations conducted from different locations within the solar system. Thus, we have identified all the Debye potentials involved in the Mie problem, namely the potential Π 0 given by (58) representing the incident EM field, the potential Π s from (90) describing the scattered EM field outside the termination shock, r > R ⋆ , and the potential Π in from (95) describing it inside the termination shock, 0 < r ≤ R ⋆ .

V. DIFFRACTION OF LIGHT OUTSIDE THE HELIOSPHERE
Solutions for the Debye potentials for the scattered wave outside and inside the termination shock, given by (90) and (95), respectively, describe the scattering of light by the extended solar corona. The presence of the Sun itself is not yet captured. For this, we need to set additional boundary conditions that describe the interaction of the Sun with the incident radiation.
A. The fully absorbing boundary conditions Boundary conditions representing the opaque Sun were introduced in [48] and were used in [11]. Here we use these conditions again. Specifically, to set the boundary conditions, we rely on the semiclassical analogy between the partial momentum, ℓ, and the impact parameter, b, that is given as ℓ = kb [27,28].
We begin with realizing that rays with impact parameter b ≤ R ⊙ are absorbed by the Sun and those with b > R ⊙ are transmitted (see discussion of that in [13]). Therefore, the fully absorbing boundary condition signifies that all the radiation intercepted by the body of the Sun is fully absorbed by it and no reflection or coherent reemission occurs. All intercepted radiation will be transformed into some other forms of energy, notably heat. Thus, we require that no scattered waves exist with impact parameter b ≪ R ⊙ or, equivalently, for ℓ ≤ kR ⊙ . It means that we need to subtract the scattered waves from the incident wave for ℓ ≤ kR ⊙ , as was discussed in [11].

B. The Debye potential for the scattered wave
To implement the boundary conditions for the EM wave outside the termination shock, we realize that the total EM field in this region is given as the sum of the incident and scattered waves, Π = Π 0 + Π s , with potentials given by (58) and (90), correspondingly. Also, using the asymptotic behavior of ζ (+) L from (64) and with the help of (84), we notice that at large distances from the Sun we can write ζ ℓ (similar to derivations of (92)). As a result, we present the Debye potential in the region r > R ⋆ in the following form: Next, relying on the representation of the Riccati-Bessel function ψ ℓ via incoming, ζ (+) ℓ , and outgoing, ζ [13] and also after (34)), we may now express the Debye potential (96) as This form of the combined Debye potential is convenient for implementing the fully absorbing boundary conditions discussed in Sec. V A. Specifically, subtracting from (97) the outgoing wave (i.e., ∝ ζ (+) ℓ ) for the impact parameters b ≤ R ⊙ or equivalently for ℓ ∈ [1, kR ⊙ ], we have or, equivalently, coming back to the form (96) This is our main result, valid for all distances and angles. It is a rather complex expression. It requires the tools of numerical analysis to fully explore the resulting EM field [49][50][51]. However, in most practically important applications we need to know the field in the forward direction. Furthermore, our main interest is to study the largest plasma impact on light propagation, which corresponds to the smallest values of the impact parameter. In this situation, we may simplify the result (99) by taking into account the asymptotic behavior of the function ζ (+) ℓ (kr), considering the field at large heliocentric distances, such that kr ≫ ℓ, where ℓ is the order of the Riccati-Bessel function (see p. 631 of [52]). For kr → ∞ and also for r ≫ r t = ℓ(ℓ + 1)/k (see [11]), such an expression is given by (B16) as which includes the contribution from the centrifugal potential in the radial equation (27) (see e.g., Appendix B or Appendix A in [13] or in [49]). As a result, expression (100) extends the argument of (64) to shorter distances, closer to the turning point of the potential (see the relevant discussion in Appendix F of [11]). By including the extended centrifugal term in (100) we can now properly describe the contribution of the plasma potential to the bending of the light ray's trajectory.
As a result, we may take the approximate behavior of ζ (+) ℓ (kr) in (100) and use (99) to present the solution for the Debye potential outside the termination shock, r > R ⋆ as The first term in (101), Π 0 (r, θ), is the Debye potential of the pure vacuum case that may be derived from the solution obtained in [11]. For practical purposes, however, it is convenient to use an exact expression for Π 0 , which may be given as [13] Π 0 (r, θ) = E 0 k 2 r sin θ e ikr cos θ − e ikr + 1 which is the Debye potential of the plasma-free EM wave. This solution is always finite and is valid for any angle θ. The second term, Π bc (r, θ), is due to the fully absorbing boundary conditions and is responsible for the geometric shadow behind the Sun. The third term, Π p (r, θ), quantifies the contribution of the solar plasma on the scattering of the EM propagating through the solar system, and evaluated at the distance r > R ⋆ . Note that because of the plasma model (1) and (3), the last sum in (101) formally extends only to ℓ = kR ⋆ corresponding to the impact parameter equal to the distance to the termination shock. As expected, for r > R ⋆ , the phase shift δ ℓ = 0 and the entire plasma-scattered wave vanishes.
With the solution for the Debye potential given in (101), and with the help of (10)-(15) (also see [11]), we may now compute the EM field in the various regions involved. Given the value of (ω p /ω) 2 (∼ 10 −2 for radio and ∼ ×10 −11 for optical wavelengths), we may neglect the effect of the solar plasma, behaving as ∝ 1/r 2 , on the amplitude of the EM wave. This is especially true at large heliocentric distances where the effect of the plasma on the amplitude of the EM wave is negligibly small. The plasma will affect the delay of the EM wave, which is fully accounted for by the solution for the Debye potentials. Thus, we can put ǫ = µ = 1 in (10)- (15) and use the following expressions to construct the EM field in the static, spherically symmetric geometry (see details in [11]): with quantities α, β and γ are computed from the known Debye potential, Π, as β(r, θ) = 1 r γ(r, θ) = 1 r sin θ ∂ rΠ ∂r + ik r ∂ rΠ ∂θ .
We may now use these equations to compute the EM field in the region beyond the termination shock.

C. EM field in the shadow region
In the shadow behind the Sun (i.e., for impact parameters b ≤ R ⊙ ) the EM field is represented by the Debye potential of the shadow, Π sh , which is given by (101) as Π sh (r, θ) = Π 0 (r, θ) + Π bc (r, θ) = Π 0 (r, θ) + e ikr r As was already discussed in [13], there is no EM field in the geometric shadow behind the Sun and, thus, no light.

D. EM field outside the shadow
In the region outside the solar shadow (i.e., for light rays with impact parameters b > R ⊙ ), the EM field is derived from the Debye potential given by the remaining terms in (101) as Expression (108) is our main result for the region outside the termination shock, r > R ⋆ . It contains all the information needed to describe the total EM field originating from a plane wave that passed through the region of the static spherically symmetric plasma of the extended solar corona, characterized by an electron number density (3) that diminishes as r −2 or faster.
To evaluate this potential, we use the exact solution for the potential Π 0 (r, θ) given by (102) and present the total solution for the Debye potential (108) in the area outside the solar shadow as follows: where the plasma scattering amplitude f p (r, θ) is given by We note that because of the contribution from the centrifugal potential in (100), the scattering amplitude, f p (r, θ), now also a function of distance, which is not the case in the typical problems describing nuclear and atomic scattering [27][28][29]53]. However, as we observed in [11,13] and as we shall see later here, in the case when we are interested in the trajectories of the particles, the presence of this term in the phase of the scattering amplitude (110) allows us to properly describe the bending of the light rays in the dispersive medium of the plasma.
As a result, the Debye potential for the plasma-scattered wave outside the termination shock takes the form with the plasma scattering amplitude f p (r, θ) given by (110). We will use (111) to derive the components of the EM field produced by this wave. For this, we use (111) in the expressions (104)-(106) to derive To establish the solution for the EM field, we use the scattering amplitude f (r, θ) from (110) to evaluate each of the expressions in (112)-(114) and substitute the results in (103). E. Solution for the function α(r, θ) and the radial components of the EM field We begin with the investigation of α(r, θ) from (104), which, in terms of the scattering amplitude, f p (r, θ), is given by (112). Thus, by using f p (r, θ) from (110), we obtain the following expression for α(r, θ): To evaluate this expression, we use the asymptotic representation for P (1) l (cos θ) from [49,54,55], valid when ℓ → ∞: This approximation can be used to transform (115) as At this point, we may replace the sum in (117) with an integral: and evaluate this integral by the method of stationary phase (see [11,48]). This method allows us to evaluate integrals of the type where the amplitude A(ℓ) is a slowly varying function of ℓ, while ϕ(ℓ) is a rapidly varying function of ℓ. The integral (119) may be replaced, to good approximation, with a sum over the points of stationary phase, ℓ 0 ∈ {ℓ 1,2,.. }, for which dϕ/dℓ = 0. Defining ϕ ′′ = d 2 ϕ/dℓ 2 , we obtain the integral Because the scattering term e i2δ * ℓ − 1 in (118) provides two contributions, each with a different phase, we will treat the integral (118) as the sum of two integrals: one with and one without the contribution from the plasma phase shift 2δ ⋆ ℓ . To demonstrate our approach, we begin with the plasma-free case.

Evaluating the plasma-free term
For the term in (118) without the plasma phase shift, the relevant ℓ-dependent part of the phase is of the form [13] ϕ [0] The phase is stationary when dϕ [0] ± /dℓ = 0, which implies therefore, we may write the solution for the points of stationary phase Note that by extending the asymptotic expansion of ζ (+) ℓ (kr) from (100) to the order of O((kr) −(2n+1) ) (i.e., using the WKB approximation as was done in Appendix B), the validity of the result (123) extends to the order of O(θ 2n+1 ).
To calculate ϕ ′′ (ℓ) to the order of O(θ 6 ) as in (124), we need to include in the phase ϕ [0] ± (ℓ) (121) another term ∝ ℓ 6 , which may be taken from (B16). This allows to compute ϕ ′′ (ℓ 0 ) The remaining integral is easy to evaluate using the method of stationary phase. Before we do that, we need to bring in the amplitude factor for the asymptotic expansion ζ (+) ℓ (kr) given by (100). This factor, which we denote by a(ℓ), is readily available from (B16) in the following form: The fact that we did not use it in (100) does not affect results of the calculations above. However, as we shall see below, its presence is needed to offset some of the terms that are present in the phase of (121). The significance of this term is realized in the fact that, for the method of stationary phase, it cancels out the contribution of the θ-dependence in (125), namely using result (126) we derive Note that by improving the asymptotic expansion of ζ (+) ℓ (kr) (100) (that, in a more complete form, is given by (B16)) to a higher order and extending the phase, from (100), to the order of O((kr) −(2n+1) ) and the amplitude, a(ℓ), from (126), to the order of O((kr) −2n ), the validity of (127) will extend to the order of O(θ 2n ). If needed, this can be achieved by following the derivations presented in Appendix B. Now, using (127), we have the amplitude of the integrand in (118), for ℓ ≫ 1, taking the form where the superscript [0] denotes the term with no plasma contribution. We can drop the 1/(ikr) term in the parentheses of this expression, as it is 1/(kr) times smaller in magnitude compared to the leading term. As a result, the plasma-free part of the expression for the δα [0] (r, θ) from (118) takes the form As before, we note that the validity of this result my be extended to the order of O(θ 2n+1 ) by extending the asymptotic expansion of ζ (+) ℓ (kr) from (100) to the order of O((kr) −(2n+1) ). This can be done by following the approach that led to (B16) and repeating the derivations presented in this section.

Evaluating the term with plasma contribution
We now turn our attention to the term in (118) that contains a plasma phase shift contribution. The relevant ℓ-dependent part of the phase is given as with the plasma contribution clearly shown. From the definition (85) and (86), this plasma phase shift is given as where we used representation of Q ⋆ βi via Euler's beta function, as shown by (54). The phase (130) is stationary when dϕ [p] ± /dℓ = 0, which, similarly to (122)-(123), implies where δθ b = d2δ * ℓ /dℓ is the semiclassical angle of deflection of a light ray [53]. This angle may be computed from (85) and (86) by taking into account the semiclassical relation ℓ = kb. As a result, the light deflection angle 5 is computed to be Note that this expression for the plasma deflection angle, δθ b , is identical to that obtained in [6,7] and used for the recent Cassini experiment [56]. As the earlier result was obtained under a different set of physical assumptions and mathematical tools, such correspondence confirms the validity of our approach, which relies on the wave-optical treatment of the problem and advocates for a direct solution of Maxwell's equations, while treating the plasma contribution to the phase of the EM wave via the eikonal approximation.
With the results given in (132), we may compute the needed expressions for the value of the phase along the path of stationary phase: and for the second derivative of the phase along the same path, simialry to (125), we have Using (133), we estimate the magnitude of the second term in this expression: Evaluating this quantity with the values from the empirical model (4), we see that for the smallest impact parameter b = R ⊙ this quantity takes the largest value of d 2 (2δ ⋆ b )/dℓ 2 = 1.57 × 10 −26 (λ/(1 µm) 3 . This results in the fact that for optical wavelengths, even at the heliocentric distance of r ≃ 6.5 × 10 3 AU, this term will be over 10 4 times smaller than the 1/(kr) term in (135), representing a small correction to ϕ ′′ (ℓ 0 ) that may be neglected for our purposes. This is equivalent to treating the deflection angle δθ b constant, which is consistent with the eikonal (or high-energy) approximation [36][37][38][39][40][41].
As a result, the expression for the second derivative of the phase from (135) takes the form The relevant, plasma-dependent part in the integral in (118) now is easy to evaluate using the method of stationary phase. Similarly to (128), we have the amplitude of the plasma-dependent term in (118) evaluated to be where the superscript [p] denotes the term due to the plasma phase shift. As it was done with (128), we can drop the 1/(ikr) term as it is much smaller compared to the leading term. Also, using the fact that sin θ ≃ b/r, we may evaluate the expression Considering (133), we see that the largest value of the bending angle, δθ b , is reached at the smallest impact parameters, b = R ⊙ , limiting the size of this angle as δθ b (R ⊙ ) ≤ 1.73 × 10 −12 λ/1 µm 2 rad. As such, for optical wavelengths the second term in (139) becomes significant only beyond 10 9 AU, which is beyond any practical significance. Therefore, we will neglect this term in our further considerations. As a result, similarly to (129), we obtain the contribution of the plasma-dependent term in (118) in the form With the expressions (129) and (140) at hand, we may now present the quantity α(r, θ) from (118) as Using the approach presented above, we may now evaluate the scattering efficiency factors β(r, θ) and γ(r, θ) needed to determine other components of the EM field. To investigate the behavior β(r, θ) from (105) (which, in terms of the scattering amplitude, f p (r, θ), is given by (113)), we will use the scattering amplitude f p (r, θ) from (110) and obtain the following expression for β(r, θ): To evaluate the magnitude of the function β(r, θ), we need to establish the asymptotic behavior of P l (cos θ)/ sin θ and ∂P (1) l (cos θ)/∂θ. For fixed θ and ℓ → ∞ this behavior is given 6 [50] as (this can be obtained directly from (116)): With these approximations, the function β(r, θ) in the region outside the geometric shadow (i.e., not on the optical axis), takes the following form: For large ℓ ≫ 1, the first term in the curly brackets in (145) dominates, so that this expression may be given as To evaluate β(r, θ) from the expression (146), we again use the method of stationary phase. For this, representing (146) in the form of an integral over ℓ, we have: As we have done with (118), we treat this integral as a sum of two integrals: a plasma-free and a plasma-dependent term. Expression (147) shows that the ℓ-dependent parts of the phase have a structure identical to (121) and (130). Therefore, the same solutions for the points of stationary phase apply. As a result, using (123) and (125), from (147) for the part of the integral that does not depend on the phase shift, δ ⋆ ℓ , we have As a result, similarly to (129), the expression for the δβ Next, using the ℓ-dependent phase (130) with the plasma phase shift included and the relevant expressions (134) and (137) we have: Thus, the term in (147) that depends on the contribution of the plasma-induced phase shift, namely δβ ± (r, θ), takes the form Using the expressions (149) and (151), we present the integral (147) as β(r, θ) = δβ To determine the remaining components of the EM field (103), we need to evaluate the behavior of the function γ(r, θ) from (106), which, in terms of the scattering amplitude f p (r, θ), is given by (114). Using f p (r, θ) from (110), we substitute it in (114) ad obtain the following expression for γ(r, θ): To evaluate this expression, we use the asymptotic behavior of P (1) l (cos θ)/ sin θ and ∂P (1) l (cos θ)/∂θ given by (143) and (144), correspondingly, and rely on the method of stationary phase. Similarly to (145), we will drop the second term in the curly brackets in (153). The remaining expression for γ(r, θ) will now be determined by evaluating the following integral: Clearly, this expression yields the same equation to determine the points of the stationary phase (121) and (130) and, thus, all the relevant results obtained in Sec. V E. Therefore, the ℓ-dependent amplitude of (154), which is independent on plasma phase shift, A [0] , is evaluated as Therefore, the plasma independent part of the function δγ Similarly, we have the following expression for the amplitude A [p] which with (134) and (137) is evaluated to be Therefore, the plasma-dependent term in (154) , δγ We may now use the expressions (156) and (158) to present the integral (154) as H. Solution for the EM field outside the termination shock At this point, we have all the necessary ingredients to present the ultimate solution for the scattered EM field in the eikonal approximation. To determine the components of the EM field, we use the expressions that we obtained for the functions α(r, θ), β(r, θ) and γ(r, θ), which are given by (141), (152) and (159), correspondingly, and substitute them in (103). As a result, we establish the solution for the scattered EM field in the region outside the solar system: At the same time, the EM field produced by the Debye potential Π 0 of the incident wave (102) has the form [13]: The total field in accord with (109) is given by the sum of (161)-(162) and (163)-(164): One can see that the total EM field behind the very large sphere, λ ≪ R ⊙ , embedded in the spherically symmetric plasma distribution, has the structure similar to the incident EM wave. However its phase and propagation direction are affected by the delay introduced by the plasma in the solar system. The EM field outside the termination shock takes a very simple form that depends on the plasma phase shift, δ * ℓ . This phase shift given by (87) is clearly showing its dependence on the ray's solar impact parameter. Eqs. (165)-(167) account for this contribution. We note that plasma defocuses light rays, introducing aberrations that are small at optical wavelengths, but are significant in the radio bands. Images taken by conventional astronomical instruments will be affected by the static plasma. Although temporal variability in the plasma may introduce additional aberrations, at optical wavelengths such effects are small and may be accounted for with standard observational techniques [10,14].

VI. DIFFRACTION OF LIGHT WITHIN THE HELIOSPHERE
A. The Debye potential inside the termination shock To establish the solution for the Debye potential inside the termination shock, we need to implement the fully absorbing boundary conditions, as we did in Sec. V. To do this, we identically rewrite (95) using a representation of the Riccati-Bessel function ψ ℓ via incoming and outgoing waves, ζ By removing from this expression the outgoing waves corresponding to the impact parameters b ≤ R ⊙ or, equivalently, for ℓ ∈ [1, kR ⊙ ], we implement the fully absorbing boundary conditions that account for the physical properties of the solar surface. The resulting expression has the form Expression (169) is our solution for the Debye potential representing the EM field that travels through the solar system in the presence of the solar plasma. Reinstating the Riccati-Bessel function ψ ℓ and, similarly to (101), taking into account the asymptotic behavior of the function ζ (+) ℓ (kr) given by (100), expression (169) representing the solution for the Debye potential within the solar system takes the form This is our main result for the Debye potential representing the EM field inside the termination shock, 0 < r ≤ R ⋆ . It describes the propagation of monochromatic EM waves on the background of a generic static spherically symmetric plasma for which the number density diminishes as r −2 or faster. The first term in (170) is the Debye potential of the incident wave modified by the plasma as the wave propagates through in the solar system. The second term is for the geometric shadow behind the Sun (similar to that discussed in [13],) also modified by the plasma. The third term represents the ongoing scattering of the EM field as it propagates through the spherically symmetric distribution of the extended solar corona, given at a particular heliocentric position within solar system, R ⊙ ≤ r ≤ R ⋆ .
Note that because of the plasma model (1) and (3), the last sum in (170) formally extends only to ℓ = kR ⋆ , corresponding to the impact parameter equal to the distance to the termination shock. For r > R ⋆ , not only does the vanishing phase shift, δ ℓ = 0, essentially eliminate this term, but this distance is also outside the boundary that characterizes the inner region, as for r > R ⋆ we enter the domain of the scattered wave discussed in Sec. V.
To discuss the diffraction of light in the solar system, we refer to the solution for the Debye potential Π given by (170). Each term in (170) has the contribution of the ongoing plasma phase shift given as δ * ℓ − δ ℓ (r), where δ * ℓ and δ ℓ (r) are given by (84) and (93), correspondingly, with eikonal phase shifts for the short-range plasma potential ξ b (r) and ξ ⋆ b are given by (65) and (85), respectively. To evaluate these terms, we derive the differential plasma-induced phase shift occurring as the wave travels through the heliocentric ranges R ⊙ ≤ r ≤ R ⋆ . Defining from (65) and (85) we compute As we can see, the quantity δξ b (r) is independent of either L or ℓ and is a function of the heliocentric distance only. Thus, in all the terms of (170) we may move the factor exp i(δ * ℓ − δ ℓ (r)) ≡ exp iδξ b (r) outside the summation over ℓ. Clearly, δξ b (r) is significant only in the immediate vicinity of the Sun, where r ≃ R ⊙ , but it falls off rapidly for larger distances. Therefore, for optical or near-IR wavelengths, any phase variability introduced by this term may be neglected for larger heliocentric distances.
With this, the solution for the Debye potential Π from (170) takes the form To investigate the solution (173) we already have all the necessary expressions. Indeed, in the shadow behind the Sun (i.e, for impact parameters b ≤ R ⊙ ), similarly to (107), the EM field is represented by the Debye potential of the shadow, Π in sh , which is given by the following plasma-corrected expression Π in sh (r, θ) = e iδξ b (r) Π 0 (r, θ) + In complete analogy with (107), we conclude that there is no EM field in the geometric shadow behind the Sun and, thus, no light (see also [13]). Next, similarly to (109)-(110), we use the exact solution for the potential Π 0 (r, θ) from (102) and present the total solution for the Debye potential in the area outside the solar geometric shadow, i.e., for R ⊙ ≤ b ≤ R ⋆ : where the scattering amplitude f in p (r, θ) is given by As we discussed in Sec. V D, the effect of plasma on the amplitude of the EM wave even in close proximity to the Sun is at most (ω p /ω) 2 ∼ 10 −2 for radio and (ω p /ω) 2 ∼ 10 −11 for optical wavelengths. Thus, we may neglect terms of order of (ω p /ω) 2 (E, H) and compute the EM field in the inner solar system by using (104)-(106), with the scattering amplitude from (176). The difference between (176) for the scattering amplitude of the wave moving inside the termination shock and a similar expression (110) for the wave's propagation in the regions outside the termination shock is in the presence of the distance-dependent phase shift, δ ℓ (r) as opposed to δ * ℓ , which depends only on the impact parameter.
Therefore, to compute the EM field for R ⊙ ≤ r ≤ R ⋆ we may use the results obtained for α(r, θ), β(r, θ) and γ(r, θ), given by (141), (152) and (159), correspondingly, where the plasma-induced phase shift δ * ℓ from (87) is replaced by δ ℓ (r) from (65) that has the form with η 2 , η 6 and η 16 given by (88). With this, all the results that we obtained in Sec. V F are directly applicable to the ranges within the termination shock, R ⊙ ≤ r ≤ R ⋆ . As a result, and similar to (160) and (162), we have the components for the scattered field. We then add the expressions for the plasma-corrected incident wave, as it propagates through the extended solar plasma, are derived from (163)-(164) by multiplying them by the factor of e iδξ b (r) . As a result, the total EM field is given by the sum of the plasma-corrected incident wave, and the wave scattered by the spherically symmetric plasma in the extended solar corona, all of which may be derived directly from (165)-(167), by replacing the δ ⋆ ℓ with δ ℓ (r) and the by multiplying each of the expressions by a factor of e iδξ b (r) . As a result, the EM field in the inner solar system is given in the following form Note that the total plasma-induced phase shift in (178)-(180) has the contribution of the eiklonal phase shift is in agreement with (56). Expressions (178)-(180) suggest that, as the EM wave travels through the solar system, its phase acquires a plasma-induced phase shift δ ⋆ ℓ + δ ℓ (r) due to the scattering of light in the spherically symmetric free electron plasma (3); its direction is also changing due to refraction in this plasma. For µ/ℓ ≪ 1, we use (94) and take µ 2 from (29), to express this phase shift in the form that shows its explicit dependence on the plasma model (3): Expression (181) is the plasma-induced phase shift of the EM wave that originated at a distant source as it travels through the solar system at heliocentric ranges R ⊙ ≤ b ≤ r ≤ R ⋆ . Although the last term in this expression depends on r, the contribution of this term is small and may be neglected for observations conducted beyond a few solar radii, Most of the phase shift is accumulated in the vicinity of the Sun, where the plasma density is dominated by contributions from the terms with higher powers of (R⊙/r) in the empirical solar plasma model (4), also shown in Fig. 2. Far outside the region of interest, caustics form, as observed in [57].
as r 5R ⊙ , as seen in Fig. 2. Therefore, for the total plasma-induced phase shift in (178)-(180), δ * ℓ + δ ℓ (r) ≃ 2δ * ℓ , with 2δ * b given by (131). Using the phenomenological model (4), we estimate the plasma-induced phase delay, d b : which, for Sun-grazing rays (i.e., b = R ⊙ ) results in the plasma phase delay of d b (R ⊙ ) = −1.2 × 10 −4 λ/1 µm 2 m, which is negligible for optical and IR wavelengths, but it is important in the microwave band. By examining (178)-(180), we see that, as the plasma-induced phase delay enters both in the phase of the EM field and it also changes the direction of the ray's propagation. At each point, the wavefront is adjusted by the amount of the plasma-induced phase shift. The angle of such a turn may be computed directly from the delay (131) as d(δ * ℓ + δ ℓ (r))/dℓ ≃ δθ b , which may be computed similarly to (133). To consider the the impact of the solar plasma on the imaging properties of the astronomical telescopes, we need to know the energy flux at the image plane situated at a particular distance from the Sun, which is given by the Poynting vector, S. Taking E and H from (178)-(180), we compute the components of the time-averaged Poynting vector S in the spherical coordinate system (r, θ, φ) as usual [11,58,59]: where the overline, a, denotes time averaging. For practical purposes, instead of the spherical coordinates (r, θ, φ), it is convenient to introduce a cylindrical coordinate system (ρ, φ, z). Following the discussion in Sec. III.F in [11], from (183), to appropriate order, we derive the following expression for the Poynting vector S in cylindrical coordinates: Note, that in the absence of plasma (i.e., when δ ℓ (r) = 0), (184) coincides with the expression for the Poynting vector of a plane wave in the Euclidean geometry, namelyS →S 0 = (c/8π)E 2 0 (0, 0, 1), corresponding to the results obtained with the geometric optics approximation [26], as well as with the wave-optical treatment of the shadow cast by a large sphere [13]. One may see that the scattering of light in the solar corona results in the refraction of light that is characterized by the angle δθ b , as expected. As shown on Fig. 3, the presence of plasma changes the shape of the shadow behind the Sun from a cylindrical to a conical shape, starting with a rotational hyperboloid region, with asymptotes characterized by the angle δθ b and eventually (far outside the heliosphere at optical wavelengths) forming caustics, as it was also observed in [57].
The result given in (184) describes the angle of deflection of the light rays as they are scattered by the solar plasma, as shown in Fig. 3. This figure illustrates schematically that the direction of light rays changes as they traverse the solar system. The diagram also shows the surfaces of equal phase, representing the delay introduced by the plasma and the resulting tilt of the EM wavefront. The conical shape of the shadow is controlled by the ray's bending angle. Using the phenomenological model (4) in (133), we estimate the plasma deflection angle, δθ b , as a function of the impact parameter and the wavelength: which suggests that for sungrazing rays (i.e., for the rays with impact parameter b ≃ R ⊙ ), the bending angle (185) reaches the value of δθ b (R ⊙ ) = 1.73 × 10 −12 λ/1 µm 2 rad, which is large for radio wavelengths, but is negligible in optical or IR bands. For typical observing situations with reasonable Sun-Earth-probe separation angles [7][8][9], expression (133) provides a good description. This, once again, justifies the application of the eikonal approximation.

VII. DISCUSSION AND CONCLUSIONS
In this paper, we studied the propagation of EM waves in the vicinity of a very large opaque sphere representing our Sun, and considered the problem of diffraction of light in the presence of the extended solar corona. We considered a free electron plasma with a generic spherically symmetric and static distribution, the electron number density of which covers the entire solar system all the way to the boundary of the heliosphere (3). We developed a wave-optical treatment of the diffraction of light in a realistic model of the plasma present in the solar system and studied the shadow cast by a large sphere in the presence of such plasma, thereby extending the results reported in [13].
In our approach we sought to solve the problem above by considering three regions that cover the solar system, namely i) the Sun itself, providing the fully absorbing boundary conditions (see Sec. V A) that let us study its shadow; ii) the region within the solar system covering heliocentric distances from the solar surface to the edge of the heliosphere, R ⊙ ≤ r ≤ R ⋆ , which yields a description of the light scattering problem that is relevant to all astronomical observations conducted within the solar system (see Sec. VI), and iii) the region that lies outside the heliosphere (see Sec. V), yielding us a description of light scattering in the vicinity of any star, which is relevant for modeling contributions of plasma in extrasolar systems to any type of astronomical observations. We were able to provide a wave-optical treatment of the EM wave propagating in all three regions. For this we used the Mie approach of finding the relevant Debye potentials, which were then used to construct the EM field. We successfully used this approach to treat similar problems of light propagating in the vicinity of a very large sphere but on different backgrounds. Thus, in [13] we studied the scattering of light in flat spacetime and in [11] we studied the same problem, but on the background of a weak gravitational field produced by the Sun. Our current paper extends these earlier results.
Solving the wave equation for the Debye potential in the presence of an arbitrary potential is a challenging problem. In fact, no closed form analytical solutions exist. However, we were able to apply the existing tools developed to address similar problems encountered in describing nuclear and atomic scattering. One of those tools-the eikonal or high-energy approximation (see Secs. IV A-IV B)-was very effective in helping us to find a solution that satisfies our objectives. Using this method, we have shown that the presence of the plasma in the solar system results in a phase shift (131) and related change in the direction of the wave's propagation through the solar system (133).
The magnitude of the phase shift (181) depends on the impact parameter (and weakly on distance). The smaller the light ray's impact parameter, the larger the phase shift and, thus, the tilting of the wavefront. The index of refraction of a free electron plasma is negative, therefore the wavefront tilts outward, away from the Sun (see Fig. 3). This effect needs to be considered when comparing observations conducted at spatially separate locations (e.g., very long baseline interferometry), especially at wavelengths longer than optical or near-IR. In practice, it would manifest itself via a differential phase shift for any two locations with two different solar impact parameters.
In fact, by tilting the wavefront, the static, spherically symmetric plasma (3) introduces optical aberrations and, thus, it affects image formation for astronomical observations. Although this effect is very small at optical wavelengths, it nevertheless may be noticeable in astronomical observations conducted in the inner solar system. Our formalism allows us to account for these effects in any field of practical astronomy. In addition, the brightness of the solar corona [16] is an obvious source of significant background noise for any astronomical observations conducted in the inner solar system and, thus, must be accounted for when discussing sensitivity in practical applications [10,14].
The approach that we developed is derived from first principles and requires no ad hoc assumptions. Therefore, it may be easily augmented, e.g., by introducing latitude-dependent corrections to the plasma model and, if needed, the rotating Sun. For this, one may rely on the formalism to describe the relativistic reference frames from [60] that may be used to treat rotation. In addition, our rigorous treatment can also serve as a basis for further study, including reliable modeling of deviations from spherical symmetry.
One of the natural next steps would be to extend the results that we obtained here to the case of the EM waves propagating in plasma on the curved spacetime background induced by the Sun's gravitational field and to extend the results reported in [11]. This work in ongoing; results will be presented elsewhere.
As a result, we obtain two approximate solutions for the partial radial function R ℓ : R ℓ (r) = c ℓ exp i kS −1 (r) + S 0 (r) + k −1 S 1 (r) + k −2 S 2 (r) + k −3 S 3 (r) + k −4 S 4 (r) + k −5 S 5 (r) + + d ℓ exp − i kS −1 (r) + S 0 (r) + k −1 S 1 (r) + k −2 S 2 (r) + k −3 S 3 (r) + k −4 S 4 (r) + k −5 S 5 (r) = = exp ℓ(ℓ + 1) 4k 2 r 2 − where c ℓ and d ℓ account for all the constants relevant to the point r 0 in (B9)-(B12). Expression (B15) is used in Section V E where we apply the method of the stationary phase to develop expressions containing the scattering amplitude. As we saw previously (e.g., [11,13]) the solution for the points of the stationary phase leads to a solution for ℓ of the form ℓ ≃ kr sin θ. This observation allow us to somewhat simplify the expressions (B15). Indeed, any term for which the exponent of ℓ in the numerator is less than the exponent of kr in the denominator is extremely small compared to the other terms. Indeed, the first and the last term in the amplitude are of order ℓ 2 /(kr) 2 and ℓ 4 /(kr) 4 , correspondingly. However, the second term is of order ℓ 2 /(kr) 4 , which is 1/(kr) 2 times smaller than the first term and 1/ℓ 2 smaller than the third term. Thus, the second term may be neglected. On the same grounds we may neglect three terms in the phase of these expressions (B15). The neglected terms may become important very close to the θ ≃ 0 direction, corresponding to extremely large heliocentric distances, which is beyond the scope of this paper.
Also, similarly to [11], we may further improve the asymptotic expression for R ℓ from (B15) by redefining the constants c ℓ and d ℓ as c ℓ → c ℓ e −i πℓ 2 , d ℓ → d ℓ e i πℓ 2 . This leads to the following expression for the radial function R ℓ : R ℓ (r) = exp ℓ(ℓ + 1) 4k 2 r 2 + [ℓ(ℓ + 1)] 2 8k 4 r 4 × × c ℓ exp i kr − πℓ 2 + ℓ(ℓ + 1) 2kr + ℓ(ℓ + 1) In [11] the asymptotic behavior of the Riccati-Bessel functions was obtained for very larger distances from the turning point for r ≫ r t ; the solution (B16) improves them by extending the argument of these functions to shorter distances, closer to the turning point. Similar expressions may be obtained from the asymptotic expansions of the Riccati-Bessel given as a finite sums which are given in [49,55].