3 Path Solver

Refer to caption — Figure 7: The path solver constructs propagation paths that connect a source and a target in a scene.

A path solver aims to determine a set of paths that connect two endpoints in a scene, as illustrated in Figure 7. From this set of paths, a time- or frequency-domain CIR can be computed. The path solver leverages SBR to efficiently compute paths in a scene. SBR operates as follows: First, rays²² 2 The processing of the rays is parallelized for accelerated computation. are traced from the source into directions specified by a lattice on the unit sphere (see Section 3.1.1). Then, an SBR loop iteratively tests the intersection of the rays with the scene geometry. At every intersection, an interaction type is selected from ${\cal I}$ , and a new ray is spawned according to the selected interaction type. The method for choosing interaction types is detailed in Section 3.1.2. The loop terminates when no more rays remain active, meaning they have either bounced out of the scene, reached a predefined maximum number of bounces, or satisfied another termination criterion. Notably, at each interaction of a ray with a surface in the scene, only one scattered ray is spawned, as spawning multiple rays would result in an exponential increase in the number of rays, quickly overwhelming the available computational resources.

When a diffuse reflection occurs during the SBR loop, the solver checks whether there is a line-of-sight (LoS) to the target, meaning that the ray segment from the interaction point to the target point is not occluded. If such a connection exists, the path is deemed valid and is returned by the solver. This process is known as next-event estimation. Subsequently, a new ray is spawned in a random direction within the hemisphere defined by the surface normal of the primitive involved. This new ray may lead to the discovery of additional paths. Next-event estimation is only possible for diffuse reflections that scatter energy into all directions. For specular reflection and refraction, only a single ray is spawned in the only possible direction. Consequently, the probability of finding paths connecting the source to the target that end with a specular reflection or refraction by SBR is effectively zero, as depicted in Figure 8(a). This is because hitting a specific point (the target) with a specular reflection or refraction requires perfect alignment of the incident ray, which has zero probability in continuous space. Likewise, for diffraction, the probability that the target is on the surface of the Keller cone is zero, as shown in Figure 8(b). When diffraction occurs, a single ray is spawned in a random direction along the Keller cone to continue the sample. Paths ending with a specular suffix can therefore not be determined solely by SBR. In particular, specular chains cannot be computed by SBR. Since specular chains typically carry a significant amount of the transported energy, SBR alone is insufficient for computing CIRs.

To compute paths ending with a specular suffix, including specular chains, the path solver uses an image method-based algorithm. However, a brute-force implementation of the image method that would test all possible combinations of primitives and wedges that constitute the scene is prohibitive even for moderately large scenes. Therefore, the path solver uses SBR as a heuristic to find path candidates. A path candidate consists of a sequence of primitives or wedges and interaction types ending with a specular reflection, refraction, or diffraction. As the path vertices determined during the SBR do not form a path that connects to the target, the image method is used to adjust the vertex locations to establish the connection. Intuitively, SBR pre-selects path candidates in the vicinity of the transmitter and receiver for the image method to compute corresponding valid paths to be returned by the solver. These steps of the Sionna RT path solver are depicted in Figure 9. Notably, only path candidates ending with a specular suffix (including specular chains) require further processing by the image method algorithm. Finally, for each valid path identified, the path solver computes the corresponding channel coefficients and propagation delays that characterize the EM wave propagation along that path.

Alternative methods for identifying paths that include specular chains or specular suffixes may involve detection spheres centered at the targets to capture paths ending with a specular suffix, and/or extending the direction of specular reflection (or refraction) to a lobe around the specular, refracted, or diffracted direction. Such methods only approximate the path vertex locations and the directions of the scattered rays, yielding inaccurate path coefficients and resulting in incorrect CIRs.

The remainder of this section follows the architecture of the path solver shown in Figure 9. First, Section 3.1 presents the SBR-based candidate generator. Section 3.2 then describes the image method-based algorithm that processes candidates to find valid paths with specular suffixes. Section 3.3 details how the solver computes the complex-valued channel coefficients and propagation delays, which characterize the radio channel. Finally, Section 3.4 discusses current limitations and potential improvements of the path solver.

3.1 Generating Candidates by Shooting and Bouncing of Rays (SBR)

The candidate generator utilizes the SBR-based method as illustrated in Figure 10. We first provide an overview of the SBR procedure, followed by a detailed explanation. We consider a scene with a single source located at $\mathbf{s}\in\mathbb{R}^{3}$ and $N_{T}$ targets at $\mathbf{t}_{k}\in\mathbb{R}^{3}$ , where $k=1,\dots,N_{T}$ . Although the path solver built into Sionna RT can handle multiple sources, we focus on a single source for the sake of clarity. The process is identical for all sources. The paths traced by the candidate generator during the SBR procedure are referred to as samples in this section to avoid confusion with the valid paths returned to the user and the candidate paths forwarded to the image method-based algorithm. A sample therefore consists of a sequence of rays connected end-to-end, traced by the candidate generator. As we will see, a single sample can result in multiple valid or candidate paths, as illustrated in Figure 11.

The SBR procedure begins by spawning a user-specified number of rays, denoted by $N_{S}$ , from the source $\mathbf{s}$ in directions $\widehat{\mathbf{k}}_{n}^{(0)}$ on a Fibonacci lattice, where $n=1,\dots,N_{S}$ . Each ray starts a sample that is traced by the candidate generator. Details on the implemented Fibonacci lattice can be found in Section 3.1.1. The SBR loop then iterates until either all samples bounce out of the scene or reach a user-specified maximum number of bounces, denoted as the maximum depth $L$ . We denote the depth of the samples by $\ell=0,\dots,L$ . The first step of the SBR loop involves testing the intersection of the samples with the scene geometry. If a sample does not intersect any object, it has bounced out of the scene, and is deactivated. When a sample $n$ with depth $\ell$ intersects the scene, we denote by $\mathbf{v}_{n}^{(\ell)}\in\mathbb{R}^{3}$ the intersection point, by $\widehat{\mathbf{n}}_{n}^{(\ell)}\in\mathbb{R}^{3}~{}\left(\left\lVert\widehat% {\mathbf{n}}_{n}^{(\ell)}\right\rVert_{2}=1\right)$ the surface normal at the intersection point oriented towards the half-space containing the incident field, by $o_{n}^{(\ell)}\in\mathbb{N}_{0}$ the index of the intersected object, and by $m_{n}^{(\ell)}\in\mathbb{N}_{0}$ the index of the specific intersected primitive within that object. The interaction type, denoted by $\chi_{n}^{(\ell)}$ , is then sampled from ${\cal I}$ , determining how the sample will continue to propagate. Details on the sampling of the interaction type are provided in Section 3.1.2. The rest of the loop depends on the sampled interaction type.

Diffuse Reflection

As previously indicated, diffuse reflection is assumed to scatter energy in all directions of the half-space containing the incident field, which enables next event estimation. Therefore, if a diffuse reflection is sampled (i.e. $\chi_{n}^{(\ell)}={\cal S}$ ), the LoS between the intersection point and each target $\mathbf{t}_{1},\cdots,\mathbf{t}_{N_{T}}$ is tested. For each target $\mathbf{t}_{k}$ with an unobstructed LoS to the intersection point, a path connecting the source to the target

p_{n}^{(\ell)}=\left(\mathbf{s},\mathbf{v}_{n}^{(1)},\cdots,\mathbf{v}_{n}^{(% \ell)},\mathbf{t}_{k}\right)

(11)

is recorded as a valid path to be returned to the user. Finally, a new ray is spawned to continue the sample path, with its direction selected uniformly at random from the hemisphere defined by $\widehat{\mathbf{n}}_{n}^{(\ell)}$ , i.e. $\widehat{\mathbf{k}}_{n}^{(\ell)}\sim{\cal U}_{2\pi}(\widehat{\mathbf{n}}_{n}^% {(\ell)})$ , where ${\cal U}_{2\pi}(\widehat{\mathbf{n}}_{n}^{(\ell)})$ represents the uniform distribution on the hemisphere.

Specular Reflection and Refraction

