Terahertz Wireless Data Center: Gaussian Beam or Airy Beam?

To address the aforementioned challenges, in this paper, we investigate Airy beam-based beamforming in a quasi-LoS environments for THz wireless data centers. We propose a novel cascaded geometric and wave channel model (CGWCM) to capture the diffraction characteristics in the blocked scenarios effectively. Furthermore, three Airy beam search schemes are developed to establish communication links under quasi-LoS. Extensive simulation results demonstrate that the proposed Airy beam search schemes outperform conventional Gaussian beam search approaches in quasi-LoS environments, offering a more practical solution for THz wireless data centers. The main contributions are summarized as follows.

• $\bullet$

  We develop a CGWCM for the quasi-LoS scenarios, which achieves a balance between accuracy and complexity. Specifically, we consider the ULA in the transceiver and present the geometrical channel model (GCM) and the wave-based channel model (WCM). In quasi-LoS, GCM exhibits low accuracy because it neglects diffraction, whereas the WCM can accurately capture the diffraction effects caused by the blockages but involves a complex integral form. To establish an accurate and low-complexity model, we develop the CGWCM, incorporating the concept of virtual arrays in blocked regions, integrating GCM and WCM, and transforming the complex integral form into a simple matrix product form. Finally, the calibration of the CGWCM is discussed to enable practical spectral efficiency evaluation in quasi-LoS scenarios.

• $\bullet$

  Motivated by the non-diffraction, self-acceleration, and self-healing characteristics, we investigate the Airy beam-based beamforming design. Three key characteristics, non-diffraction, self-acceleration, and self-healing, are introduced, which are helpful for the transmission of the Airy beams in the blocked scenarios. We further explore the analog precoder of Airy beams by imposing a cubic phase profile on a Gaussian beam, which can be practically implemented by adjusting the phase shifts in THz hybrid beamforming architectures. Then, the digital precoder in the Tx and combiner in the Rx are designed with the Airy beam-based analog precoder and the available channel information.

• $\bullet$

  We propose Airy beam search schemes for the communication link establishment in the quasi-LoS scenarios without prior knowledge of blockage position. We first derive the sampling interval by analyzing the correlation between arbitrary Airy beams and obtain an exhaustive search scheme for optimal beam with high training overhead. To reduce this overhead, a hierarchical focusing-Airy beam search is proposed, which achieves near-optimal beam by first identifies the best focusing beam and then refining the curving coefficient. To further reduce overhead and accelerate link establishment, a low-complexity beam search is proposed, which uses prior position information to steer the beam toward the receiver and obtain a sub-optimal Airy beam.

• $\bullet$

  We evaluate the performance of the CGWCM and proposed beam search schemes based on real channel data measured in a data center. Results demonstrate that the CGWCM remains more accurate than GCM in quasi-LoS scenarios while having a more simplified form than WCM. The proposed hierarchical focusing-Airy beam search scheme outperforms the Gaussian beam search schemes and achieves performance comparable to exhaustive search schemes but with reduced overhead. Meanwhile, the proposed low-complexity beam search scheme can rapidly establish communication links while maintaining good performance with significantly lower overhead.

The remainder of the paper is organized as follows. In Sec. II, the system model of THz UM-MIMO systems is provided. In Sec. III, we introduce the GCM and WGM. By combining these two models, the low-complexity CGWCM model is proposed. Then, we analyze the characteristics and Airy beam-based beamforming design in Sec. IV. Furthermore, in Sec. V, the beam correlation is analyzed and beam search schemes are proposed. After an in-depth analysis and numerical evaluation of the proposed CGWCM and beam search schemes in Sec. VI, the paper is summarized in Sec. VII.

The GCM models electromagnetic wave propagation as straight-line transmission, neglecting wave effects, and is widely adopted in the literature [27]. As shown in Fig. 2(a), the ULAs are deployed at the Tx and Rx, where $y_{ti}$ and $y_{rj}$ denote the vertical coordinates of the $i^{\mathrm{th}}$ transmitting and $j^{\mathrm{th}}$ receiving antennas, respectively, with $i=1,\cdots,N_{t}$ and $j=1,\cdots,N_{r}$. To simplify analysis, the blockage is modeled as a rectangle of width $W$ and length $T=T_{1}+T_{2}$. $D$ denotes the distance between the Tx and Rx and $L$ is the distance from the Tx to the blockage. The four lines extending from the $i^{\mathrm{th}}$ transmitting antenna to the vertices of the blockage are labeled as $l_{1}$, $l_{2}$, $l_{3}$, and $l_{4}$, intersecting the Rx array plane at points $b_{1}$, $b_{2}$, $b_{3}$, and $b_{4}$, which can be expressed as

