Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Identify the limits of geometric optics ray tracing by numerically solving the vector Kirchhoff integral

Open Access Open Access

Abstract

The properties of a pencil of light as defined approximately in the geometric optics ray tracing method are investigated. The vector Kirchhoff integral is utilized to accurately compute the electromagnetic near field in and around the pencil of light with various beam base sizes, shapes, propagation directions and medium refractive indices. If a pencil of light has geometric mean cross section size of the order p times the wavelength, it can propagate independently to a distance p2 times the wavelength, where most of the beam energy diffuses out of the beam region. This is consistent with a statement that van de Hulst made in a classical text on light scattering. The electromagnetic near fields in the pencil of light are not uniform, have complicated patterns within short distances from the beam base, and the fields tend to converge to Fraunhofer diffraction fields far away from the base.

© 2020 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

1. Introduction

The geometric optics method is an approximate approach that is extensively employed to compute the single-scattering properties of dielectric non-spherical particles (e.g. [18]) with sizes much larger than the incident wavelength, where numerically accurate methods (e.g. [913]) are not computationally feasible. It has been well known that geometric optics methods only give results with acceptable accuracy when the particle size parameter (defined as the particle characteristic radius, divided by the incident wavelength λ, and multiplied by 2π) is much larger than unity (e.g. [1,1415]).

In the geometric optics method, the scattering field is decomposed into the diffracted field, and the fields resulting from reflection and refraction external and internal to the particle. The reflection and refraction contributions are computed by the ray tracing technique, in which electromagnetic waves are approximated as a number of rays. A ray is defined as the path of light within the geometric optics framework. The ray is a geometric line for which the ray cross section is not a concern. A bundle of ray has a finite cross section in a certain shape, which is also called beam or a pencil of light. Hereafter, ‘a pencil of light’ and ‘beam’ are used interchangeably. The electromagnetic field in a beam is assumed to be a plane wave. Individual beams are traced independently within the particle based on Snell’s law and the Fresnel equations. The beams emerging from the particle surface are assumed to propagate to an infinitely far distance to generate the far field. Indeed, the ray tracing technique only considers the wave properties of reflection, refraction, and polarization. A beam only contributes to second moments of the far field along its propagation direction, and the interferences in a single beam, among beams and between a beam and the diffracted field are neglected.

Some efforts have been made (e.g. [1621]) to improve the conventional geometric optics method by applying physical optics based modifications. For example, rigorous electromagnetic surface or volume integrals are applied to compute far fields induced by the near fields obtained with the ray tracing technique inside the particle or on the particle surface [22]. With this modification, interferences among the plane waves representing individual beams and with the diffracted field are implicitly incorporated for computing the far field, and a beam contributes to the far field in all directions. Comparisons with numerically accurate approaches show that this physical optics modification can significantly improve the accuracy of the geometric optics method [17,22]. However, the remaining errors attributed to the near field obtained by the geometric optics method inside the particle or on the particle surface still exist.

van de Hulst [1] states that if the base of a pencil of light has a width with the order , the pencil of light can exist independently over a length of the order p2λ. This statement can be viewed as a corollary of the Fresnel–Huygens principle, which is an approximate mathematical description of the diffraction of a wave.

Consider the propagation of a plane wave in free space with two parallel equal-phase planes P and Q, as shown in Fig. 1. According to the Fresnel–Huygens principle, all fields on plane P serve as point sources of spherical waves. The superposition of the spherical waves gives the field on plane Q, and the envelope of the spherical wave fronts constitutes the equal-phase plane Q [23]. The field at Q1 is expressed by

$${u_{{\textrm{Q}_1}}} = \int\limits_\textrm{P} {{u_\textrm{P}}\frac{{i\exp (ikr)}}{{r\lambda }}K(\chi )dS} ,$$
where uP is the field amplitude on plane P, and can be either the electric or magnetic field component, k is the (real-valued) wavenumber, r is the distance between P0 and Q1 (note that r >> λ), and K(χ) is a direction factor. The angle χ is between P0Q0 and P0Q1, as shown in Fig. 1. The integral kernel in Eq. (1) is a spherical wave emitted by the field on a surface element of P. Because of plane wave propagation, the field amplitude at any point of plane Q is
$${u_\textrm{Q}} = \exp (ikl){u_\textrm{P}},$$
where l is the distance between P (P0) and Q (Q0). If χ = 0, then K(χ) = 1 [23], and Eq. (1) becomes
$${u_{{\textrm{Q}_1}}} = \int\limits_\textrm{P} {{u_\textrm{P}}\frac{{i\exp (ikl)}}{{l\lambda }}dS} .$$
Comparing Eqs. (2) and (3), if the area of dS is equal to , and up is nonzero only on the area of dS, Eqs. (2) and (3) only differ by the imaginary unity i or exp(/2). The field intensity at Q1 can be viewed as contributed to only by the field on the surface element dS with area . The fields on the plane P outside dS tend to cancel out each other at Q1, and contribute to a phase difference of π/2 [1].

 figure: Fig. 1.