When a specular reflection or refraction is sampled (i.e. $\chi_{n}^{(\ell)}={\cal R}$ or $\chi_{n}^{(\ell)}={\cal T}$ ), the path (11) is a valid path for target $\mathbf{t}_{k}$ , where $k\in\left\{1,\dots,N_{T}\right\}$ , only if

\frac{\mathbf{t}_{k}-\mathbf{v}_{n}^{(\ell)}}{\left\lVert\mathbf{t}_{k}-% \mathbf{v}_{n}^{(\ell)}\right\rVert_{2}}=\widehat{\mathbf{k}}_{n}^{(\ell)},

(12)

where

\widehat{\mathbf{k}}_{n}^{(\ell)}=\begin{cases}\frac{\widehat{\mathbf{k}}_{n}^% {(\ell-1)}-2\left(\left(\widehat{\mathbf{k}}_{n}^{(\ell-1)}\right)^{\mathsf{T}% }\widehat{\mathbf{n}}_{n}^{(\ell)}\right)\widehat{\mathbf{n}}_{n}^{(\ell)}}{% \left\lVert\widehat{\mathbf{k}}_{n}^{(\ell-1)}-2\left(\left(\widehat{\mathbf{k% }}_{n}^{(\ell-1)}\right)^{\mathsf{T}}\widehat{\mathbf{n}}_{n}^{(\ell)}\right)% \widehat{\mathbf{n}}_{n}^{(\ell)}\right\rVert_{2}}&\text{if }\chi_{n}^{(\ell)}% ={\cal R}\\ \widehat{\mathbf{k}}_{n}^{(\ell-1)}&\text{if }\chi_{n}^{(\ell)}={\cal T}\end{cases}

(13)

and where $\mathbf{u}^{\mathsf{T}}$ denotes the transpose of $\mathbf{u}$ . Since targets are modeled as points, condition (12) effectively occurs with zero probability, and, consequently, the path (11) is not a valid path. However, the specular suffix $p_{n}^{(\ell_{d}:)}$ , where $\ell_{d}$ is defined as in (10), can be further processed to determine vertices $\left(\widetilde{\mathbf{v}}_{n}^{(\ell_{d}+1)},\cdots,\widetilde{\mathbf{v}}_% {n}^{(\ell)}\right)$ such that

\widetilde{p}_{n}^{(\ell)}=\left(\mathbf{s},\mathbf{v}_{n}^{(1)},\cdots,% \underbrace{\widetilde{\mathbf{v}}_{n}^{(\ell_{d}+1)},\cdots,\widetilde{% \mathbf{v}}_{n}^{(\ell)},\mathbf{t}_{k}}_{\text{Processed specular suffix}}\right)

(14)

is a valid path, as illustrated in Figure 12. This additional processing is performed using the image method in the algorithm described in Section 3.2.

Importantly, for a sequence of primitives and interaction types

\left(\mathbf{v}_{n}^{(\ell_{d})},\left(o_{n}^{(\ell_{d}+1)},m_{n}^{(\ell_{d}+% 1)},\chi_{n}^{(\ell_{d}+1)}\right),\cdots,\left(o_{n}^{(\ell)},m_{n}^{(\ell)},% \chi_{n}^{(\ell)}\right),\mathbf{t}_{i}\right)

(15)

where $\mathbf{v}_{n}^{(0)}=\mathbf{s}$ , there exists at most a single set of vertices $\widetilde{\mathbf{v}}_{n}^{(\ell_{d}+1)},\cdots,\widetilde{\mathbf{v}}_{n}^{(% \ell)}$ such that $\widetilde{p}_{n}^{(\ell)}$ is a valid path. Consequently, if $\ell_{d}=0$ , i.e. $p_{n}^{(\ell)}$ is a specular chain, then further processing by the image method of multiple identical candidates (15) would result in multiple identical valid paths, which would lead to erroneous CIRs. Therefore, if the sample $p_{n}^{(\ell)}$ is a specular chain, a hashing-based de-duplication mechanism, detailed in Section 3.1.3, ensures that specular chain candidates are only evaluated once. Notably, the path hashing mechanism discards redundant candidates “on-the-fly”, significantly reducing memory and compute requirements. However, when $\ell_{d}>0$ , then $\mathbf{v}_{n}^{(\ell_{d})}$ corresponds to a diffuse reflection point that is the result of a randomly sampled ray direction (either at the source or the last previous diffuse reflection point). Since the probability of two different paths having identical randomly sampled diffuse reflection points is effectively zero, path uniqueness is guaranteed in these cases without requiring the de-duplication mechanism.

As shown in Figure 10, a LoS occlusion test to the target is used as a heuristic to discard path candidates that have little chance of resulting in a valid path after processing by the image method-based algorithm. Finally, a new ray is traced into the direction (13) to continue the sample path.

Diffraction

When tracing rays, the probability that a ray exactly intersects the edge of a wedge is zero. Therefore, when a ray intersects a primitive and diffraction is sampled (i.e., $\chi_{n}^{(\ell)}={\cal D}$ ), the intersection point between the ray $\left(\mathbf{v}_{n}^{(\ell-1)},\widehat{\mathbf{k}}_{n}^{(\ell-1)}\right)$ and the primitive $\left(o_{n}^{(\ell)},m_{n}^{(\ell)}\right)$ is projected onto one of the primitive’s edges, yielding the vertex $\mathbf{v}_{n}^{(\ell)}$ , as illustrated in Figure 13. The index of the edge onto which the intersection point is projected is denoted by $\iota^{(\ell)}_{n}$ .

The path (11) is a valid path to a target $\mathbf{t}_{k}$ (with $k\in\{1,\dots,N_{T}\}$ ) only if there exists an angle $\phi\in[0,n\pi]$ , where $n\pi$ is the exterior angle of the wedge, such that

\frac{\mathbf{t}_{k}-\mathbf{v}_{n}^{(\ell)}}{\left\lVert\mathbf{t}_{k}-% \mathbf{v}_{n}^{(\ell)}\right\rVert_{2}}=\widehat{\mathbf{k}}_{n}^{(\ell)},

(16)

with the diffraction direction given by

\widehat{\mathbf{k}}_{n}^{(\ell)}=\sin{\left(\beta_{0}\right)}\cos{\phi}% \widehat{\mathbf{t}}_{n}^{(\ell)}+\sin{\left(\beta_{0}\right)}\sin{\phi}% \widehat{\mathbf{n}}_{n}^{(\ell)}+\cos{\left(\beta_{0}\right)}\widehat{\mathbf% {e}}_{n}^{(\ell)},

(17)

where $\beta_{0}$ is the angle between the incident ray direction $\widehat{\mathbf{k}}_{n}^{(\ell-1)}$ and the wedge vector $\widehat{\mathbf{e}}_{n}^{(\ell)}$ , and $\left(\widehat{\mathbf{t}}_{n}^{(\ell)},\widehat{\mathbf{n}}_{n}^{(\ell)},% \widehat{\mathbf{e}}_{n}^{(\ell)}\right)$ is an orthonormal basis such that $\widehat{\mathbf{t}}_{n}^{(\ell)}=\widehat{\mathbf{n}}_{n}^{(\ell)}\times% \widehat{\mathbf{e}}_{n}^{(\ell)}$ (see Figure 13). Note that $\widehat{\mathbf{t}}_{n}^{(\ell)}$ is tangent to the wedge face. In (17), ensuring that the angle between the diffracted ray and the edge, denoted by $\beta_{0}$ , matches the angle between the incident ray and the edge, places the ray on the Keller cone. As with specular reflection and refraction, since targets are modeled as points, condition (16) is satisfied with zero probability for all $\phi\in[0,n\pi]$ , rendering the sampled path invalid. Therefore, candidate paths ending with a diffraction are further refined using the image method-based algorithm described in Section 3.2, just as for paths ending with specular reflection or refraction. As shown in Figure 10, the processing of diffraction thus follows the same steps as for specular reflection and refraction: a hashing-based de-duplication mechanism ensures uniqueness of specular chains, and an occlusion test to the target is used as a heuristic to eliminate candidates unlikely to produce valid paths after refinement with the image method. Finally, a new ray is traced in the direction (17) to continue the sample path, with $\phi$ sampled uniformly from $[0,n\pi]$ .