In GCM, the complex path gain from the $i^{\mathrm{th}}$ transmitting antenna to the $j^{\mathrm{th}}$ receiving antenna is denoted as $\alpha^{ij}_{g}$ and the complex gain of all antenna pairs between the Tx and Rx is computed individually. Since GCM models electromagnetic wave propagation as a straight line, antenna pairs with an unblocked path have a non-zero complex gain, while those with a blocked path have zero gain, meaning the transmitted signal is completely obstructed. In Fig. 2(a), for a fixed $y_{ti}$, the propagation path remains unblocked if $y_{rj}>\text{max}(b_{1},b_{2})$ or $y_{rj}<\text{min}(b_{3},b_{4})$. Otherwise, the propagation path is blocked. Consequently, the quasi-LoS channel response between the $i^{\mathrm{th}}$ transmitting antenna and the $j^{\mathrm{th}}$ receiving antenna can be expressed as

where $D^{ij}$ represents the communication distance from the $i^{\mathrm{th}}$ transmitting antenna to the $j^{\mathrm{th}}$ receiving antenna and $\mathbf{H}_{\mathrm{G}}$ denotes the $N_{r}\times N_{t}$-dimensional geometrical channel matrix. Considering the free spreading loss, the path gain $|\alpha^{ij}_{g}|$ in (3) can be explicitly written as $|\alpha^{ij}_{g}|=\frac{c}{4\pi fD^{ij}}$, where $c$ is the speed of light and $f$ is the carrier frequency.

In exchange for its low computational complexity of $\mathcal{O}(N_{t}N_{r})$, GCM has low accuracy in quasi-LoS conditions as it neglects electromagnetic wave diffraction.

The WCM is derived from the Huygens-Fresnel principle, which accounts for Electromagnetic (EM) wave diffraction and effectively characterizes wave behavior, particularly in the quasi-LoS [23]. According to this principle, each point on the wavefront acts as a secondary source emitting spherical waves, and the superposition of these secondary waves determines the resulting field at the receiver. To describe wavefront propagation on a plane, the Rayleigh–Sommerfeld (RS) integral can be employed [23], enabling the calculation of the electric field at any point based on a known input field. Let the $x_{0}$-plane be defined at $x=x_{0}$, and $\mathbf{E}(x_{0},y_{0})$ denote the electric field at position $(x_{0},y_{0})$. The electric field on the $x_{0}$-plane is represented as a column vector $\mathbf{E}(x_{0},\colon)$, which contains values at all positions $(x_{0},y)$ on the plane. Based on the input field $\mathbf{E}(x_{\mathrm{in}},\colon)$, the electric field at an arbitrary location $(x_{0},y_{0})$ can be expressed as

where $Q_{\mathrm{w}}$ is a scalar parameter, $x_{\mathrm{in}}$ is the propagation coordinate of input plane, $y^{\prime}$ is the transverse coordinate of input electric field, $r=\sqrt{(x_{0}-x_{in})^{2}+(y_{0}-y^{\prime})^{2}}$ and $k=\frac{2\pi}{\lambda}$.

As shown in Fig. 2(b), the $i^{\mathrm{th}}$ transmitting antenna emits a point electric field $\mathbf{E}(0,y_{ti})$ and forms an electric field column vector $\mathbf{E}(L,\colon)$ on the L-plane, which can be expressed as

where $y_{pn}$ is denoted as the transverse coordinate of $n^{\mathrm{th}}$ point on the plane and $r_{n}=\sqrt{L^{2}+(y_{ti}-y_{pn})^{2}}$. However, because of the presence of blockage, a part of the electric field is blocked and cannot transmit continuously. In [23], a binary matrix $\mathbf{B}(x,y)$ is introduced, which takes the value of 0 if the blocker includes point $(x,y)$ and takes 1 otherwise. Therefore, the blocked electric field $\mathbf{E}_{\mathrm{B}}$ on the L-plane can be written as

We iteratively apply the RS integral to the planes within the blockage region, with spacing $\Delta x$. After obtaining the electric field $\mathbf{E}_{\mathrm{B}}(L+W,\colon)$ of $(L+W)$-plane, we can obtain the quasi-LoS channel response between the $i^{\mathrm{th}}$ transmitting antenna and the $j^{\mathrm{th}}$ receiving antenna as