Fig. 1. Illustration of Huygens–Fresnel principle.

Download Full Size | PDF

In other words, the field on a region with the area , can propagate independently to a distance l, which is the van de Hulst statement mentioned above. The validity of the van de Hulst statement still needs to be examined using a rigorous method, and quantified to help evaluate the geometric optics method performance.

In this study, we compute the electromagnetic near field of beams by numerically evaluating the vector Kirchhoff integral for various cross sections, and analyze van de Hulst’s statement about the “existence” of a pencil of light in geometric optics ray tracing. Our discussion is limited to a non-absorbing isotropic dielectric medium. The method and computation setup are described in Section 2. Section 3 gives results and analysis. Section 4 concludes the study.

2. Methodology

During ray tracing, the electromagnetic field in a beam is assumed to be a monochromatic time-harmonic plane wave. Its electric field vector E is expressed in complex exponential form as

$${\textbf {E}}({{\textbf {r}}}) = {{\textbf {E}}_0}\textrm{exp}[{i({{{\textbf {k}}} \cdot {{\textbf {r}}} - \omega t} )} ],$$
where E0 is the electric filed amplitude vector that is a constant in the beam, r is a position vector pointing from the origin to the observation point, k is the wavevector, t is time and ω is the angular frequency. The beam propagates along the k direction. The magnitude of k is the wavenumber k = 2πm/λ0, in which m is the refractive index, and λ0 is the wavelength in a vacuum. In this paper we use the Gaussian unit system. The magnetic field H and magnetic induction B are related to the electric field by
$${{\textbf {H}}}({{\textbf {r}}}) = {{\textbf {B}}}({{\textbf {r}}}) = \frac{{\nabla \times {\textbf {E}}({{\textbf {r}}})}}{{i{k_0}}},$$
where k0 is the wavenumber in a vacuum. For a plane wave, Eq. (5) can be simplified to
$${{\textbf {H}}}({{\textbf {r}}}) = {{\textbf {B}}}({{\textbf {r}}}) = \frac{{{{\textbf {k}}} \times {\textbf {E}}({{\textbf {r}}})}}{{{k_0}}}.$$

In this study, we focus on the time-averaged Poynting vector of the electromagnetic field associated with a beam. The time-averaged Poynting vector is defined as

$$\bar{{{\textbf {S}}}} = \frac{c}{{8\pi }}{\mathop{\rm Re}\nolimits} ({{\textbf {E}} \times {{{\textbf {H}}}^\ast }} ),$$
where c is the speed of light. The time-averaged Poynting vector of a plane wave is obtained by combining Eqs. (6) and (7),
$$\bar{{{\textbf {S}}}} = \frac{c}{{8\pi {k_0}}}{\mathop{\rm Re}\nolimits} ({|{\textbf {E}}{|^2}} ){{\textbf {k}}}.$$

The direction of a Poynting vector is the energy propagation direction, and the Poynting vector magnitude represents the energy of the electromagnetic fields. Obviously, the energy of a plane wave propagates along k, which is also the beam propagation direction under the ray tracing scheme. Moreover, in the ray tracing scheme, the Poynting vectors across a beam are uniform.

An exact plane wave only exists in an infinite medium. In ray tracing, a beam originates from an area (or beam base) with a finite size. The electromagnetic field in a beam is induced by the field on the beam base. The geometric optics approximation assumes that the field on the beam base continues to propagate as a plane wave along the beam direction. To evaluate the geometric optics approximation, a rigorous method should be utilized to compute the electromagnetic fields in a beam induced by the field on the beam base.

Here we use the vector Kirchhoff integral to compute the induced electromagnetic field. The vector Kirchhoff integral is expressed as [24]

$${\textbf {E}}({{\textbf {r}}}) = \int_S {\{{i{k_0}[{{{\textbf {n}}^{\prime}} \times {{\textbf {B}}}({{\textbf {r}}^{\prime}})} ]G + [{{{\textbf {n}}^{\prime}} \times {\textbf {E}}({{\textbf {r}}^{\prime}})} ]\times \nabla^{\prime}G + [{{{\textbf {n}}^{\prime}} \cdot {\textbf {E}}({{\textbf {r}}^{\prime}})} ]\nabla^{\prime}G} \}d{{{{\textbf {r}}^{\prime}}}^2}} ,$$
in which S is the area where the beam emerges, and n′ is the unit vector normal to the area. n′ is pointing away from the region the beam propagates to. G is the Green function for an outgoing wave in an infinite space, and is expressed as
$$G({{\textbf {r}}},{{\textbf {r}}^{\prime}}) = \frac{{\exp ({ik|{{\textbf {r}}} - {{\textbf {r}}^{\prime}}|} )}}{{4\pi |{{\textbf {r}}} - {{\textbf {r}}^{\prime}}|}}.$$