As illustrated in Figure 13, projecting the intersection point onto an edge alters the scattered ray direction $\widehat{\mathbf{k}}_{n}^{(\ell-1)}$ departing from $\mathbf{v}_{n}^{(\ell-1)}$ . This modification invalidates any preceding specular reflections or refractions, meaning that path segments that are not suffixes would require refinement if a diffraction is followed by a diffuse reflection. To avoid significantly increased computational complexity from refining non-suffix path segments, the path solver prohibits sampling paths that contain both diffraction and diffuse reflection. Such paths typically contribute negligibly to the total transported energy and therefore have little impact on the accuracy of the CIRs. Moreover, Sionna RT currently supports only first-order diffraction, meaning that paths may contain at most one diffraction event.

The candidate generator depicted in Figure 10 thus yields a set of valid paths connecting the source to all targets that do not require further processing (except for the computation of the corresponding path coefficients and delays) and end with a diffuse reflection. Additionally, it produces a set of candidate paths ending with specular suffixes, which require further processing to be either discarded or refined into valid paths.

3.1.1 Sampling Initial Ray Directions

The SBR procedure starts by spawning rays from the source with directions obtained from a spherical Fibonacci lattice with $N_{S}$ points, which is defined as follows:

\widehat{\mathbf{k}}_{n}^{(0)}=\begin{bmatrix}\sin{\theta_{n}}\cos{\phi_{n}}\\ \sin{\theta_{n}}\sin{\phi_{n}}\\ \cos{\theta_{n}}\end{bmatrix}

(18)

where

	$\displaystyle\theta_{n}$	$\displaystyle=$	$\displaystyle\arccos{\frac{2n}{N_{S}}}$
	$\displaystyle\phi_{n}$	$\displaystyle=$	$\displaystyle 2\pi\frac{n}{\varphi}$

and $-\left\lfloor\frac{N_{S}}{2}\right\rfloor\leq n\leq\left\lceil\frac{N_{S}}{2}% \right\rceil-1$ , $\varphi=\frac{1+\sqrt{5}}{2}$ is the golden ratio, and $\left\lfloor\cdot\right\rfloor$ and $\left\lceil\cdot\right\rceil$ the floor and ceiling function, respectively. Note that when defining the Fibonacci lattice, the index $n$ takes values in the range $-\left\lfloor\frac{N_{S}}{2}\right\rfloor,\dots,\left\lceil\frac{N_{S}}{2}% \right\rceil-1$ instead of the usual $1,\dots,N_{S}$ . Sampling rays from the Fibonacci lattice results in ray tubes sharing approximately the same opening angle $\frac{4\pi}{N_{S}}$ .

The Fibonacci lattice was selected for its simple construction and its benefits over other lattices for numerical integration [13, 14]. Figure 14 shows the spherical Fibonacci lattice with 100, 1000, and 10000 points. The Sionna RT path solver uses $10^{6}$ points as the default value for $N_{S}$ , and higher values will ensure that more valid paths are found.

3.1.2 Sampling Interaction Types

Within the SBR loop, when a ray intersects the scene geometry, the algorithm generates at most one scattered ray. This constraint is necessary because generating multiple rays at each intersection, such as both reflected and refracted rays, would lead to an exponential increase in computational complexity and memory usage with the depth of the path. Instead, a single interaction type is randomly chosen based on a distribution:

{\cal Q}=\left\{q\left({\cal R}\right),q\left({\cal S}\right),q\left({\cal T}% \right),q\left({\cal D}\right)\right\},\qquad\text{where }q\left({\cal R}% \right)+q\left({\cal S}\right)+q\left({\cal T}\right)+q\left({\cal D}\right)=1% \text{ and }q\left(\cdot\right)\geq 0.

(19)

Here, $q\left(\chi\right)$ represents the probability of selecting the interaction type $\chi\in{\cal I}$ . Importantly, the choice of distribution ${\cal Q}$ does not affect the calculation of the paths’ electric fields, which are solely dependent on the propagation environment. However, the choice of distribution ${\cal Q}$ significantly impacts the sample efficiency of the ray tracer, which affects the number of samples $N_{S}$ needed to compute paths that connect a source to a target and capture the majority of the transported energy. Sampling only one interaction type per intersection may appear limiting, but this is offset by the large number of samples generated from each source. In most scenarios, a high sample count allows identifying all significant paths.

For each sample $n$ at depth $\ell$ , the interaction type is sampled independently for every interaction along every path. For specular reflection, diffuse reflection, and refraction, the probability of selecting each interaction type is proportional to the squared amplitudes of the corresponding reflection and refraction coefficients. Diffraction, on the other hand, is sampled using a fixed probability $q\left({\cal D}\right)\in[0,1]$ . This choice is motivated by the fact that diffraction typically carries much less energy compared to reflection and refraction, which would result in very low probabilities for sampling diffraction if its probability were also based on the diffraction coefficients. The resulting distribution for the interaction types is:

$\displaystyle q\left({\cal R}\right)$	$\displaystyle=\left(1-q\left({\cal D}\right)\right)R^{2}\frac{\left\|r_{\perp}% \right\|^{2}+\left\|r_{\parallel}\right\|^{2}}{\left\|r_{\perp}\right\|^{2}+\left\|r% _{\parallel}\right\|^{2}+\left\|t_{\perp}\right\|^{2}+\left\|t_{\parallel}\right\|^% {2}}$	(20)
$\displaystyle q\left({\cal S}\right)$	$\displaystyle=\left(1-q\left({\cal D}\right)\right)S^{2}\frac{\left\|r_{\perp}% \right\|^{2}+\left\|r_{\parallel}\right\|^{2}}{\left\|r_{\perp}\right\|^{2}+\left\|r% _{\parallel}\right\|^{2}+\left\|t_{\perp}\right\|^{2}+\left\|t_{\parallel}\right\|^% {2}}$
$\displaystyle q\left({\cal T}\right)$	$\displaystyle=\left(1-q\left({\cal D}\right)\right)\frac{\left\|t_{\perp}\right% \|^{2}+\left\|t_{\parallel}\right\|^{2}}{\left\|r_{\perp}\right\|^{2}+\left\|r_{% \parallel}\right\|^{2}+\left\|t_{\perp}\right\|^{2}+\left\|t_{\parallel}\right\|^{2}}$
$\displaystyle q\left({\cal D}\right)$	$\displaystyle\in[0,1]$

In Sionna RT, $q\left({\cal D}\right)$ is set to $0.2$ . Here, $r_{\perp}$ , $r_{\parallel}$ , $t_{\perp}$ , and $t_{\parallel}$ are the Fresnel coefficients (131), $S$ is the scattering coefficient (see Section A.9), and $R$ is the reflection reduction factor (164). This approach ensures that interaction types responsible for scattering more energy are sampled more frequently—-a technique known as importance sampling [15, 16]. This is similar to importance sampling techniques for light transport paths commonly used in computer graphics. A brief overview of importance sampling in the context of ray tracing for the simulation of radio wave propagation is provided in Appendix B.

To illustrate the importance of optimized sampling strategies, we compare uniform sampling, where specular reflection, diffuse reflection, and refraction are selected with equal probability, to importance sampling using the scene shown in Figure 15(a). The scene consists of a metallic box containing a glass screen, with a source and target positioned on opposite sides. Diffuse reflection, specular reflection, refraction, and diffraction are enabled, and the maximum depth $L$ is set to $5$ . As shown in Figure 15(b), using importance sampling leads to better sample efficiency when compared to uniform sampling, as convergence requires approximately $100\times$ fewer samples to capture most of the received energy. This demonstrates the benefit of selecting interactions according to the material properties of the scatterers.

3.1.3 Specular Chain Candidates Deduplication

As mentioned in Section 3.1 and Figure 10, to avoid storing multiple identical specular chain candidates, a path hashing method is used to discard duplicates of previously found candidates “on-the-fly”. A specular chain candidate is defined by a sequence of interactions with the scene, each characterized by the intersected object $o_{n}^{(\ell)}\in\mathbb{N}_{0}$ , the corresponding intersected primitive $m_{n}^{(\ell)}\in\mathbb{N}_{0}$ , the diffracted edge $\iota_{n}^{(\ell)}\in\{0,1,2\}$ (if $\chi_{n}^{(\ell)}={\cal D}$ ), and the interaction type $\chi_{n}^{(\ell)}\in\left\{{\cal R},{\cal T},{\cal D}\right\}$ with $\ell_{d}=0$ , i.e.