Following the above procedure, the complete channel matrix accounting for the diffraction effect can be obtained.

Although the WCM provides an accurate characterization of EM wave propagation in the presence of blockages, its integral form leads to high computational complexity, particularly when modeling UM-MIMO channels with a large number of elements.

Inspired by both the GCM and WCM, we propose the CGWCM with a simplified form and high accuracy. As shown in Fig. 2(c), we place $M$ virtual arrays at regular intervals of $\Delta x$ within the blockage region to account for diffraction effects based on the Huygens-Fresnel principle, as in the WCM. The $m^{\mathrm{th}}$ virtual array is denoted as $P_{m}$, and the transverse coordinate of $n^{\mathrm{th}}$ antenna in the virtual array is denoted as $y_{pn}$. Each virtual array consists of $N$ antennas with a spacing of $d=\lambda/2$, and the channel between adjacent arrays is modeled using the GCM. The channel between the Tx and $P_{1}$ is represented by $\mathbf{H}_{\mathrm{TP}}\in\mathbb{C}^{N\times N_{t}}$, the channel between $P_{m}$ and $P_{m+1}$ is denoted by $\mathbf{H}_{m(m+1)}\in\mathbb{C}^{N\times N}$ and the channel between $P_{M}$ and the Rx is $\mathbf{H}_{\mathrm{PR}}\in\mathbb{C}^{N_{r}\times N}$. Given the transmitting signal $\mathbf{x}\in\mathbb{C}^{N_{t}\times 1}$, the receiving signal at $P_{1}$ is $\mathbf{y}_{1}=\mathbf{H}_{\mathrm{TP}}\mathbf{x}$. Because of the blockage, the antennas within the blocked region cannot transmit the signal. Let $\mathbf{B}_{\mathrm{m}}$ be the binary matrix for the $m^{\mathrm{th}}$ virtual array $P_{m}$, where entries from row indices $b_{m1}$ to $b_{m2}$ are set to zero, i.e., $\mathbf{B}_{\mathrm{m}}(b_{m1}:b_{m2},\colon)$ = 0, and the remaining entries are one. Accordingly, the blocked received signal at $P_{1}$ is given by

where the dimensional of $\mathbf{B}_{1}$ is the same as $\mathbf{H}_{\mathrm{TP}}$, i.e., $\mathbf{B}_{1}\in\mathbb{C}^{N\times N_{t}}$. Then, $P_{1}$ transmits the signal to the next virtual array and repeats the similar processing. After cascading through the virtual arrays, the final received signal $\mathbf{y}_{b}$ can be expressed as

where $Q_{c}$ is a scalar parameter and $\mathbf{B}_{2}$, $\cdots$, $\mathbf{B}_{\mathrm{M}}$ are $N\times N$ dimensional. Therefore, the quasi-LoS channel matrix under the proposed CGWCM framework can be given by

Compared with WCM, which computes RS integrals individually for each antenna pair with a complexity of $\mathcal{O}(N_{t}N_{r}MN^{2})$, CGWCM uses matrix operations and achieves a more efficient and parallelizable formulation with a total complexity of $\mathcal{O}(MN^{3}+N^{2}N_{t}+NN_{t}N_{r})$, making it well-suited for large-scale THz UM-MIMO systems.

In LoS scenarios, the GCM provides highly accurate channel modeling. However, in quasi-LoS scenarios, CGWCM introduces virtual arrays to account for diffraction, which leads to inaccuracies in the complex path gain. To address this, we propose a calibration method that corrects the path gain of CGWCM by leveraging the accurate gain provided by the GCM in the LoS case. Specifically, we introduce virtual arrays at the same locations as those used in the quasi-LoS CGWCM model and compute the LoS channel using CGWCM, denoted as $\mathbf{H}_{\mathrm{C}}^{\mathrm{LoS}}$. We then compute the LoS channel using the GCM, denoted as $\mathbf{H}_{\mathrm{G}}^{\mathrm{LoS}}$, which serves as a reference for accurate complex path gain. By comparing these two channels, we derive the amplitude and phase calibration parameters as

Then the calibrated LoS and quasi-LoS channel can be given by $\mathbf{H}_{\mathrm{Gc}}^{LoS}=\alpha_{cal}e^{j\phi_{cal}}\mathbf{H}_{\mathrm{G}}^{\mathrm{LoS}}$ and $\mathbf{H}_{\mathrm{Gc}}^{qLoS}=\alpha_{cal}e^{j\phi_{cal}}\mathbf{H}_{\mathrm{G}}^{qLoS}$. This calibration effectively compensates for the gain mismatch introduced by virtual arrays and ensures we can practically evaluate the spectral efficiency in the quasi-LoS scenarios.