The vector Kirchhoff integral is derived from Maxwell’s equations and the Green theorem, and is a rigorous formula to compute the electromagnetic field (both near and far fields) induced by the field on the surface integral area. We assume that the field on area S is a plane electromagnetic wave, which is identical to the geometric optics approximation, even though in a realistic light scattering event, the field on a small area of a particle may not be plane wave. With the plane wave assumption, we can focus on the properties of a single beam by excluding other geometric optics approximations. The induced electromagnetic near fields computed by Eqs. (5) and (9) are compared with the fields obtained with the ray tracing technique.

Because the induced electric field computed by Eq. (9) may not be a plane wave, Eq. (8) cannot be applied to calculate Poynting vectors. Equation (5) is used to derive the analytical expression of the induced magnetic field. We use the two-dimensional Simpson’s rule to numerically evaluate the surface integral in Eq. (9) and the derived magnetic field counterpart. The Poynting vector is computed by Eq. (7). Our numerical tests show that five integration points per wavelength distance are sufficient to get convergent results in the numerical integration. Square and circular shapes of the beam cross sections are considered in the calculations. We also perform the computations for various beam cross section sizes, beam propagation directions and medium refractive indices.

Without loss of generality, in the calculations, the electric field on the beam base area is set to

$${\textbf {E}}({{\textbf {r}}^{\prime}}) = \sqrt {\frac{{8\pi {k_0}}}{{mc}}} \textrm{exp}[{i({{{\textbf {k}}r^{\prime}} - \omega t} )} ]({1,\;\;0,\;\;0} ),$$
and the wavelength is set to 2π so the numerical value of the plane wave Poynting vector is unity in units of erg s-1 cm-2.

3. Results and analysis

In the following context, all lengths are expressed in terms of number of wavelengths, including the dimension of the beam cross section and the distance from the beam base. Cartesian coordinates are used to define the positions and directions. The coordinate origin is on the geometric center of the beam base, the z direction is opposite to the beam base normal vector, and the xoz plane is parallel to the paper.

We use the vector Kirchhoff integral to compute the electromagnetic near fields in and around the beams, and compute the Poynting vectors by Eq. (7). Figures 2 and 3 show the magnitudes and directions of the Poynting vectors on the xoz plane. Unlike the ray tracing approximations, the electromagnetic fields are not uniform in the beam along and perpendicular to propagation direction and the fields extend outside the beam boundary. The Poynting vectors have components both along and perpendicular to the beam propagation directions. The propagation-direction Poynting vector components represent the energy of the wave if it were transverse, while the perpendicular-direction components include the energy of the longitudinal wave. The fields that are farther from the beam base have smaller perpendicular-direction components. This implies that the field at an infinite distance from the beam base is transverse, which is consistent with the scattered far field [25]. It is also notable that the fields in the beam tend to become more uniform in the perpendicular direction, although their magnitudes keep decreasing along the progagation direction.

 figure: Fig. 2.

Fig. 2. The near field Poynting vectors on the xoz plane for beam base sizes d = 5, 10 and 15 times the wavelength. The medium refractive index is m = 1. The propagation direction of the field on the beam base is along the z axis. The contour plots show the Poynting vector magnitude $|\bar{{{\textbf {S}}}}|$, and the arrows show the Poynting vector directions. The lengths of the arrows are proportional to the square root of the corresponding Poynting vector magnitude. The red dashed lines indicate the distance from the beam base l = d2 times the wavelength. Upper plot: The beam base shape is a square with side width d; Lower plot: The beam base shape is a circle with diameter d.

Download Full Size | PDF

 figure: Fig. 3.

Fig. 3. The near field Poynting vectors on the xoz plane for beam base sizes d = 10 times the wavelength. The medium refractive index is m = 1. The propagation directions of the field on the beam base have angles θ = 0°, 30°, and 60° relative to the z axis. The contour plots show the Poynting vector magnitude $|\bar{{{\textbf {S}}}}|$, and the arrows show the Poynting vector directions. The lengths of the arrows are proportional to the square root of the corresponding Poynting vector magnitude. The red solid lines indicate the distance from the beam base l = d2 times the wavelength, and the red dashed lines indicate the distance from the beam base l′ = d2cosθ times the wavelength. Left plot: The beam base shape is a square with side width d; Right plot: The beam base shape is a circle with diameter d.

Download Full Size | PDF

The field patterns for a square and circle beam base are different, in particular in the region close to the beam base, where the field distributions on the beam cross section are complicated due to interference. The beams with the same beam base shape but various beam base sizes d have similar field patterns, but the scales of the patterns are different as shown in Fig. 2. At distance l = d2 times the wavelength, the Poynting vector magnitude has decreased substantially in the beams, but the magnitudes are similar for different beam base sizes.