c^{(\ell)}_{n}=\left(\mathbf{s},\left(o_{n}^{(1)},m_{n}^{(1)},\iota_{n}^{(1)},% \chi_{n}^{(1)}\right),\cdots,\left(o_{n}^{(\ell)},m_{n}^{(\ell)},\iota_{n}^{(% \ell)},\chi_{n}^{(\ell)}\right)\right).

(21)

For specular chain candidates, we do not need to store the path vertices. If the candidate leads to a valid path, these vertices will be computed. Invalid candidates are discarded.

The hashing process guarantees that each specular chain candidate uniquely represents the geometry of a path after refinement with the image method (see Section 3.2). This is achieved by computing a hash value that uniquely identify the sequence of primitives constituting the path. Importantly, primitives that are coplanar (i.e., lie in the same plane) are treated as equivalent during hashing. This is because the image method assumes all reflection surfaces extend infinitely, making the exact in-plane position of a primitive irrelevant to the path geometry. Furthermore, refraction can be omitted from the hashing procedure since it does not change the ray direction and therefore does not alter the geometry of the path, as refraction is modeled without angular deflection. Figure 16 illustrates this concept: two candidate paths, $p_{1}$ and $p_{2}$ , yield to the same refined path $\widetilde{p}$ , and thus they should share the same hash value. This holds even though only one candidate is refracted by primitive $(a)$ , and both are reflected by different but coplanar primitives $(b)$ and $(c)$ , respectively.

\texttt{hash}_{n}^{(0)}\leftarrow 0

\ell\leftarrow 1

3: while ¡Stop condition¿ do

\triangleright

SBR loop

\cdots

\texttt{inter\_hash}_{n}^{(\ell)}\leftarrow\textsc{hash\_intersection}(o_{n}^{% (\ell)},m_{n}^{(\ell)},\iota_{n}^{(\ell)},\chi_{n}^{(\ell)})

\texttt{hash}_{n}^{(\ell)}\leftarrow\textsc{hash\_update}(\texttt{hash}_{n}^{(% \ell-1)},\texttt{inter\_hash}_{n}^{(\ell)})

\cdots

\ell\leftarrow\ell+1

9: end while

10:

11: procedure hash_intersection(

o,m,\iota,\chi

)

12: if

\chi={\cal R}

then return

\textsc{hash\_plane}(o,m)

13: else if

\chi={\cal D}

then return

\textsc{hash\_edge}(o,m,\iota)

14: end if

15: end procedure

16:

17: procedure hash_plane(

o, m

)

18:

\mathbf{n}\leftarrow\textsc{get\_normal}(o,m)

\triangleright

Read the primitive normal

19:

\mathbf{p}\leftarrow\textsc{get\_vertex}(o,m)

\triangleright

Read a primitive vertex

20:

d\leftarrow\mathbf{n}^{\mathsf{T}}\mathbf{p}

\triangleright

Offset of the plane to the origin

21:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(\mathbf{n}_{x}),\texttt{H})

22:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(\mathbf{n}_{y}),h)

23:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(\mathbf{n}_{z}),h)

24:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(d),h)

25: return

h

26: end procedure

27:

28: procedure hash_edge(

o,m,\iota

)

29:

\mathbf{p}_{1},\mathbf{p}_{2}\leftarrow\textsc{get\_endpoints}(o,m,\iota)

\triangleright

Read the edge endpoints

30:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(\mathbf{p}_{1,x}),\texttt{H})

31:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(\mathbf{p}_{1,y}),h)

32:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(\mathbf{p}_{1,z}),h)

33:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(\mathbf{p}_{2,x}),h)

34:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(\mathbf{p}_{2,y}),h)

35:

h\leftarrow\textsc{FNV1a}(\texttt{Q}(\mathbf{p}_{2,z}),h)

36: return

h

37: end procedure

38:

39: procedure hash_update(

\texttt{curr\_hash},\texttt{inter\_hash}

)

40: return

(1373\cdot\texttt{curr\_hash}+\texttt{inter\_hash})

\triangleright

Update rolling hash

41: end procedure

Algorithm 1 Hashing of specular chains

Algorithm 1 depicts how path hashing works during ray tracing to avoid computing duplicate specular chain candidates. It starts by initializing a hash value $\texttt{hash}_{n}^{(0)}$ to $0$ for each sample $n=1,\dots,N_{S}$ being traced. Then, as the sample interacts with objects in the scene, the algorithm computes a hash value $\texttt{inter\_hash}_{n}^{(\ell)}$ for each intersection, and updates the path hash $\texttt{hash}_{n}^{(\ell)}$ by combining it with the hash of the intersection, creating a unique identifier for the entire sample up to that point. The computation of the interaction hash depend on the interaction type. For specular reflection, hashing is performed based on the primitive normal vector and an offset positioning the plane with respect to the origin. An integer hash value is then computed by sequentially combining the normal components and the offset using the FNV-1a hash function [17]. As numerical imprecision can lead to slight variations in the normal vectors of coplanar primitives, the real-valued normal components are quantized to a fixed precision and then converted to integers using a function $\texttt{Q}(\cdot)$ prior to hashing. Discussion on the choice of the quantization function $\texttt{Q}(\cdot)$ is provided below. For diffraction, hashing is performed on the two edge endpoints by applying the FNV-1a hash function sequentially to the components of the endpoints. In both $\textsc{hash\_plane}($ ) and $\textsc{hash\_edge}($ ), the first call to the FNV-1a hash function uses a large constant H to seed the hash function.

Once the interaction hash is computed, the hash_update function updates the sample hash value $\texttt{hash}_{n}^{(\ell)}$ through a polynomial rolling hash function [18] with a prime base of $1373$ . For each new intersection, it updates the running hash by multiplying the current hash value by the prime number and adding the new intersection hash. When a specular chain candidate is considered for a target $k\in\left\{1,\dots,N_{T}\right\}$ , the hash value of the candidate is paired with the target index $k$ to form a new hash value

\texttt{hash}_{n,k}^{(\ell)}=\textsc{FNV1a}(\texttt{hash}_{n}^{(\ell)},k).

(22)

This ensures that the specular chain candidate is unique for each target.

Practical Consideration: Hash Array

To track information about previously encountered specular chains, an array of integers of size $N_{H}$ is allocated in memory and initialized to zero, referred to as the hash array. When considering a specular chain candidate for target $k$ with corresponding hash value $\texttt{hash}_{n,k}^{(\ell)}$ , the corresponding index $i_{n,k}^{(\ell)}$ in the hash array is computed as

i_{n,k}^{(\ell)}=\texttt{hash}_{n,k}^{(\ell)}\mod N_{H}.

(23)

If the hash array entry at index $i_{n,k}^{(\ell)}$ is greater than zero, it indicates that this path has already been encountered and should be discarded. Otherwise, the entry is incremented by one and the candidate is retained for further processing by the image method. It is important to note that the index $i_{n,k}^{(\ell)}$ is not guaranteed to be unique for different specular chains. When two different specular chains share the same index value despite being different, the path solver will discard one of them. Such an event is referred to as a collision, and can result in the discarding of valid paths. However, when two candidates are identical, their indices values will also be identical.

Valid paths and candidates are stored in a buffer of size $N_{B}$ . This buffer holds, for each path, the sequence of interaction types, vertices, intersected objects, primitives, and more. When a valid path is found or a candidate is identified as unique, indicated by its corresponding entry $i_{n,k}^{(\ell)}$ (23) in the hash array being set to zero, the path data are stored in the buffer. It is important to note that the index $i_{n,k}^{(\ell)}$ is solely used for indexing the hash array and not for storing path data in the buffer. Instead, path data are stored contiguously in the next available slot of the path buffer, as shown in Figure 17. If the buffer becomes full, any newly found paths or candidates are discarded. Discarding paths rather than overwriting existing ones is based on the reason that paths with lower depth, which are added to the buffer first as the SBR loop advances from $\ell=0$ to $\ell=L$ , typically carry more energy and should not be replaced by higher depth paths when the buffer becomes full. Consequently, the buffer size $N_{B}$ determines the maximum number of paths and candidates that can be stored.