As demonstrated in Sec. IV, the phase vector of an arbitrary Airy beam in a ULA system is characterized by three fundamental parameters: the cubic coefficient $a$, the distance $r$ and the angle $\theta$. These parameters collectively determine the Airy beam and can be expressed as

Then, the correlation of two arbitrary Airy beams with distinct cubic coefficients $a_{m}$ and $a_{n}$, focusing on polar coordinates $(r_{m},\theta_{m})$ and $(r_{n},\theta_{n})$, respectively, can be formulated as

where $y_{ti}$ is the vertical coordinate of $i^{\mathrm{th}}$ transmitting antenna defined as in (2).

We consider a symmetric antenna array configuration about the origin with an inter-element spacing of $d$. The antenna positions can be expressed as $y_{ti}=(i-\frac{N_{t}+1}{2})d$ for $i=1,2,\cdots,N_{t}$ and the array spans from $-\frac{N_{t}-1}{2}d$ to $\frac{N_{t}-1}{2}d$. Then, the correlation expression can be rewritten as

The discrete summation can be approximated as a continuous integral when the number of antennas $N_{t}$ is sufficiently large. This approximation, based on the Riemann sum concept, becomes increasingly accurate as the inter-element spacing $d$ decreases relative to the array size, which can be expressed as

It appears challenging to directly derive the sampling scheme for the cubic coefficient of the curving characteristics as well as the distance and angle parameters of focusing characteristics from (V-A). To simplify the analysis, we can focus on each part individually and sequentially.

For Airy beam search, an optimal beam vector is selected from a pre-designed beam codebook through a training procedure between the Tx and Rx. In conventional far field beam searching and near field beam focusing searching, the codewords in their codebook are determined by a unique spatial angle or point location. Similarly, for the Airy beam codebook, we aim to identify a distinctive feature that can uniquely determine an Airy beam codeword. Specifically, since the Airy beam can only exist in near field region, as the analysis in Sec. IV, the Airy beam can be defined through a point location information due to its focusing property and curving coefficient due to its self-accelerating trajectory. Thus, an Airy beam codeword can be uniquely defined by a tuple of parameters, i.e., $(a,r,\theta)$. Each codeword determines a unique three-dimensional point location, specified by the two-dimensional focusing point location and the curving coefficient. The codebook is designed to cover the desired range of these parameters, ensuring comprehensive representation of the target Airy beams within the specified region.

Based on the beam correlation analysis, the beam sample interval can be derived to design the codebook. Specifically, we assign the values $C_{a_{m},a_{n}}(\overline{\alpha})$, $C_{r_{m},r_{n}}(\overline{\beta})$ and $C_{\theta_{m},\theta_{n}}(\overline{\gamma})$ as $\xi_{a}$, $\xi_{r}$ and $\xi_{\theta}$, respectively, and derive the corresponding $\overline{\alpha}_{\xi}$, $\overline{\beta}_{\xi}$ and $\overline{\gamma}_{\xi}$. Subsequently, the sample intervals can be expressed as

where in (30b), $\theta$ is a fixed parameter and in (30c) $u=1,2,\cdots,N_{t}$. We define $[a_{\mathrm{min}},a_{\mathrm{max}}]$ as the potential range of the curving coefficient. Specifically, since the location of the blockage is unknown, the possibility of concave and convex trajectory should be covered in $[a_{\mathrm{min}},a_{\mathrm{max}}]$. Thus, we let $a_{\mathrm{min}}<0$, $a_{\mathrm{max}}>0$ and $|a_{\mathrm{min}}|=|a_{\mathrm{max}}|$. We denote $J$ as the number of sample curving coefficient and the $j^{\mathrm{th}}$ sample value can be expressed as

where $j=1,2,\cdots,J$ and $a_{\mathrm{max}}=a_{\mathrm{min}}+(J-1)s_{a}$. Moreover, we define $r_{\mathrm{max}}$ as the distance between the Tx and Rx, i.e., $r_{\mathrm{max}}=D$. Therefore, in the fixed angle $\theta$, $\frac{\cos^{2}{\theta}}{r_{\mathrm{max}}}$ is lowest. We denote $K$ as the number of sample distance and the $k^{\mathrm{th}}$ sample distance can be derived as