Figure 3 shows the Poynting vectors for various beam directions. If a beam direction deviates from the normal direction of the beam base, the energy in the beam decays faster along the propagation direction than the beam along the normal direction of the beam base. The energy is also asymmetric on the beam cross sections. As shown in Fig. 3, at distance l = d2 times the wavelength, energy in the three beams with the same base shape and size but different propagation directions are clearly different. However, at distance l′ = d2cosθ, where θ is the angle between the propagation direction and z direction, the beams seem to have similar energy if they have the same base shape and size. $\sqrt {{d^2}\cos \theta } $ can be viewed as the geometric mean of the beam cross section.

Figures 4 and 5 show the Poynting vectors on various beam cross sections perpendicular to the propagation directions. The beams have base size d = 15 times the wavelength. The evolution of the field patterns along the beam propagation direction is clear in Figs. 4 and 5. The beams with square and circle bases have substantially different energy distributions on their cross sections at distances less than about 160 times the wavelength from the beam bases. However, at distance farther than about 160λ, their field patterns look similar. The field distribution on the beam cross section tends to become more and more uniform, which is consistent with the phenomenon observed in the xoz plane (Figs. 2 and 3). At distance l = d2 from the beam bases (z′ = 225.0λ plots in Figs. 4 and 5), the field distributions for square and circle beam bases resemble each other.

 figure: Fig. 4.

Fig. 4. The near field Poynting vectors on the beam cross section planes at various distances from the beam base. The cross section planes are perpendicular to the propagation directions of the field on the beam base. The beam base is a square, and has side width d = 15 times the wavelength. The medium refractive index is m = 1. The propagation direction of the field on the beam base is along the z axis. The primed x, y and z are the Cartesian coordinates with z′ along the beam propagation direction. The contour plots show the Poynting vector magnitude $|\bar{{{\textbf {S}}}}|$, and the arrows show the Poynting vector directions. The lengths of the arrows are proportional to the square root of the corresponding Poynting vector magnitude.

Download Full Size | PDF

 figure: Fig. 5.

Fig. 5. Same as Fig. 4, except the beam base shape is a circle with diameter d = 15 times the wavelength.

Download Full Size | PDF

We also consider the energy flux crossing the beam cross sections,

$$P = \int\limits_{S^{\prime}} {\bar{{{\textbf {S}}}}({{\textbf {r}}}) \cdot \tilde{k}{d^2}{{\textbf {r}}}} ,$$
where $\tilde{k}$ is a unit vector along the beam propagation direction, and S′ is the cross section region. In the ray tracing scheme, Eq. (12) is equal to the area of the cross section multiplied by the Poynting vector component along $\tilde{k}$, and is conservative along the propagation direction if the medium is non-absorptive. We use the vector Kirchhoff integral calculation results to compute the accurate beam energy flux Eq. (12).

Figure 6 shows the beam energy flux computed by the vector Kirchhoff integral relative to the geometric optics ray tracing counterpart, which is equal to the energy flux at the beam base. The beam energy flux decreases along the propagation directions. For the same beam base sizes, the beam energy flux decreases faster if the beam direction deviates from the beam base normal direction. At distance l = d2cosθ, all of the beam energy fluxes decay more than 80%. With the same beam direction, the beam energy fluxes decay to almost the same fraction at distance l = d2cosθ times the wavelength from the beam bases for various beam base sizes. This is consistent with Fig. 3, where at distance l = d2cosθ, the Poynting vector magnitudes are similar for different beam base sizes.

 figure: Fig. 6.

Fig. 6. The ratios between the energy fluxes crossing the beam cross section computed by the vector Kirchhoff integral (PK) and ray tracing (PR). The beam bases are circles, and have diameters d = 5, 10 and 15 times the wavelength. The propagation directions of the field on the beam base have angles θ = 0° and 60° relative to the z axis. The medium refractive index is m = 1. The black dots on the curves indicate the distance from the beam base l = d2cosθ times the wavelength.

Download Full Size | PDF

We also compute the energy fluxes for beams in media with various refractive indices. The beam bases are square with side length d = 10 times the wavelength in a vacuum, and the propagation directions are along z. The results are shown in Fig. 7. At distance l = md2, the energy fluxes decay to the same fraction for various medium refractive indices. These results suggest that a beam can propagate further in a medium with larger refractive index. Note that 0 is equal to mdλ, where λ is the wavelength in the medium λ0/m, and md2λ0 is equal to (md)2λ. Thus, the phenomenon is still consistent with van de Hulst’s statement, and the wavelength in that statement should be the wavelength in the medium rather than in a vacuum.

 figure: Fig. 7.