Due to the substantial amount of data stored for each path, setting $N_{B}$ to excessively high values can lead to memory issues. This necessitates separating the hash array size $N_{H}$ from the buffer size $N_{B}$ . This is because a small hash array size can increase the probability of collisions. To reduce the risk of collisions, $N_{H}$ is set to

N_{H}=\max\left(N_{B},10^{6}\right)

(24)

ensuring a minimum array size of $10^{6}$ elements. In the interface of the Sionna RT built-in path solver, $N_{B}$ corresponds to the samples_per_src parameter, and $N_{H}$ corresponds to max_num_paths_per_src.

Practical Consideration: Quantization Function

An important practical consideration is the design of the quantization function $\texttt{Q}(\cdot)$ . When quantizing floating-point values to integers, values that are very close but lie on opposite sides of a quantization boundary can be mapped to different integers and therefore yield different hash values. For example, with rounding as the quantization function, $\texttt{Q}(1.4999)=1$ and $\texttt{Q}(1.5001)=2$ . This can cause identical candidates, which differ only due to numerical imprecision, to be treated as distinct, resulting in duplicate paths. To mitigate this issue, multiple hash values are computed for each sample, each using a different quantization function. A path is only considered unique if all of its corresponding hash values are unique. In the current implementation of Sionna RT, two quantization functions are used: rounding and floor, leading to two hash values being computed per candidate. This approach can be interpreted as a form of locality-sensitive hashing [19], as it aims to group together candidates that differ only due to numerical imprecision by using multiple hash functions to reduce the risk of misclassification.

3.2 Image Method-based Candidate Processing

The image method processes path candidates to determine valid paths that are specular chains, or end with specular suffixes. Consider a candidate path

c_{n}^{(\ell)}=\left(\mathbf{s},\left(o_{n}^{(1)},m_{n}^{(1)},\iota_{n}^{(1)},% \chi_{n}^{(1)}\right),\cdots,\underbrace{\left(o_{n}^{(\ell_{d}+1)},m_{n}^{(% \ell_{d}+1)},\iota_{n}^{(\ell_{d}+1)},\chi_{n}^{(\ell_{d}+1)}\right),\cdots,% \left(o_{n}^{(\ell)},m_{n}^{(\ell)},\iota_{n}^{(\ell)},\chi_{n}^{(\ell)}\right% ),\mathbf{t}}_{\text{Specular suffix}}\right)

(25)

where $\ell$ denotes the path depth, $\ell_{d}$ is the specular suffix index defined as in (10), and $\mathbf{t}$ is the target of the path. Recall that the sequence of interactions $\boldsymbol{\chi}^{(\ell)}=\left(\chi_{n}^{(1)},\cdots,\chi_{n}^{(\ell)}\right)$ was determined by the SBR-based candidate generator. For path candidates that end with a specular suffix but are not specular chains (i.e. $\ell_{d}>0$ ), the image method processes only the specular suffix. In these cases, the vertex of the last diffuse reflection $\mathbf{v}_{n}^{(\ell_{d})}$ serves as the source point for the image method calculations. Since the processing algorithm is identical for both specular chains and specular suffixes, we treat them equivalently in the following discussion. For simplicity of notation, we will refer to both types as specular chains and assume $\ell_{d}=0$ throughout the remainder of this section.

Paths without Diffraction

We begin by describing the image method process for candidates that contain no diffraction; as illustrated in Figure 18, this process consists of two main steps. First, the method iteratively computes images of the source by reflecting it across each plane at which a specular reflection occurs. Since refracted paths are modeled without angular deflection, incorporating refraction into the image method is straightforward: refractions may simply be skipped when computing the source images, as they do not alter the ray’s direction. Then, the algorithm backtracks from the target to determine the actual path vertices by intersecting lines between successive image sources.

Formally, the image method iteratively computes images of the source $\underline{\mathbf{s}}_{n}^{(i)}$ for $i=0,\dots,\ell$ , starting from the original source $\mathbf{s}$ :

\underline{\mathbf{s}}_{n}^{(i)}=\begin{cases}\mathbf{s}&\text{if }i=0\\ \underline{\mathbf{s}}_{n}^{(i-1)}-2\left(\left(\underline{\mathbf{s}}_{n}^{(i% -1)}-\mathbf{p}_{n}^{(i)}\right)^{\mathsf{T}}\widehat{\mathbf{n}}_{n}^{(i)}% \right)\widehat{\mathbf{n}}_{n}^{(i)}&\text{if }\chi_{n}^{(i)}={\cal R}\\ \underline{\mathbf{s}}_{n}^{(i-1)}&\text{if }\chi_{n}^{(i)}={\cal T}\end{cases}

(26)

Here, $\mathbf{p}_{n}^{(i)}$ represents any point on the plane containing the primitive $\left(o_{n}^{(i)},m_{n}^{(i)}\right)$ , and $\widehat{\mathbf{n}}_{n}^{(i)}$ denotes the surface’s normal vector. For specular reflections ( $\chi_{n}^{(i)}={\cal R}$ ), the image is computed by reflecting the previous image across the surface plane, while for refractions ( $\chi_{n}^{(i)}={\cal T}$ ), the image position remains unchanged.

During backtracking, the algorithm traces rays iteratively from the target $\mathbf{t}$ back to the source $\mathbf{s}$ by intersecting lines between successive image points. The intersection points $\widetilde{\mathbf{v}}_{n}^{(i)}$ are the refined vertices of the path. Specifically, a sequence of $\ell$ rays $\left(\widetilde{\mathbf{v}}_{n}^{(i+1)},\widehat{\mathbf{d}}_{n}^{(i)}\right)$ , $i=\ell,\dots,1$ , is traced, where $\widetilde{\mathbf{v}}_{n}^{(\ell+1)}=\mathbf{t}$ and

	$\displaystyle\widehat{\mathbf{d}}_{n}^{(i)}$	$\displaystyle=\frac{\underline{\mathbf{s}}_{n}^{(i)}-\widetilde{\mathbf{v}}_{n% }^{(i+1)}}{\left\lVert\underline{\mathbf{s}}_{n}^{(i)}-\widetilde{\mathbf{v}}_% {n}^{(i+1)}\right\rVert_{2}}$		(27)
	$\displaystyle\widetilde{\mathbf{v}}_{n}^{(i)}$	$\displaystyle=\texttt{intersection}\left(\left(\widetilde{\mathbf{v}}_{n}^{(i+% 1)},\widehat{\mathbf{d}}_{n}^{(i)}\right),\left(o_{n}^{(i)},m_{n}^{(i)}\right)% \right),\text{only if }i>0$		(27)

where $\widetilde{\mathbf{v}}_{n}^{(i)}$ is the intersection point between the ray $\left(\widetilde{\mathbf{v}}_{n}^{(i+1)},\widehat{\mathbf{d}}_{n}^{(i)}\right)$ and the primitive $\left(o_{n}^{(i)},m_{n}^{(i)}\right)$ , as shown in Figure 18(b). The path is considered valid only if two conditions are met: (i) each ray $\left(\widetilde{\mathbf{v}}_{n}^{(i+1)},\widehat{\mathbf{d}}_{n}^{(i)}\right)$ successfully intersects a primitive that is coplanar with the primitive $\left(o_{n}^{(i)},m_{n}^{(i)}\right)$ , and (ii) the line segment between $\widetilde{\mathbf{v}}_{n}^{(i+1)}$ and $\widetilde{\mathbf{v}}_{n}^{(i)}$ is unobstructed by any other scene geometry. If either condition fails, then the path is discarded. The case $i=0$ corresponds to the final segment of the path, which connects to the source $\mathbf{s}$ . If the path is valid, then the returned path will be (14).

Paths with Diffraction