where $k=1,2,\cdots,K$. For the angle sample, we define $\theta_{\mathrm{min}}$ and $\theta_{\mathrm{max}}$ as $-\frac{\pi}{2}$ and $\frac{\pi}{2}$, respectively. We denote $V$ as the number of sample angle and the $v^{\mathrm{th}}$ sample angle can be derived as

where $v=1,2,\cdots,V$ and $u$ determines the sample interval of angle domain. As shown in Fig. 6(a), for the exhaustive beam search, the codebook should contain all the sample points determined by $(a,r,\theta)$ in three-dimensional space. Therefore, the total number of codewords in the exhaustive beam search codebook is given by $JKV$, and the corresponding codebook can be expressed as

where $A,R,\Theta$ are the set of three parameters, respectively.

In the Tx, a codeword is selected from the codebook for transmission in each time slot. The Rx employs an omnidirectional analog combiner and a unit matrix digital combiner to process the received signal. The received signal in $t^{\mathrm{th}}$ time slot with $t^{\mathrm{th}}$ codeword $f_{t}$ can be denoted as

where $\mathbf{W}$ is combiner in the Rx and $\mathbf{H}_{b}$ is the blocked channel. Since the number of codewords in the exhaustive codebook is $JKV$, the total number of time slots is $T_{e}=JKV$. After $T_{e}$ time slots, the receiver can select an optimal beam which can be described as

Finally, an Airy beam vector $\mathbf{f}_{\mathrm{e}}=\mathbf{a}(a_{\mathrm{e}},r_{\mathrm{e}},\theta_{\mathrm{e}})$ can be used to communicate with the Rx in the quasi-LoS scenario.

To reduce the overhead of exhaustive beam search, we propose a hierarchical Focusing-Airy beam search scheme that narrows the spatial search space and performs a two-stage refinement, where the first stage searches over focusing beams to determine distance and angle and the second stage refines the curving coefficient. As shown in Fig. 6(b), firstly, we narrow the scope of the beam search in the spatial domain to reduce the number of codewords. Since the LoS is blocked, the beam search can be initiated in the LoS region within the spatial domain between the Tx and Rx. Suppose that the numbers of antennas in the Tx and Rx are the same and the length of the ULA is $L_{a}$. The polar coordinate of the points located in the LoS range can be expressed as

The number of $\text{point}_{\mathrm{LoS}}$ is denoted as $T_{f}^{1}$. After we determine the scope of beam search in the spatial domain, we propose a two-stage beam search scheme. In the first stage, we set $a$ equal to 0, which means that all the beams in the codebook of the first stage are focusing beams. The codebook for the first stage can be denoted as

After $T_{f}^{1}$ time slots, the optimal focusing beam $\mathbf{f}_{f}=\mathbf{a}(0,r_{f},\theta_{f})$ can be selected from $\mathcal{F}_{f}$. In the second stage, we fix the distance coefficient $r_{f}$ and angle coefficient $\theta_{f}$ and design a codebook $\mathcal{F}_{a}$ to search the Airy beam with different curving coefficient, which can be expressed as

After $T_{f}^{2}$ time slots where $T_{f}^{2}$ is the total number of codewords in $\mathcal{F}_{f}^{2}$, a near-optimal Airy beam vector can be obtained to keep good performance in quasi-LoS communication scenario.

Instead of exhaustively searching through all possible Airy beams in the scenario, the hierarchical focusing-Airy beam search scheme reduces the search space by narrowing down the potential beam directions and conducting the hierarchical search, which significantly reduces overhead and speeds up the search, making it more practical for real-time applications.

To enable fast and reliable link recovery in quasi-LoS scenarios, we propose a low-complexity beam search scheme that leverages prior Rx position knowledge to directly focus the beam on the Rx and refine the curving coefficient with significantly reduced overhead. Specifically, as shown in Fig. 6(c), We choose the focusing point on the curve $\frac{\cos{\theta}}{r}=\frac{1}{D}$. Similar to the hierarchical focusing-Airy beam search, we also choose the point located in the LoS range. We define the angle range$[\theta_{l}^{min},\theta_{l}^{max}]$ where $\theta^{max}_{l}$ and $\theta^{min}_{l}$ are set as $\arctan{\frac{L_{a}}{2D}}$ and $-\arctan{\frac{L_{a}}{2D}}$. We set the curving coefficient $a$ as 0, the angle sample interval as $s_{\theta}$ and design a codebook $\mathcal{F}_{l}^{1}$ for searching an appropriate focusing beam, which can be denoted as