Fig. 7. The ratios between the power crossing the beam cross section computed by the vector Kirchhoff integral (PK) and ray tracing (PR). The beam bases are square, and have side width d = 10 times the vacuum wavelength. The propagation direction of the field on the beam base is along the z axis. Four medium refractive indices m (=1.0, 1.3, 1.5, and 2.0) are considered. The black dots on the curves indicate the distance from the beam base l = md2 times the wavelength in a vacuum.

Download Full Size | PDF

We find that the accuracy of applying physical-geometric optics method (PGOM) [20] to compute particle single-scattering properties is dependent on the particle refractive index. The PGOM utilizes the geometric ray tracing technique to obtain the near field in or around the particle, and applies electromagnetic volume (PGOMV) or surface (PGOMS) integrals to obtain scattering fields in the far-field limit. Figure 8 shows two examples of PGOMS and PGOMV computational results of scattering phase functions compared with the numerically accurate invariant imbedding T-matrix (II-TM) method [13]. The scattering phase functions in Fig. 8 are for the same particle shape. If the particle refractive index is 1.308, the PGOMS and PGOMV are both almost consistent with the II-TM result. If the refractive index is 1.05, PGOMS performs worse than PGOMV comparing with the II-TM result. The PGOMS is actually always less accurate if the real part of the refractive index of the particle is only slightly higher than 1 according to our previous experiments.

 figure: Fig. 8.

Fig. 8. The scattering phase functions of a hexagonal column with aspect ratio (2a/H) 1 and size parameter (ka) 100 computed by the II-TM, PGOMS, and PGOMV methods. Left: refractive index m=1.308; Right: refractive index m=1.05.

Download Full Size | PDF

It has been proved that the electromagnetic volume and surface integrals are equivalent in computing the scattering field [22] so the difference between PGOMS and PGOMV should be attributed to the ray tracing-computed near field. The PGOMS uses the electromagnetic surface integral, which computes the far field using the near fields only on the particle external surface. According to the current study, the fields on the particle external surface computed by ray tracing are less accurate if the particle has a smaller refractive index, which contributes to the inaccuracy of the far field. In comparison, the PGOMV uses the electromagnetic volume integral, which uses all fields within the particle to compute the far fields so the errors caused by ray tracing may be reduced. Thus, the PGOMV is more accurate than PGOMS when the real part of the particle refractive index is close to 1.

As shown in Figs. 4 and 5, the field distribution in a beam becomes more and more uniform after the beam propagates away from its base. We examine the Poynting vector component along the propagation direction ($\bar{{{\textbf {S}}}} \cdot \tilde{k}$) on the beam cross section at l = d2cosθ. $\bar{{{\textbf {S}}}} \cdot \tilde{k}$ at y = 0 is shown in Fig. 9. Beams with circle and square bases both have bell shape distributions of $\bar{{{\textbf {S}}}} \cdot \tilde{k}$. The Fraunhofer (far-field) diffraction through a square aperture generates an energy distribution with the form of [sin(η)/η]4, in which η is proportional to the diffracted field position. If the aperture is a circle, the Fraunhofer diffraction generates an energy distribution with the form of [J1(η)/η]2, where J1 is the first order Bessel function of the first kind [23]. We fit $\bar{{{\textbf {S}}}} \cdot \tilde{k}$ numerical results by

$$\bar{{{\textbf {S}}}} \cdot \tilde{k} = A{\left( {\frac{{\sin \kappa x}}{{\kappa x}}} \right)^4} + B,$$
for the square beam base, and
$$\bar{{{\textbf {S}}}} \cdot \tilde{k} = A{\left[ {\frac{{{J_1}(\kappa x)}}{{\kappa x}}} \right]^2} + B,$$
for the circle beam base, where A, B and κ are fitting parameters to be determined. The fitting results are almost perfectly consistent with the numerical results, even though distance l = d2cosθ may not be far enough for the Fraunhofer diffraction approximation.

 figure: Fig. 9.

Fig. 9. The Poynting vector component along the propagation direction of the field on the beam base at the distance l=225 times the wavelength from the beam base, and y = 0. The beam base size is 15 times the wavelength. The medium refractive index is m = 1. The numerical results for the square and circle beam bases are fitted by two analytical equations.

Download Full Size | PDF

In Fig. 10, we show the computed energy distributions for beam bases with size d = 5 at distance 1000 times wavelength from the bases, where far field diffraction approximation applies. The energy distribution at 1000 times wavelength for the square beam base is consistent with the Fraunhofer diffraction patterns for a square aperture. The maxima and minima of energy distribution are obvious. For the circular beam base, the energy distribution at 1000 times wavelength is an Airy disk, which is Fraunhofer diffraction pattern of a circular aperture. The consistencies with Fraunhofer diffraction patterns suggest that the numerical results converge to Fraunhofer diffraction in the far field, where the beam as defined in the geometric optics method definitely does not exist.

 figure: Fig. 10.