To handle paths involving diffraction, we adopt an approach inspired by [20]. Since Sionna RT currently only supports first-order diffraction, each sample $n$ contains at most a single diffracting edge, whose position within the path is denoted by $i_{{\cal D}}\in\left\{1,\dots,\ell\right\}$ . For notational simplicity, we omit the explicit dependence of $i_{{\cal D}}$ on $n$ .

The core idea is to use image theory to compute the diffraction point in image space and subsequently map it back to the original space. The motivation for this approach is that, in the image space, all specular reflections and refractions are effectively “removed”, leaving only the diffraction event. This reduces the problem to a simpler form: finding the appropriate diffraction point for a path of the shape source $\rightarrow$ diffraction $\rightarrow$ target. More precisely, for every sample $n$ , the diffracting edge at depth $i_{{\cal D}}$ , which is characterized by an endpoint $\mathbf{o}_{n}$ and an edge vector $\widehat{\mathbf{e}}_{n}$ , has its images computed by reflecting its endpoint and edge vector across each plane at which a specular reflection occurs over the following interactions $i=i_{{\cal D}}+1,\dots,\ell$ :

	$\displaystyle\underline{\mathbf{o}}_{n}^{(i)}$	$\displaystyle=\begin{cases}\mathbf{o}_{n}&\text{if }i=i_{{\cal D}}\\ \underline{\mathbf{o}}_{n}^{(i-1)}&\text{if }\chi_{n}^{(i)}\neq{\cal R}\\ \underline{\mathbf{o}}_{n}^{(i-1)}-2\left(\left(\underline{\mathbf{o}}_{n}^{(i% -1)}-\mathbf{p}_{n}^{(i)}\right)^{\mathsf{T}}\widehat{\mathbf{n}}_{n}^{(i)}% \right)\widehat{\mathbf{n}}_{n}^{(i)}&\text{if }\chi_{n}^{(i)}={\cal R}\end{cases}$		(28)
	$\displaystyle\widehat{\underline{\mathbf{e}}}_{n}^{(i)}$	$\displaystyle=\begin{cases}\widehat{\mathbf{e}}_{n}&\text{if }i=i_{{\cal D}}\\ \widehat{\underline{\mathbf{e}}}_{n}^{(i-1)}&\text{if }\chi_{n}^{(i)}\neq{\cal R% }\\ \widehat{\underline{\mathbf{e}}}_{n}^{(i-1)}-2\left(\left(\widehat{\underline{% \mathbf{e}}}_{n}^{(i-1)}\right)^{\mathsf{T}}\cdot\widehat{\mathbf{n}}_{n}^{(i)% }\right)\widehat{\mathbf{n}}_{n}^{(i)}&\text{if }\chi_{n}^{(i)}={\cal R}\end{cases}$		(28)

where $\mathbf{p}_{n}^{(i)}$ represents any point on the plane containing the primitive $\left(o_{n}^{(i)},m_{n}^{(i)}\right)$ , and $\widehat{\mathbf{n}}_{n}^{(i)}$ denotes the surface’s normal vector. Applying image theory, specular reflections and refractions can be ignored when searching for the diffraction point in image space. The problem thus reduces to finding the diffraction point in image space, $\underline{\widetilde{\mathbf{v}}}_{n}^{(i_{{\cal D}},\ell)}$ , of the path $\underline{\mathbf{s}}_{n}^{(\ell)}\rightarrow\underline{\widetilde{\mathbf{v}% }}_{n}^{(i_{{\cal D}},\ell)}\rightarrow\mathbf{t}$ . According to Fermat’s principle, the diffraction point minimizes the total path length, denoted by ${\cal L}(x_{n})$ . Specifically,

	$\displaystyle{\cal L}(x_{n})$	$\displaystyle=\left\lVert\mathbf{t}-\underline{\widetilde{\mathbf{v}}}_{n}^{(i% _{{\cal D}},\ell)}\right\rVert_{2}+\left\lVert\underline{\mathbf{s}}_{n}^{(% \ell)}-\underline{\widetilde{\mathbf{v}}}_{n}^{(i_{{\cal D}},\ell)}\right% \rVert_{2}$		(29)
	$\displaystyle\underline{\widetilde{\mathbf{v}}}_{n}^{(i_{{\cal D}},\ell)}$	$\displaystyle=\underline{\mathbf{o}}_{n}^{(\ell)}+x_{n}\widehat{\underline{% \mathbf{e}}}_{n}^{(\ell)}$		(30)

where $x_{n}\in\mathbb{R}$ . For first-order diffraction, this minimization has a closed-form solution (see Annex C).

Once the diffraction point is found in the image space, the refined path vertices $\widetilde{\mathbf{v}}_{n}^{(i)}$ , $i=1,\dots,\ell$ , are determined by backtracking as in the case without diffraction. However, for depths $i=i_{{\cal D}}+1,\dots,\ell$ , the algorithm uses the intermediate images of the diffraction point:

\widetilde{\underline{\mathbf{v}}}_{n}^{(i_{{\cal D}},i)}=\underline{\mathbf{o% }}_{n}^{(i)}+x_{n}\widehat{\underline{\mathbf{e}}}_{n}^{(i)}

(31)

in place of $\underline{\mathbf{s}}_{n}^{(i)}$ in (27). For $i=1,\dots,i_{{\cal D}}-1$ , $\underline{\mathbf{s}}_{n}^{(i)}$ is used as in the non-diffraction case. The diffraction point is simply $\widetilde{\mathbf{v}}_{n}^{(i_{{\cal D}})}=\mathbf{o}_{n}+x_{n}\widehat{% \mathbf{e}}_{n}$ . This process is depicted in Figure 19.

Typically, most of the candidates are discarded by the image method. This is illustrated in Figure 20, which shows the number of candidates discarded as a function of the path depth for a single source and target. The considered scene is shown in Figure 20(a). The number of samples $N_{S}$ is $10^{6}$ , diffuse reflection is disabled, and diffraction is enabled. As seen in Figure 20(b), most of the candidates are discarded by the image method because only specular chains are generated by the candidate generator, and most of them do not result in valid paths.

3.3 Channel Coefficients, Delays, and Doppler Shifts Computation

The final step of the path solver computes complex-valued channel coefficients and real-valued delays for each valid path (see Figure 9). For clarity, we will focus our discussion on a single valid path with index $n$ and depth $\ell$ , though the same procedure applies to all valid paths computed by both the SBR-based candidate generator and image method. We denote the path vertices as $\mathbf{v}_{n}^{(i)},~{}i=1,\dots,\ell$ , which may correspond to vertices previously computed by the image method (denoted earlier as $\widetilde{\mathbf{v}}_{n}^{(i)}$ ). The directions of propagation between vertices are represented by unit vectors $\widehat{\mathbf{k}}_{n}^{(i)},~{}i=0,\dots,\ell$ , defined as

\widehat{\mathbf{k}}_{n}^{(i)}=\frac{\mathbf{v}_{n}^{(i+1)}-\mathbf{v}_{n}^{(i% )}}{\left\lVert\mathbf{v}_{n}^{(i+1)}-\mathbf{v}_{n}^{(i)}\right\rVert_{2}},% \qquad i=0,\dots,\ell

(32)

where $\mathbf{v}_{n}^{(0)}=\mathbf{s}$ is the source and $\mathbf{v}_{n}^{(\ell+1)}=\mathbf{t}$ is the target. Note that these propagation directions may differ from those used during the initial SBR loop described in Section 3.1 due to the refinement by the image method.

Ensuring Energy Conservation

As explained in Section 3.1.2, at each interaction point, an interaction type is randomly chosen based on the probability distribution ${\cal Q}$ (19). Initially, the energy emitted by the source is distributed across $N_{S}$ samples according to the transmit antenna pattern (see Section A.3). Since only one interaction type is selected at each interaction point, ray tubes that do not match the chosen interaction type are discarded, as illustrated in Figure 21. The remaining ray tubes must therefore compensate for the discarded ones to maintain energy conservation. Given that there is a finite number of specular chains with a specific depth, the channel gain (5) observed at the target due to specular chains is accurate if all specular chains are identified by the path solver, which is assumed to be the case. This assumption will be reasonable if the number of samples $N_{S}$ is sufficiently large.

