Intensity, Wavefront Curvature, and Polarization
The following describes the algorithm used to compute the ray intensity and polarization when Intensity computation is set to Compute intensity or Compute intensity and power in the physics interface Intensity Computation section. In this algorithm, each ray is treated as a propagating wavefront subtending a small solid angle. The algorithm is only valid when the rays propagate in homogeneous media. The following auxiliary dependent variables are defined for each ray:
In 3D:
The initial ray intensity I0 and normalized Stokes parameters sn1, sn2, and sn3. By allocating four degrees of freedom in this way it is possible to characterize any intensity magnitude and polarization state.
The integral A of the attenuation coefficient along the ray path.
Two principal radii of curvature, r1 and r2, which represent the maximum and minimum radii of curvature of the intersection of the wavefront with an arbitrary plane.
Three components of a unit vector e1 in the direction corresponding to the first principal radius of curvature. This information is used to reinitialize the principal radii of curvature at curved boundaries.
In 2D:
The initial ray intensity I0 and normalized Stokes parameters sn1, sn2, and sn3.
The integral A of the attenuation coefficient along the ray path.
One radius of curvature r1 of the wavefront. All wavefronts are assumed to be converging or diverging cylindrical waves, so it is not necessary to define a second radius of curvature.
A different algorithm is used to compute the intensity and polarization of each ray if the Intensity computation is instead set to Compute intensity in graded media or Compute intensity and power in graded media. This algorithm accounts for the effect of graded media on the ray intensity, but it is also slower and less accurate for homogeneous media. The following auxiliary dependent variables are defined for each ray:
In 3D:
The intensity help variable Γ and the normalized Stokes parameters sn01, sn02, and sn03, which characterize the intensity and polarization of the ray. They are reinitialized at material discontinuities and walls.
The integral A of the attenuation coefficient along the ray path.
Two principal curvature calculation help variables α1 and α2 and the rotation angle φ which indicates the orientation of the principal curvature directions. An additional help variable is used internally to detect poles in the local coordinate system definition and to redefine the local coordinate system accordingly.
In 2D:
The intensity help variable Γ and the normalized Stokes parameters sn01, sn02, and sn03, which characterize the intensity and polarization of the ray. They are reinitialized at material discontinuities and walls.
The integral A of the attenuation coefficient along the ray path.
For the options Compute Power, Compute intensity and power, and Compute intensity and power in graded media an additional auxiliary dependent variable is defined to indicate the total power transmitted by each ray. Like intensity, power is affected by absorbing media and by reflection or refraction at boundaries. Unlike intensity, power does not change due to the focusing or divergence of a ray bundle.
The Stokes Parameters
The four Stokes parameters completely characterize the intensity and polarization of a fully polarized, partially polarized, or unpolarized ray. The Stokes parameters can be interpreted as indicators of the ray intensity that would be measured when sending a ray through various arrangements of polarizers and wave retarders.
The Stokes parameters of a ray are defined as in Ref. 9 as follows. Consider a ray propagating along the z-axis of a Cartesian coordinate system, with electric field components polarized within the xy-plane. In 3D models the x-axis is parallel to the first principal curvature direction e1. In 2D models the x-axis is always oriented in the out-of-plane direction.
Let I(θ,ε) be the intensity of radiation polarized in the xy-plane at an angle θ counterclockwise from the x-axis, when a phase angle of ε is introduced between the x- and y-components. For example, I(0,0) is the intensity that would be measured after sending a ray through a linear polarizer with a transmission axis parallel to the x-axis, and I(π/4, π/2) is the intensity that would be measured after sending a ray through a system of polarizers and wave retarders that only transmits radiation with right-handed circular polarization. The four Stokes parameters are then defined as follows:
is the sum of the intensity of linearly polarized light in the x and y directions. It represents the total intensity of the ray.
is the difference in intensity between light that is linearly polarized in the x direction and light that is linearly polarized in the y direction.
is the difference in intensity between linearly polarized light in the direction of the line y = x and linearly polarized light in the direction of the line y = −x.
For example, a ray with linear polarization parallel to the x-axis has Stokes parameters s0 = s1 = I and s2 = s3 = 0. A fully polarized ray is characterized by the relation
A partially polarized or unpolarized ray is characterized by the relation
The degree of polarization P is defined as
A degree of polarization P = 1 corresponds to fully polarized radiation, whereas a degree of polarization P = 0 corresponds to unpolarized radiation.
In the COMSOL implementation of ray intensity calculation, the degrees of freedom are not the Stokes parameters themselves, but rather normalized Stokes parameters sn1, sn2, and sn3. Then the actual Stokes parameters are recovered through the relations s1 = Isn1, s2 = Isn2, and s3 = Isn3. The reason for this decoupling is that the normalized Stokes parameters remain bounded, enabling accurate reinitialization of ray power, even if a ray is reflected or refracted very close to a caustic, where the intensity (in the geometrical optics limit) can become arbitrarily large.
Principal Radii of Curvature
In 3D models, each ray is treated as a wavefront for which two principal radii of curvature, r1 and r2, are defined. In addition, the values of the two principal radii of curvature are stored as r1,i and r2,i whenever the ray reaches a boundary.
Within domains, the equations
are solved. Negative radii of curvature indicate that the wavefront is expanding as the ray propagates, while positive radii of curvature indicate that the wavefront is converging. A continuous locus of points at which either of the principal radii of curvature equals zero is called a caustic surface. The unit vector e1 is defined so that r1 is the radius of curvature of the intersection of the wavefront with the plane tangent to e1 and the wave vector k. Similarly, if a vector e2 is defined so that
Then r2 is the radius of curvature of the intersection of the wavefront with a plane tangent to e2 and k. The components of e1 are stored as auxiliary dependent variables for each ray. The components of e2 can then be derived from e1 and k at any time.
The principal radii of curvature are reinitialized at material discontinuities, and the orientation of the unit vectors e1 and e2 may change. If the unit vector in the direction of the incident ray ni is not parallel to the surface normal ns, then it is possible to define a unique tangent plane, called the plane of incidence, that contains the incident, reflected, and refracted rays. The unit vector normal to this plane, u0, is defined as
If the ray is normal to the surface, then the incident ray is parallel to the refracted ray and antiparallel to the reflected ray, and u0 can be any arbitrary unit vector orthogonal to ni. In addition to the unit vector normal to the plane of incidence, the following unit vectors tangent to the incident wavefront, refracted wavefront, reflected wavefront, and surface normal are defined:
where the subscripts i, t, and r denote the incident, refracted, and reflected wavefronts, respectively. For a wavefront propagating in a direction n, with principal curvatures k1 = 1/r1 and k2 = 1/r2 defined for directions e1 and e2, the principal curvatures in two other orthogonal directions e1' and e2' (both orthogonal to n) are
(3-19)
where θ is the angle of rotation about n which transforms the vectors e1 and e2 to e1' and e2', respectively. Because k1' and k2' are not principal curvatures, it is necessary to include off-diagonal elements of the wavefront curvature tensor equal to k12'.
The following algorithm is used to reinitialize the principal radii of curvature of the wavefront and their orientations. The reinitialization of the curvature variables follows the method of Stavroudis in Ref. 11.
1
Given ni and ns, compute unit vectors in the directions of the reflected and refracted rays, nr and nt.
2
Compute the vectors u0, ui, ut, ur, and us.
3
Compute the angle of rotation θ(i) needed to transform the local coordinate system with axes parallel to e1, e2, and ni to a local coordinate system with axes parallel to u0, ui, and ni. (Note: the superscript (i) is used to distinguish this coordinate system rotation angle from the angle of incidence of the ray, often denoted θi.)
4
5
Given the two principal curvatures of the surface, k1,s and k2,s with directions e1,s and e2,s, transform the curvature variables to a local coordinate system with axes parallel to u0, us, and ns. Let the new curvature variables be denoted by k1,s', k2,s', and k12,s'.
6
The curvatures of the wavefront and the surface have now been defined in coordinate systems that share the axis u0 and only differ by a rotation by the angle of incidence θi about u0. Defining the variables η and γ as in Material Discontinuity Theory, compute the curvature variables of the refracted ray in a coordinate system defined by u0, ut, and nt using the equations
7
Obtain the principal curvatures of the refracted ray by rotating the coordinate system defined by u0, ut, and nt by an angle θ(t) about nt. The angle θ(t) is defined as
8
9
Invert the principal curvatures to obtain the principal radii of curvature of the refracted ray. Rotate u0 about nt by the angle θ(t) to obtain the reinitialized principal curvature direction e1,t.
10
11
Obtain the principal curvatures of the reflected ray by rotating the coordinate system defined by u0, ur, and nr by an angle θ(r) about nr. The angle θ(r) is defined as
12
13
Invert the principal curvatures to obtain the principal radii of curvature of the reflected ray. Rotate u0 about nr by the angle θ(r) to obtain the initialized principal curvature direction e1,r.
In 2D components, only one principal radius of curvature is computed; each ray is treated as a cylindrical wave. In addition, no auxiliary degrees of freedom for the principal curvature direction are required because the out-of-plane direction can always be treated as one of the axes of the local coordinate system that defines the orientation of the wavefront.
In 2D axisymmetric model components, radii of curvature are computed for the in-plane direction, that is the rz-plane, as well as the out-of-plane or azimuthal directions. During ray-boundary interactions, these radii of curvature are reinitialized as if the ray interacted with a 3D surface of revolution, which may have finite radii of curvature in both the in-plane and out-of-plane directions.
Stokes Vector Calculation
The values of the intensity I and normalized Stokes parameters are stored as the auxiliary dependent variables I0, sn1, sn2, and sn3 when a ray is released. These auxiliary variables are updated when the ray hits a boundary. At any point along the ray’s trajectory, each Stokes parameter is equal to
In 2D, the second principal radius of curvature r2 is treated as an arbitrarily large distance that remains constant for each ray. Because s1 and s2 are based on the differences in intensity between orthogonal polarizations of radiation, a local coordinate system must be defined for the ray. Because auxiliary degrees of freedom have already been allocated for e1 as explained in Principal Radii of Curvature, a local coordinate system is defined with axes parallel to e1, e2, and k. The axes parallel to e1 and e2 then function as the x- and y-axes in The Stokes Parameters, respectively.
Stokes Vector Reinitialization
Reflection and Refraction of S- and P-polarized Rays
Recall from The Fresnel Equations that the ratios of the electric field amplitudes of the incident, reflected, and refracted rays can be expressed as a set of Fresnel coefficients,
where the subscripts s and p denote s- and p-polarized rays, or rays with linear polarizations perpendicular to and parallel to the plane of incidence, respectively. The incident ray is assumed to move from a region of refractive index n1 toward a region of refractive index n2. The angles θi and θt are the angle of incidence and angle of refraction, respectively.
If n1 > n2 and the angle of incidence exceeds the critical angle θ = asin(n2/n1), the incident ray undergoes total internal reflection, and the reflected ray has intensity equal to that of the incident ray.
The intensity of the transmitted and reflected waves are related to the intensity of the incident wave Ii by the equations
for p-polarized rays, and
for s-polarized rays. The coefficients being multiplied by the incident ray intensity are called the reflectance R and transmittance T:
Phase Shift Calculation for Reflected and Refracted Rays
Any fully polarized ray can be resolved into a pair of s- and p-polarized rays, with a phase shift δ between them. The reflected s- and p-polarized rays at a material discontinuity may then undergo different phase shifts; this may, for example, cause a linearly polarized incident ray to be yield a reflected ray with elliptical polarization. In addition, a phase shift may be applied to the transmitted ray if one or more Thin Dielectric Film subnodes are added to a Material Discontinuity node. The phase shifts of s- and p-polarized rays are computed from the complex-valued Fresnel coefficients:
Unpolarized and Partially Polarized Rays
An unpolarized ray has the following properties:
An unpolarized ray has Stokes parameters s1 = s2 = s3 = 0.
A partially polarized ray is neither completely deterministic nor completely random. It is characterized by a degree of polarization between 0 and 1.
Often, the ray being reflected or refracted at a boundary is unpolarized or partially polarized; or the polarization direction neither lies in the plane of incidence nor perpendicular to it. In this case, a more general formulation for the new values of the Stokes parameters is needed. The Stokes parameters of the reflected and refracted rays can be expressed in terms of the Stokes parameters of the incident ray, the transmittance or reflectance at the boundary, and the phase delay introduced by the boundary (if any) (Ref. 15):
where the subscripts r and t represent reflection and refraction, respectively; and the subscripts s and p indicate the transmittance or reflectance for s- and p-polarized rays.