Fig. 10. The far field Poynting vector magnitude $|\bar{{{\textbf {S}}}}|$ on the beam cross section planes at 1000 times wavelength from beam base. For viewing purpose, the contour plots are the quarter root of $|\bar{{{\textbf {S}}}}|$. The cross section planes are perpendicular to the propagation directions of the field on the beam base. The medium refractive index is m = 1. The primed x, y and z are the Cartesian coordinates with z′ along the beam propagation direction. The propagation direction of the field on the beam base is along the z axis. Left: The beam base is a square, and has side width d = 5 times the wavelength; Right: The beam base shape is a circle with diameter d = 15 times the wavelength.

Download Full Size | PDF

4. Summary and conclusion

We use the rigorous vector Kirchhoff integral to numerically compute the electromagnetic near field in a beam (i.e. a pencil of light, a bundle of rays) defined in the geometric optics ray tracing scheme. The numerical results are used to examine a statement made by van de Hulst about the existence of a beam. The statement is that, if a beam has a cross section with the size of the order , the beam can have independent existence over a length with the order of p2λ. The statement is shown to be a corollary of the Huygens–Fresnel principle. The accurate computational results show that the fields within a beam are not uniform, and have complicated structures in the region close to the beam base. The fields in the beams are different for various beam base shapes and propagation directions. At a distance l = d2cosθ times the wavelength, beam fields with different base sizes show similar patterns and energy fluxes. Based on these findings, the statement about beam existence can be made more specific: if λ is the wavelength in the beam propagation medium, a beam having the geometric mean width of the order can propagate independently to a distance p2λ away from the beam base. The electromagnetic fields in the beam tend to be uniform, but the energy flux diverges during its propagation. The beam loses most of its energy at p2λ. For example, a beam with a circular cross section loses more than 80% of its energy at the distance p2λ. The fractions of energy flux loss at p2λ are similar for different p values. At p2λ, the energy distribution in the beam is virtually the same as the Fraunhofer diffraction pattern.

Understanding the existence of a beam and the electromagnetic fields in the beam is required to evaluate the performance of the geometric optics ray tracing technique in computing particle single-scattering properties. The beam energy flux decays more slowly during propagation if the beam has larger base size. In geometric optics ray tracing approximation, the beam energy flux is assumed to be conserved during propagation in a non-absorbing medium. Thus, if the exact beam energy flux decays more slowly during propagation, the geometric optics approximation tends to be more accurate. Using geometric optics ray tracing to compute light scattering by a particle will produce more accurate results if the beam base size is chosen to be larger, and it is required that the particle must be much larger than the incident wavelength. The PGOM ray tracing must be applied with the beam as large as the particle dimension, so the geometric optics ray approximation error is minimized.

In this study, the existence of any single beam is investigated independently of any other beams. Interference among different beams is not considered in the study. Neglecting the interference among beams in geometric optics ray tracing certainly introduces errors. However, it is complicated to account for the interference among beams during ray tracing. Instead, according to the current study, if a beam has a larger base size, less energy will “leak out” of the beam so the interference among beams has less effect on the calculation results. Thus, the error of neglecting the interference among beams is also reduced if we assume a beam with a larger base size.

Funding

Texas A & M University T3 Grant (246280-00000); National Science Foundation (AGS-1826936); National Aeronautics and Space Administration Radiation Sciences Program.

Acknowledgement

The computations were conducted at the Texas A&M University Supercomputing Facility.

Disclosures

The authors declare no conflicts of interest.

References

1. H. C. van de Hulst, Light Scattering by Small Particles (John Wiley and Sons, 1957).

2. P. Wendling, R. Wendling, and H. K. Weickmann, “Scattering of solar radiation by hexagonal ice crystals,” Appl. Opt. 18(15), 2663 (1979). [CrossRef]  

3. Q. Cai and K. N. Liou, “Polarized light scattering by hexagonal ice crystals: theory,” Appl. Opt. 21(19), 3569–3580 (1982). [CrossRef]  

4. A. Macke, “Scattering of light by polyhedral ice crystals,” Appl. Opt. 32(15), 2780–2788 (1993). [CrossRef]  

5. J. A. Lock, “Ray scattering by an arbitrarily oriented spheroid. I. Diffraction and specular reflection,” Appl. Opt. 35(3), 500 (1996). [CrossRef]  

6. J. A. Lock, “Ray scattering by an arbitrarily oriented spheroid. II. Transmission and cross-polarization effects,” Appl. Opt. 35(3), 515–531 (1996). [CrossRef]  

7. P. Yang and K. N. Liou, “Geometric-optics–integral-equation method for light scattering by nonspherical ice crystals,” Appl. Opt. 35(33), 6568 (1996). [CrossRef]  

8. A. G. Borovoi and I. A. Grishin, “Scattering matrices for large ice crystal particles,” J. Opt. Soc. Am. A 20(11), 2071 (2003). [CrossRef]  