In the context of diffuse reflections, an incident ray tube scatters an amount of energy that is determined by its footprint on the diffusely reflective surface (see Section A.9). Due to the random sampling of interaction types, only a portion of the ray tubes that should reach a diffusely reflective surface are actually traced (see Figure 21), which requires compensating for the untraced ones. For a given path $p^{(\ell)}$ that includes a diffuse reflection and has a specular suffix index $\ell_{d}$ (10), we denote by $\Pr\left(\boldsymbol{\chi}_{n}^{(\ell_{d})}\right)=\prod_{\ell^{\prime}=1}^{% \ell_{d}}q\left(\chi_{n}^{(\ell^{\prime})}\right)$ the probability of sampling the prefix leading to the last diffuse reflection point, where $q\left(\cdot\right)$ is defined in (19). Consequently, on average, only a fraction $\Pr\left(\boldsymbol{\chi}_{n}^{(\ell_{d})}\right)$ of the ray tubes that should reach the diffusely reflective surface are actually traced. To account for this, the electric field is scaled by $\nicefrac{{1}}{{\sqrt{\Pr\left(\boldsymbol{\chi}_{n}^{(\ell_{d})}\right)}}}$ at the point of diffuse reflection.

Algorithm 2 Electric field computation

\mathbf{E}_{n}\leftarrow\textsc{tx\_pattern}(\widehat{\mathbf{k}}^{(0)}_{n})

\triangleright

Electric field vector

\tau_{n}\leftarrow 0

\triangleright

Cumulative delay

r_{n}\leftarrow 0

\triangleright

Length of the ray tube

\gamma_{n}\leftarrow 1

\triangleright

Cumulative path probability

s_{n}\leftarrow 0

\triangleright

Length of the incident ray tube on the diffracting edge

6: for

\ell^{\prime}=1,\dots,\ell

\mathbf{T}_{n}^{(\ell^{\prime})}\leftarrow\textsc{Evaluate\_Material}(o_{n}^{(% \ell^{\prime})},\chi_{n}^{(\ell^{\prime})},\dots)

\mathbf{E}_{n}\leftarrow\mathbf{T}_{n}^{(\ell^{\prime})}\mathbf{E}_{n}

\gamma_{n}\leftarrow\gamma_{n}\cdot q\left(\chi_{n}^{(\ell^{\prime})}\right)

10:

r_{n}\leftarrow r_{n}+\left\lVert\mathbf{v}_{n}^{(\ell^{\prime})}-\mathbf{v}_{% n}^{(\ell^{\prime}-1)}\right\rVert_{2}

11: if

\chi_{n}^{(\ell^{\prime})}={\cal S}

then

12:

\mathbf{E}_{n}\leftarrow\frac{\mathbf{E}_{n}}{\sqrt{\gamma_{n}}}

13:

\gamma_{n}\leftarrow 1

14:

r_{n}\leftarrow 0

15: else if

\chi_{n}^{(\ell^{\prime})}={\cal D}

then

16:

s_{n}\leftarrow r_{n}

17:

r_{n}\leftarrow 0

18: end if

19:

\tau_{n}\leftarrow\tau_{n}+\frac{\left\lVert\mathbf{v}_{n}^{(\ell^{\prime})}-% \mathbf{v}_{n}^{(\ell^{\prime}-1)}\right\rVert_{2}}{c}

20: end for

21:

r_{n}\leftarrow r_{n}+\left\lVert\mathbf{t}-\mathbf{v}_{n}^{(\ell)}\right% \rVert_{2}

22:

\tau_{n}\leftarrow\tau_{n}+\frac{\left\lVert\mathbf{t}-\mathbf{v}_{n}^{(\ell)}% \right\rVert_{2}}{c}

23:

\mathbf{E}_{n}^{\text{rx}}\leftarrow\textsc{rx\_pattern}(-\widehat{\mathbf{k}}% _{n}^{(\ell)})

24:

a_{n}\leftarrow\frac{\lambda}{4\pi}\left(\mathbf{E}_{n}^{\mathsf{T}}\mathbf{E}% _{n}^{\text{rx}}\right)

25: if Diffracted path then

26:

a_{n}\leftarrow\frac{a_{n}}{\sqrt{s_{n}r_{n}(s_{n}+r_{n})}}

27: else

28:

a_{n}\leftarrow\frac{a_{n}}{r_{n}}

29: end if

30: return

a_{n},~{}\tau_{n}

Path Coefficient and Delay Computation

Algorithm 2 outlines the computation of the electric field along a path, based on the EM background in Appendix A. The algorithm calculates the electric field vector for each path by iterating through its sequence of interactions. Computation begins by initializing the electric field $\mathbf{E}_{n}$ with the transmit antenna pattern, evaluated in the initial departure direction $\widehat{\mathbf{k}}_{n}^{(0)}$ . For dual polarization, two independent electric field vectors are tracked per path. The variable $\tau_{n}$ accumulates the total propagation delay, and $r_{n}$ tracks the ray tube length; both start at zero. The cumulative path probability $\gamma_{n}$ , used for energy normalization, is initialized to one (see previous paragraph) and updated at each interaction. The variable $s_{n}$ is only used for diffracted paths. For these paths, it accumulates the length of the incident ray tube on the diffracting edge, which is required to compute the diffraction spreading factor, i.e., the field intensity decrease factor due to free space propagation (see Section A.8).

At each interaction, the electric field is updated by the transfer matrix $\mathbf{T}_{n}^{(\ell^{\prime})}$ , which depends on the intersected object $o_{n}^{(\ell^{\prime})}$ , the interaction type $\chi_{n}^{(\ell^{\prime})}$ , and other parameters (omitted for brevity). The path probability $\gamma_{n}$ is multiplied by the current interaction probability $q\left(\chi_{n}^{(\ell^{\prime})}\right)$ . For diffuse reflections, the electric field is divided by $\sqrt{\gamma_{n}}$ to ensure energy conservation and $\gamma_{n}$ and $r_{n}$ are reset to one and zero, respectively, as diffuse reflection spawns a new ray tube. For diffractions, the length of the incident ray tube on the diffracting edge $r_{n}$ is stored in $s_{n}$ and $r_{n}$ is reset to zero. Each segment’s travel time is added to $\tau_{n}$ , with $c$ the speed of light.

After processing all interactions, the algorithm adds the final segment to the target to both $\tau_{n}$ and $r_{n}$ . The receive antenna pattern is then evaluated for the incoming direction, and the channel coefficient is computed following Section A.6. Finally, the spreading factor (lines 25–29) is computed and applied to the channel coefficient depending on whether the path involves diffraction or not.

Remark

The transfer matrices $\mathbf{T}_{n}^{(\ell^{\prime})}$ and the electric field vectors $\mathbf{E}_{n}$ are complex-valued, but are implemented using their real-valued representations defined as

\begin{bmatrix}\Re\left(\mathbf{T}_{n}^{(\ell^{\prime})}\right)&-\Im\left(% \mathbf{T}_{n}^{(\ell^{\prime})}\right)\\ \Im\left(\mathbf{T}_{n}^{(\ell^{\prime})}\right)&\Re\left(\mathbf{T}_{n}^{(% \ell^{\prime})}\right)\end{bmatrix}

(33)

for matrices and

\begin{bmatrix}\Re\left(\mathbf{E}_{n}\right)\\ \Im\left(\mathbf{E}_{n}\right)\end{bmatrix}

(34)

for vectors, where $\Re{\left(\cdot\right)}$ and $\Im{\left(\cdot\right)}$ denote the real and imaginary parts, respectively.

3.3.1 Handling Multiple Antennas at the Transmitter and Receiver

As explained in the introduction (Section 1), Sionna RT supports antenna arrays of any size at both the transmitter and receiver, utilizing either synthetic or non-synthetic arrays. In the case of non-synthetic arrays, each transmit antenna is treated as a source and each receive antenna as a target, meaning paths are traced from every transmit antenna, and the solver calculates paths for each transmit-receive antenna pair. For synthetic arrays, paths originate from a single “virtual” source at the center of the transmitter array, and a single “virtual” target is considered at the center of every receiver array. Channel coefficients are then computed by applying phase shifts synthetically. This method allows for efficient scaling with the number of antennas, but it is an approximation that is valid only when the transmitter and receiver array sizes are small relative to the distance from the radio devices to each other and to the scatterers. Formally, for synthetic arrays, the channel matrix for the $n$ -th path is expressed as:

\mathbf{H}_{\text{array},n}=a_{n}\mathbf{u}_{A}^{\text{rx}}\left(\mathbf{u}_{A% }^{\text{tx}}\right)^{\mathsf{T}}

(35)

where $a_{n}$ is the channel coefficient for the $n$ -th path, calculated using Algorithm 2. The array response vectors for the transmit ( $\mathbf{u}_{A}^{\text{tx}}$ ) and receive ( $\mathbf{u}_{A}^{\text{rx}}$ ) antenna arrays are defined as:

	$\displaystyle\mathbf{u}_{A}^{\text{tx}}$	$\displaystyle=\left[e^{j\frac{2\pi}{\lambda}\left(\widehat{\mathbf{k}}_{n}^{(0% )}\right)^{\mathsf{T}}\mathbf{d}_{1}^{\text{tx}}},\cdots,e^{j\frac{2\pi}{% \lambda}\left(\widehat{\mathbf{k}}_{n}^{(0)}\right)^{\mathsf{T}}\mathbf{d}_{N_% {A}^{\text{tx}}}^{\text{tx}}}\right]^{\mathsf{T}},$		(36)
	$\displaystyle\mathbf{u}_{A}^{\text{rx}}$	$\displaystyle=\left[e^{j\frac{2\pi}{\lambda}\left(-\widehat{\mathbf{k}}_{n}^{(% \ell)}\right)^{\mathsf{T}}\mathbf{d}_{1}^{\text{rx}}},\cdots,e^{j\frac{2\pi}{% \lambda}\left(-\widehat{\mathbf{k}}_{n}^{(\ell)}\right)^{\mathsf{T}}\mathbf{d}% _{N_{A}^{\text{rx}}}^{\text{rx}}}\right]^{\mathsf{T}}.$		(36)

Here, $\widehat{\mathbf{k}}_{n}^{(0)}$ ( $\widehat{\mathbf{k}}_{n}^{(\ell)}$ ) represents the departure (arrival) direction of the $n$ -th path, $N_{A}^{\text{tx}}$ ( $N_{A}^{\text{rx}}$ ) denotes the number of antennas in the transmit (receive) array, and $\mathbf{d}_{i}^{\text{tx}}$ ( $\mathbf{d}_{i}^{\text{rx}}$ ) indicates the relative position of the $i$ -th antenna element to the virtual source (target) in the transmit (receive) array.

3.3.2 Doppler Shifts

The impact of moving scene objects on the paths can be simulated in two ways. The first approach involves moving objects in small incremental steps along their trajectories and recomputing the propagation paths at each step. While this method provides the highest accuracy, it is computationally expensive. However, when trajectory lengths are small (not exceeding a few wavelengths), a more efficient approach may be used. In this alternative method, we assign a velocity vector $\boldsymbol{\varv}\left(o\right)\in\mathbb{R}^{3}$ to each object $o$ and compute the time evolution by accumulating Doppler shifts along the propagation paths.

Figure 22 illustrates this approach. When replaying the valid paths using Algorithm 2, the Doppler shift for each path $n$ is computed by accumulating the contributions from individual objects:

\nu_{n}^{(\ell^{\prime})}=\begin{cases}0&\text{if }\ell^{\prime}=0\\ \nu_{n}^{(\ell^{\prime}-1)}+\frac{\boldsymbol{\varv}\left(o_{n}^{(\ell^{\prime% })}\right)^{\mathsf{T}}\left(\widehat{\mathbf{k}}_{n}^{(\ell^{\prime})}-% \widehat{\mathbf{k}}_{n}^{(\ell^{\prime}-1)}\right)}{\lambda}&\text{if }1\leq% \ell^{\prime}<\ell\\ \end{cases}

(37)

where $\lambda$ is the wavelength. This results in the total Doppler shift $\nu_{n}^{(\ell)}$ caused by object movements along the path. The final Doppler shift observed at the target combines the accumulated shift from moving objects with the shifts caused by source and target motion:

\nu_{n}=\nu_{n}^{(\ell)}+\frac{\boldsymbol{\varv}\left(\mathbf{s}\right)^{% \mathsf{T}}\widehat{\mathbf{k}}_{n}^{(0)}}{\lambda}-\frac{\boldsymbol{\varv}% \left(\mathbf{t}\right)^{\mathsf{T}}\widehat{\mathbf{k}}_{n}^{(\ell)}}{\lambda}.

(38)

where $\boldsymbol{\varv}\left(\mathbf{s}\right)\in\mathbb{R}^{3}$ and $\boldsymbol{\varv}\left(\mathbf{t}\right)\in\mathbb{R}^{3}$ denote the velocity vectors of the source and target, respectively.

3.4 Current Limitations and Potential Improvements

We discuss in this section some current limitations and potential improvements of the path solver.

Sampling of Ray Directions

Improved sample efficiency could be achieved by sampling initial ray directions and scattered ray directions based on the transmit antenna pattern and the scattering pattern, respectively. This approach is motivated by the substantial sample efficiency gains shown in Section 3.1.2 through the use of importance sampling. By applying optimized sampling strategies for ray directions, further improvements may be realized, as illustrated in Figure 23.

Handling Large Numbers of Targets and Sources

The current method for handling multiple targets involves iterating over all targets during each SBR loop iteration. Figure 24(b) illustrates this, showing the compute time of the candidate generator and the average number of candidates per target found as the number of targets varies. In this experiment, the built-in Sionna RT scene “simple street canyon” depicted in Figure 24(a), was used. Targets were sampled uniformly at random on a plane parallel to the ground at an elevation of 1.5 meters. The maximum depth $L$ was set to 5, diffuse reflections were disabled, and the number of samples $N_{S}$ was set to $10^{6}$ . As observed, the compute time increases with the number of targets due to the iteration over all targets in each SBR loop iteration. However, the number of candidates per target remains stable, although with variations due to hash collisions. Note that, as discussed in Section 3.2, most candidates do not result in valid paths. This approach may be optimized by selecting only targets near intersection points or the source, which will enhance scalability when dealing with a large number of targets.

Scaling with the number of sources presents its own set of challenges. When $N_{S}$ samples are generated per source, the total number of samples produced by the candidate generator becomes $N_{S}\cdot N_{O}$ , where $N_{O}$ denotes the number of sources. This results in increased computation time as the number of sources grows, as illustrated in Figure 24(c). However, the number of candidates per source remains stable. In this experiment, a single target is positioned at the center of the scene at an elevation of 1.5 meters, while sources are randomly sampled on a plane parallel to the ground at an elevation of 70 meters, above all buildings. Recall that variations in the number of candidates per source may occur due to collisions. For a comprehensive view, Figure 24(d) presents results where the total number of samples is kept constant at $10^{6}$ , resulting in $N_{S}=\left\lfloor\frac{10^{6}}{N_{O}}\right\rfloor$ samples per source. In this scenario, the computation time increases at a much slower rate with the number of sources, but the number of candidates found per source decreases rapidly as the number of sources increases.

Higher-Order Diffraction and Diffuse Reflection

The path solver supports only first-order diffraction and does not handle paths that contain both diffraction and diffuse reflection. Enabling support for such combined paths would require the refinement of non-suffix path segments, which is computationally expensive and thus not implemented at present. However, these mixed paths typically transport a negligible amount of energy and therefore have little impact on the accuracy of the computed CIRs.

Regarding high-order diffraction, there are algorithms for determining the diffraction vertices [22]. However, calculating the electric field along these high-order diffraction paths remains an open problem, and as a result, high-order diffraction is not currently supported.

Collisions of Specular Chain Candidates

The path solver does not handle collisions of specular chain candidates. More precisely, when two different specular chains share the same index value in the hash array despite having a different hash due to the modulo operation (23), the path solver discards one of them. This could be avoided by using a collision resolution mechanism, and by dynamically adjusting the hash array size to reduce the probability of collisions. Currently, it is assumed that the hash array size is set to a sufficiently large value to make collisions negligible.

NVIDIA Sionna RT