9. M. I. Mishchenko, “Light scattering by randomly oriented axially symmetric particles,” J. Opt. Soc. Am. A 8(6), 871 (1991). [CrossRef]  

10. P. Yang and K. N. Liou, “Finite-difference time domain method for light scattering by small ice crystals in three-dimensional space,” J. Opt. Soc. Am. A 13(10), 2072 (1996). [CrossRef]  

11. M. A. Yurkin and A. G. Hoekstra, “The discrete dipole approximation: An overview and recent developments,” J. Quant. Spectrosc. Radiat. Transf. 106(1-3), 558–589 (2007). [CrossRef]  

12. C. Liu, R. Lee Panetta, and P. Yang, “Application of the pseudo-spectral time domain method to compute particle single-scattering properties for size parameters up to 200,” J. Quant. Spectrosc. Radiat. Transf. 113(13), 1728–1740 (2012). [CrossRef]  

13. L. Bi and P. Yang, “Accurate simulation of the optical properties of atmospheric ice crystals with the invariant imbedding T-matrix method,” J. Quant. Spectrosc. Radiat. Transf. 138, 17–35 (2014). [CrossRef]  

14. A. Macke, M. I. Mishchenko, K. Muinonen, and B. E. Carlson, “Scattering of light by large nonspherical particles: ray-tracing approximation versus T-matrix method,” Opt. Lett. 20(19), 1934 (1995). [CrossRef]  

15. P. Yang and K. N. Liou, “Light scattering by hexagonal ice crystals: comparison of finite-difference time domain and geometric optics models,” J. Opt. Soc. Am. A 12(1), 162 (1995). [CrossRef]  

16. K. Muinonen, “Scattering of light by crystals: a modified Kirchhoff approximation,” Appl. Opt. 28(15), 3044–3050 (1989). [CrossRef]  

17. P. Yang and K. N. Liou, “Light scattering by hexagonal ice crystals: solutions by a ray-by-ray integration algorithm,” J. Opt. Soc. Am. A 14(9), 2278 (1997). [CrossRef]  

18. L. Bi, P. Yang, G. W. Kattawar, Y. Hu, and B. A. Baum, “Scattering and absorption of light by ice particles: Solution by a new physical-geometric optics hybrid method,” J. Quant. Spectrosc. Radiat. Transf. 112(9), 1492–1508 (2011). [CrossRef]  

19. A. Borovoi, A. Konoshonkin, and N. Kustova, “The physical-optics approximation and its application to light backscattering by hexagonal ice crystals,” J. Quant. Spectrosc. Radiat. Transf. 146, 181–189 (2014). [CrossRef]  

20. B. Sun, P. Yang, G. W. Kattawar, and X. Zhang, “Physical-geometric optics method for large size faceted particles,” Opt. Express 25(20), 24044 (2017). [CrossRef]  

21. A. Konoshonkin, A. Borovoi, N. Kustova, H. Okamoto, H. Ishimoto, Y. Grynko, and J. Förstner, “Light scattering by ice crystals of cirrus clouds: From exact numerical methods to physical-optics approximation,” J. Quant. Spectrosc. Radiat. Transf. 195, 132–140 (2017). [CrossRef]  

22. P. Yang, J. Ding, R. L. Panetta, K.-N. Liou, G. W. Kattawar, and M. I. Mishchenko, “On the convergence of numerical computations for both exact and approximate solutions for electromagnetic scattering by nonspherical dielectric particles (invited review),” Prog. Electromagn. Res. 164, 27–61 (2019). [CrossRef]  

23. M. Born, E. Wolf, A. B. Bhatia, P. C. Clemmow, D. Gabor, A. R. Stokes, A. M. Taylor, P. A. Wayman, and W. L. Wilcock, Principles of Optics (Cambridge University Press, 1999).

24. J. D. Jackson, Classical Electrodynamics, 2nd ed. (John Wiley & Sons, Inc., 1975).

25. M. Mischenko, L. Travis, and A. Lacis, Scattering, Absorption, and Emission of Light by Small Particles (Cambridge University Press, 2002).

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (10)

Fig. 1.
Fig. 1. Illustration of Huygens–Fresnel principle.
Fig. 2.
Fig. 2. The near field Poynting vectors on the xoz plane for beam base sizes d = 5, 10 and 15 times the wavelength. The medium refractive index is m = 1. The propagation direction of the field on the beam base is along the z axis. The contour plots show the Poynting vector magnitude $|\bar{{{\textbf {S}}}}|$ , and the arrows show the Poynting vector directions. The lengths of the arrows are proportional to the square root of the corresponding Poynting vector magnitude. The red dashed lines indicate the distance from the beam base l = d2 times the wavelength. Upper plot: The beam base shape is a square with side width d; Lower plot: The beam base shape is a circle with diameter d.
Fig. 3.
Fig. 3. The near field Poynting vectors on the xoz plane for beam base sizes d = 10 times the wavelength. The medium refractive index is m = 1. The propagation directions of the field on the beam base have angles θ = 0°, 30°, and 60° relative to the z axis. The contour plots show the Poynting vector magnitude $|\bar{{{\textbf {S}}}}|$ , and the arrows show the Poynting vector directions. The lengths of the arrows are proportional to the square root of the corresponding Poynting vector magnitude. The red solid lines indicate the distance from the beam base l = d2 times the wavelength, and the red dashed lines indicate the distance from the beam base l′ = d2cosθ times the wavelength. Left plot: The beam base shape is a square with side width d; Right plot: The beam base shape is a circle with diameter d.
Fig. 4.
Fig. 4. The near field Poynting vectors on the beam cross section planes at various distances from the beam base. The cross section planes are perpendicular to the propagation directions of the field on the beam base. The beam base is a square, and has side width d = 15 times the wavelength. The medium refractive index is m = 1. The propagation direction of the field on the beam base is along the z axis. The primed x, y and z are the Cartesian coordinates with z′ along the beam propagation direction. The contour plots show the Poynting vector magnitude $|\bar{{{\textbf {S}}}}|$ , and the arrows show the Poynting vector directions. The lengths of the arrows are proportional to the square root of the corresponding Poynting vector magnitude.
Fig. 5.
Fig. 5. Same as Fig. 4, except the beam base shape is a circle with diameter d = 15 times the wavelength.
Fig. 6.
Fig. 6. The ratios between the energy fluxes crossing the beam cross section computed by the vector Kirchhoff integral (PK) and ray tracing (PR). The beam bases are circles, and have diameters d = 5, 10 and 15 times the wavelength. The propagation directions of the field on the beam base have angles θ = 0° and 60° relative to the z axis. The medium refractive index is m = 1. The black dots on the curves indicate the distance from the beam base l = d2cosθ times the wavelength.
Fig. 7.
Fig. 7. The ratios between the power crossing the beam cross section computed by the vector Kirchhoff integral (PK) and ray tracing (PR). The beam bases are square, and have side width d = 10 times the vacuum wavelength. The propagation direction of the field on the beam base is along the z axis. Four medium refractive indices m (=1.0, 1.3, 1.5, and 2.0) are considered. The black dots on the curves indicate the distance from the beam base l = md2 times the wavelength in a vacuum.
Fig. 8.
Fig. 8. The scattering phase functions of a hexagonal column with aspect ratio (2a/H) 1 and size parameter (ka) 100 computed by the II-TM, PGOMS, and PGOMV methods. Left: refractive index m=1.308; Right: refractive index m=1.05.
Fig. 9.
Fig. 9. The Poynting vector component along the propagation direction of the field on the beam base at the distance l=225 times the wavelength from the beam base, and y = 0. The beam base size is 15 times the wavelength. The medium refractive index is m = 1. The numerical results for the square and circle beam bases are fitted by two analytical equations.
Fig. 10.
Fig. 10. The far field Poynting vector magnitude $|\bar{{{\textbf {S}}}}|$ on the beam cross section planes at 1000 times wavelength from beam base. For viewing purpose, the contour plots are the quarter root of $|\bar{{{\textbf {S}}}}|$ . The cross section planes are perpendicular to the propagation directions of the field on the beam base. The medium refractive index is m = 1. The primed x, y and z are the Cartesian coordinates with z′ along the beam propagation direction. The propagation direction of the field on the beam base is along the z axis. Left: The beam base is a square, and has side width d = 5 times the wavelength; Right: The beam base shape is a circle with diameter d = 15 times the wavelength.

Equations (14)

Equations on this page are rendered with MathJax. Learn more.

u Q 1 = P u P i exp ( i k r ) r λ K ( χ ) d S ,
u Q = exp ( i k l ) u P ,
u Q 1 = P u P i exp ( i k l ) l λ d S .
E ( r ) = E 0 exp [ i ( k r ω t ) ] ,
H ( r ) = B ( r ) = × E ( r ) i k 0 ,
H ( r ) = B ( r ) = k × E ( r ) k 0 .
S ¯ = c 8 π Re ( E × H ) ,
S ¯ = c 8 π k 0 Re ( | E | 2 ) k .
E ( r ) = S { i k 0 [ n × B ( r ) ] G + [ n × E ( r ) ] × G + [ n E ( r ) ] G } d r 2 ,
G ( r , r ) = exp ( i k | r r | ) 4 π | r r | .
E ( r ) = 8 π k 0 m c exp [ i ( k r ω t ) ] ( 1 , 0 , 0 ) ,
P = S S ¯ ( r ) k ~ d 2 r ,
S ¯ k ~ = A ( sin κ x κ x ) 4 + B ,
S ¯ k ~ = A [ J 1 ( κ x ) κ x ] 2 + B ,
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.