Switch-based Hybrid Beamforming for Wideband MIMO-OFDM Systems with Beam Squint
Mengyuan Ma, Student Member, IEEE, Nhan Thanh Nguyen, Member, IEEE, and Markku Juntti, Fellow, IEEE

This paper was presented in part in Proc. ITG Workshop Smart Antennas 2021 <cit.>.

This work has been supported in part by Academy of Finland under 6G Flagship (Grant No. 318927) and EERA Project (Grant No. 332362).
Mengyuan Ma,  Nhan Thanh Nguyen, and Markku Juntti are with the Centre for Wireless Communications, University of Oulu, 90014 Oulu, Finland (e-mail:{mengyuan.ma, nhan.nguyen, markku.juntti}@oulu.fi).

Switch-based hybrid beamforming (SW-HBF) architecture is a promising implementation to realize massive multiple-input multiple-output (mMIMO) systems because it can reduce the complexity and the power consumption of the phase-shifter-based hybrid beamforming (PS-HBF) scheme. Recent studies reveal that the performance of PS-HBF is detrimentally affected by beam squint existing in the wideband communication systems while the performance of SW-HBF remains unclear. In this paper, we consider the beam squint effect and SW-HBF design in wideband mMIMO orthogonal frequency-division multiplexing systems. We first derive closed-form expressions of the expected array gains and beam squint ratio (BSR) of both SW-HBF and PS-HBF schemes. The results show that compared to the latter, the former suffers from less severe beam squint effect. Motivated by this, we focus on the problem of SW-HBF design to maximize average spectral efficiency. The formulated problem is a non-convex mixed-integer programming which is significantly challenging to solve optimally. To tackle it, we propose low-complexity algorithms based on the tabu search scheme and projected gradient ascend method by judiciously employing the problem structure. Then, given the analog beamformers, the optimal digital ones are obtained with closed-form solutions. Our extensive numerical results demonstrate that the proposed SW-HBF schemes can approach the performance of the exhaustive method with much lower complexity. Besides, the proposed SW-HBF schemes can outperform the PS-HBF schemes in terms of both spectral and energy efficiency in wideband mMIMO system with severe beam squint effect.

Switch-based Hybrid beamforming, wideband, mMIMO, spectral efficiency, energy efficiency, beam squint effect.

§ INTRODUCTION

Wideband communication systems with large utilizable spectrums are promising to meet the ever-increasing demand for ultra-high-speed data rates of future wireless networks <cit.>. The millimeter-wave (mmWave) band (30 GHz – 300 GHz) and terahertz (THz) band (0.1 THz – 10 THz) communications thereby have been considered essential for the 5G, and 6G of wireless communications <cit.>. In mmWave and THz band, the short wavelength of the carrier frequency allows a large number of antennas to be deployed compactly, enabling the implementation of massive multiple-input multiple-output (mMIMO) systems <cit.>. However, the mMIMO system deploying the conventional fully digital beamforming (DBF) requires an excessively large number of power-hungry radio frequency (RF) chains, which causes prohibitively power consumption and hardware costs to the system. Therefore, hybrid beamforming (HBF) architectures have been proposed as crucial enablers for mmWave mMIMO systems. In such architectures, the overall beamformer is divided into a high-dimension analog beamformer, and a low-dimension digital beamformer <cit.>, which significantly reduces the number of RF chains in the transceiver. Moreover, HBF can attain the multiplexing gain of DBF thanks to the channel sparsity of wideband mMIMO systems <cit.>.

§.§ Prior Works

In recent years, tremendous efforts have been devoted to developing efficient HBF schemes. Driven by the tradeoff between spectral efficiency (SE) and energy efficiency (EE), various HBF architectures have emerged, based on the fact that the analog beamformer can be realized by a network of phase shifters (PSs) <cit.>, switches <cit.> or both <cit.>. Specifically, it is proved in <cit.> that the PS-based HBF (PS-HBF) scheme can realize the optimal digital beamformer in narrowband MIMO systems as long as the number of RF chains is not less than twice the number of data streams.  By comparing the performance of PS-HBF and switch-based HBF (SW-HBF) architectures, Méndez-Rial et al. <cit.> show that the latter schemes can achieve a higher EE than the former, especially in large systems, but at the price of a lower SE. In <cit.>, HBF architectures deploying both PSs and switches are investigated. Such architectures can attain a better tradeoff between the SE and EE compared to the purely SW-HBF or PS-HBF schemes.

All the aforementioned works mainly focus on narrowband systems wherein the effect of frequency-selective channels is negligible. In contrast, HBF design for wideband multi-carrier systems is more challenging because the analog beamformer must be shared across all the subcarriers. To tackle the challenge, efficient PS-HBF algorithms are proposed in <cit.> by searching the eigenvector spaces and constructing the common analog beamformer across all frequencies. However, in those HBF strategies, the beam squint issue <cit.> has been ignored, which may cause significant achievable SE loss <cit.>. In conventional narrowband systems, the analog beamformer aligned with the physical direction at the center frequency can achieve the maximum array gain for signals at other frequencies. However, in wideband mMIMO systems, due to frequency-dependent PSs, the beam aligned with the carrier frequency squints to the incident angles of the incoming signals at other frequencies, resulting in non-negligible array gain loss, termed beam squint effect <cit.>. Furthermore, the beam squint effect is severer for systems with wider bandwidth and larger numbers of antennas.

To mitigate the beam squint effect, several approaches have been proposed <cit.>. Specifically, wide beams can compensate for the loss in array gain induced by the beam squint effect, as shown in <cit.>. Alternatively, high-resolution true-time-delayers (TTD) can be embedded into in the RF front-end to solve the beam squint issue because they can generate frequency-dependent phase shifts based on the signal traveling time. By deploying one TTD to control each antenna, the full array gain can be attained <cit.>. To improve the system EE, TTD-aided HBF architectures with a reduced number of TDDs are proposed in <cit.>. However, the algorithms developed in <cit.> heavily depend on the perfect path angles of the channel, which have to be estimated in practice. Furthermore, the TTD at THz band typically exhibits high power consumption, insertion loss, and hardeare complexity <cit.>. In addition, the practical implementation of the PSs for high frequencies is not a trivial task <cit.>. The latency caused by configuring the phase may lead to performance erosion if the channel is changes rapidly <cit.>. Moreover, the myriad of high-resolution PSs and TTDs employed in the mMIMO systems may cause excessively high power consumption and hardware cost <cit.>, and thereby significantly degrade the system EE.

§.§ Contributions

Compared to PS-HBF architectures with or without TTDs, SW-HBF exhibits several advantages. First of all, it requires lower power consumption as well as hardware cost compared to a phase shifter or a TTD <cit.>. For example, the power consumption of a switch is only ... mW [xxx], while those of a phase shifter and a TTD are ...mW [yyy], respectively. Second, employing only switches simplifies the RF chain design and allows fast adjustment to rapid channel variations <cit.>. Moreover, as shown in <cit.>, SW-HBF can alleviate the beam squint effect in wideband mMIMO systems without inducing additional power consumption or increasing the area and cost of the processing circuit as TTDs do. These motivate us to thoroughly investigate the potentials of SW-HBF architectures in wideband mMIMO Orthogonal frequency-division multiplexing (MIMO-OFDM) systems under significant beam squinting. According to our best knowledge, this has not been well explored in the literature. We first evaluate and compare the impacts of the beam squint on the PS-HBF and SW-HBF schemes by studying their array gains. Our analytical results show that the SW-HBF architectures are more robust to beam squint because it does not form highly-directional beams relying on the coherent phase as the PS-HBF techniques do. In this light, we focus on the design of SW-HBF schemes and propose efficient solutions to the analog precoding/combining matrices. Our specific contributions can be summarized as follows:

* First, we evaluate the relationship between the beam squint effect and the system parameters. <cit.> presents the beam squint ratio (BSR) that quantifies the degree of the beam squint effect. Different from <cit.>, we further take the antenna spacing distance into consideration and obtain a succinct closed-form expression of BSR. The results show how serious the beam squint effect is caused by the large number of antennas, antenna spacing distance, and wide bandwidth in wideband mMIMO systems. Furthermore, the findings reveal that spatial oversampling could be an additional method to reduce the beam squint effect. Besides, comparing the expected array gain of the PS-HBF and SW-HBF schemes indicates that the latter is more robust to the beam squint effect.

* We then focus on the design of SW-HBF scheme. We propose several SW-HBF schemes for wideband mMIMO systems by solving the problem of joint designing the precoders and combiners, which is formulated as maximizing the average SE with the consideration of the transmit power budget, rank constraints, and the hardware limit of switches. The problem is a challenging non-convex mixed-integer programming. To make it tractable, we decouple the design of analog and digital beamformer as two sub-problems. We first propose a tabu search (TS) based algorithm and then refine the TS-based algorithm with the aid of good initialization acquired by the projected gradient ascent (PGA) algorithm. Furthermore, we design a low-complexity search algorithm employing the tabu and repeated sampling procedure. Finally, the digital beamformers are obtained with closed-form solutions.

* Finally, extensive numerical simulations are provided to validate the conclusions of the analysis about beam squint effect, and to evaluate the performance of our proposed SW-HBF algorithms. We show the accuracy of the derived closed-form BSR and how the system parameters impact the beam squint effect. Furthermore, the analysis of the expected array gain for PS-HBF and SW-HBF schemes is verified with experiments. Finally, we present numerical simulations to evaluate the convergence, the SE and EE performance, and the complexity of the proposed SW-HBF algorithms. The results demonstrate that under the severe beam squint effect, the proposed SW-HBF schemes outperform the state-of-art PS-HBF algorithm in terms of both the system EE and SE. Moreover, the developed SW-HBF schemes can obtain near-optimal results as the exhaustive search do but with only a few orders of complexity higher than the PS-HBF algorithm, implying their efficiency.

We note that compared to our prior work <cit.>, we present in this work the detailed analysis of beam squint effects for the PS-based array and the switch-based array. Moreover, rather than only focusing on hybrid combining as in [xxx], we herein consider both hybrid precoding and combing design and present more meticulous algorithm design, computational complexity analysis, and extensive numerical discussions. To focus on the analysis and design of SW-HBF schemes, we assume the availability of full channel state information (CSI) at both the Tx and Rx <cit.>. In practice, CSI can be estimated with compressed sensing method <cit.>, and atomic norm minimization approach <cit.>.

§.§ Organization and notation

The rest of paper is organized as follows. The channel model and signal model are presented in Section <ref>. In Section <ref>, we analyze the beam squint effect and the expected array gain of PS-based array and switch-based array. Then, the SW-HBF schemes are proposed in Section <ref>. In Section <ref>, we shows simulation results, and finally, conclusions are drawn in Section <ref>.

Throughout this paper, scalars, vectors, and matrices are denoted by the lower-case, boldface low-case, and boldface upper-case letters, respectively; Superscripts $\left(\cdot\right)^T$, $\left(\cdot\right)^H$, $\left(\cdot\right)^{-1}$, and  $\left(\cdot\right)^{\dagger}$ stand for transpose, conjugate transpose,  matrix inversion, and  Moore–Penrose inverse, respectively; $\left\|\cdot\right\|_p$ and $\left\|\cdot\right\|_F$ denotes the $p$-norm and Frobenius norm, respectively; Furthermore, the trace and expectation operators are represented by $\tr\left(\cdot \right)$ and $\Es\left(\cdot\right)$, respectively; We use $\Ac$ to denote a set and $\left|\Ac \right|$ to denote the cardinality of $\Ac$; $\left|a \right|$, $\left|\ab \right|$ and $\left|\Ab \right|$ indicates the absolute value of $a$, $\ab$ and determinant of the square matrix $\Ab$, respectively; Moreover, the $m$ by $n$ dimensional complex space and binary space are represented by $\Cs^{m\times n}$ and $\Bc^{m\times n}$, respectively; $\Cc\Nc(\boldsymbol{\mu}, \Cb)$ denotes the distribution of a circularly symmetric complex Gaussian (CSCG) random vector variable with mean $\boldsymbol{\mu}$ an covariance matrix $\Cb$; In addition, $\Ib_N$ signifies a $N\times N$ identity matrix; $\mathbf{1}_N$ represents a $N\times 1$ vector whose elements are $1$; Finally, the element of $\Ab$ at $i$-th row and $j$-th column is denoted by $\Ab(i,j)$.

§ SYSTEM MODEL

§.§ Signal Model

SW-HBF transceiver structure in a wideband MIMO-OFDM system

We consider a point-to-point mmWave MIMO-OFDM system illustrated in Fig. <ref>. The transmitter (Tx) equipped with $N_{\rm T}$ antennas simultaneously transmits $N_{\rm{s}}$ data streams to the receiver (Rx) equipped with $N_{\rm R}$ antennas. For ease of exposition and without loss of generality, we assume that Tx and Rx are equipped with the same number of RF chains, denoted by $N_{\rm RF}$, which satisfies $N_{\rm{S}} \leq N_{\rm RF} \ll \min(N_{\rm T},N_{\rm R})$. The HBF architectures equipped at both the Tx and Rx can utilize a disjoint set of switches or PSs to realize analog beamforming. In this work, we focus on the switch-based analog beamformer due to the reasons discussed earlier. The switches can be configured at high speeds as nanosecond or even sub-nanosecond when the coherence time for a carrier frequency of $f_ {\rm c} = 60 {\rm GHz}$ at an $ 800 {\rm Hz}$ Doppler is around 1 ms <cit.>. As a result, the switching rates and settling times of the currently existing RF switches are enough for wideband mmWave systems. Additionally, the insertion loss of witches, which is only about 1dB <cit.>, can be negligible in practical mmWave systems.

Denote by $K$ the number of subcarriers and by $\ssb_k \in \Cs^{N_{\rm S} \times 1}$ the signal vector transmitted via the $k$-th subcarrier, $\Es\left[\ssb_k\ssb^H_k\right]=\Ib_{N_{\rm S}},\; k=1,\cdots,K$.
At the Tx, $\ssb_k$ is first precoded by a digital precoder $\Fb_{\rm BB}[k] \in \Cs ^{N_{\rm RF} \times N_{\rm S}}$ in the frequency domain, which is then transformed into time domain by the inverse fast Fourier transforms (IFFTs), and up-converted to the RF domain using the RF chains. The RF signals are then further precoded by an analog precoder $\Fb_{\rm RF}\in \Bc ^{N_{\rm T} \times N_{\rm RF}}$ with $\Bc=\{0,1 \}$. Note that the  analog precoder $\Fb_{\rm RF}$ is frequency-flat. More specifically, it is used to precode all the $K$ vectors $\{ \Fb_{\rm BB}[k] \ssb_k \}_{k=1}^K$ at all the $K$ subcarriers.

At the Rx, the received signal vector is first combined by a frequency-flat analog combiner $\Wb_{\rm RF}\in \Bc ^{N_{\rm R} \times N_{\rm RF}}$. After performing fast Fourier transforms (FFTs), the combined signal is further combined in frequency domain by a baseband combiner $\Wb_{\rm BB}[k] \in \Cs ^{N_{\rm RF} \times N_{\rm S}}$ for each subcarrier. Let $\Hb_k\in \Cs^{N_{\rm R}\times N_{\rm T}}$ be the channel at the $k$-th subcarrier. The received signal at the $k$-th subcarrier is given as

\begin{equation}\label{eq:received signal}
\yb_k=\Wb_{BB}[k]^H \Wb_{\rm RF}^H  \Hb_k \Fb_{\rm RF}\Fb_{\rm BB}[k] \ssb_k + \Wb_{\rm BB}[k]^H \Wb_{\rm RF}^H \nb_k,
\end{equation}

where $\nb_k\sim \Cc\Nc(\boldsymbol{0}, \sigma_{\rm n}^2\Ib_{N_{\rm R}})$ is the additive white Gaussian noise vector. By letting $\Fb_k = \Fb_{\rm RF}\Fb_{\rm BB}[k]$, $\| \Fb_k \|^2_F \leq P_{\rm b}, \forall k$, with $P_{\rm b}$ being the power budget for each subcarrier, and letting $\Wb_k = \Wb_{\rm RF}\Wb_{\rm BB}[k]$, we can rewrite (<ref>) in a more compact form as

\begin{equation}\label{eq:transmitted signal}
\yb_k= \Wb_k^H  \Hb_k \Fb_k \ssb_k + \Wb_k^H \nb_k.
\end{equation}

§.§ Wideband Channel Model

We adopt the geometric channel model for wideband mMIMO systems, which can be described as sum of the contributions of multiple propagation paths <cit.>. Assuming that uniform linear arrays (ULAs) are employed[In this paper, we consider the ULA for simplicity, but the analysis of the beam squint effect can be easily extended to the uniform plannar array (UPA), which has a similar channel structure as ULA <cit.>. Besides, the proposed SW-HBF schemes are applicable to arbitrary antenna arrays.], the $d$-th tap of the channel at frequency $f$ can be modeled as <cit.>

\begin{equation}\label{eq:d delay channel}
\Hb_f[d]=\sum\limits_{l=1}^{L_{\rm p}} \alpha_l p\left(d T_{\rm s}-\tau_l\right) \ab_{\rm r}\left( \theta_l^{\rm r}, f\right) \ab_{\rm t}^H\left( \theta_l^{\rm t}, f\right),
\end{equation}

where $L_{\rm p}$ is the number of distinct scattering paths, and $\alpha_l \sim \Cc\Nc(0,1)$ and $\tau_l$ are the complex channel gain and the delay of the $l$-th cluster, respectively; $\theta_l^r$ and $\theta_l^t$ represent the angles of arrival/departure (AoA/AoD) of the $l$-th path, respectively; $T_{\rm s}$ is the the sampling period; $p(\tau)$ denotes the pulse-shaping filter for $T_s$-spaced signaling evaluated at $\tau$ seconds <cit.>. In (<ref>), $\ab_{\rm r}\left( \theta_l^{\rm r}, f\right) $ and $\ab_{\rm t}\left( \theta_l^{\rm t}, f\right)$ denote the transmit and receive array response vectors at frequency $f$ for the $l$-th path with AoA $\theta_l^{\rm r}$ and AoD $\theta_l^{\rm t}$, respectively. The array response vector $\ab(\theta,f) \in \{\ab_{\rm r}\left( \theta_l^{\rm r}, f\right), \ab_{\rm t}\left( \theta_l^{\rm t}, f\right)\}$ is given as

\begin{equation}\label{eq:array steer vector for MIMO}
\ab(\theta,f)=\frac{1}{\sqrt{N}}\left[1,e^{-j 2 \pi \Delta \frac{f}{f_{\rm c}}\sin(\theta) },\cdots, e^{-j 2 \pi (N-1)  \Delta \frac{f}{f_{\rm c}}\sin(\theta) }  \right]^T,
\end{equation}

where $N \in \{N_{\rm T}, N_{\rm R} \}$,  $\theta \in \{ \theta_l^{\rm r}, \theta_l^{\rm t} \}$, and $\Delta$ denotes the antenna spacing normalized with the wavelength.
The frequency-domain channel at the $k$-th subcarrier can be expressed as <cit.>

\begin{equation}\label{eq:subcarrier channel}
\Hb_k=\sum_{d=0}^{D-1} \Hb_{f_k}[d] e^{-j\frac{2\pi k}{K}d},
\end{equation}

where $f_{k}$ denotes the frequency of the $k$-th subcarrier given by

\begin{equation}\label{eq:subcarrier difinition}
f_{k}=f_{\rm c}+\left(k-\frac{K+1}{2} \right)\frac{B}{K},\ \forall k,
\end{equation}

with $B$ denoting the system bandwidth.

§ ANALYSIS OF THE BEAM SQUINT EFFECT

Generally, analog beamformers are designed to generate beams oriented toward the physical direction of the channel path components to achieve the maximum array gain, which can compensate for severe path loss <cit.>. Then, digital precoders are designed based on the analog beamformers to achieve spatial multiplexing gain. Therefore, the array gain provided by the analog beamformer has a critical impact on the achievable rate performance <cit.>. In this section, we first obtain a closed-form expression of the beam squint ratio, which is defined to quantify the degree of beam squint effect. To measure the impact of beam squint on the arrray gain of all subcarriers, we further analyze the expected gain for both PS-based array and switch-based array.

§.§ Beam Squint Ratio

The beam squit ratio (BSR) is defined in <cit.> as the ratio between the physical direction deviation and half of the beamwidths for all subcarrier frequencies and physical directions. Specifically, let $N$ denote the number of antennas in a ULA, and let $\Bar{\theta}\triangleq \sin \theta$ and $\xi_k\triangleq \frac{f_k}{f_{\rm c}}$. Because physical direction deviation for the $k$-th subcarrier and the half of the beamwidth can be given by $\left|\left(1-\xi_k\right)\Bar{\theta} \right|$ and $\frac{1}{N\Delta}$ <cit.>, respectively. The BSR is expressed as <cit.>

\begin{equation}\label{eq:BSR definition}
{\rm BSR}\triangleq \frac{1}{2} \int^{1}_{-1} \frac{1}{K} \sum\limits_{k=1}^K \frac{\left|\left(1-\xi_k\right)\Bar{\theta} \right|}{\frac{1}{N\Delta}} d\Bar{\theta}.
\end{equation}

We note that, unlike the definition of BSR in <cit.> for half-wavelength antenna spacing, we take the normalized antenna spacing $\Delta$ into consideration in (<ref>). Furthermore, we prove that $\Delta$ has significant impact on the beam squint effect. Because $N\Delta\left|\left(1-\xi_k\right)\Bar{\theta} \right|>1$ represents that $\xi_k\Bar{\theta}$ has moved out of the beam mainlobe aligned with carrier frequency, ${\rm BSR} >1$ implies the beam squint effect is severe and able to cause serious performance degradation. Nevertheless, the integral in (<ref>) limits the insights of the BSR. In the following proposition, we present a closed-form expression of the BSR.

When the number of subcarriers $K$ is large, the BSR can be approximately given as

\begin{equation}\label{eq:BSR closed-form}
{\rm BSR}\left(N,b,\Delta\right) {\approx} \frac{N b\Delta}{8},
\end{equation}

where $N$, $b$, and $\Delta$ denote the number of antennas, fractional bandwidth, i.e., $b= \frac{B}{f_{\rm c}}$, and normalized antenna spacing, respectively. A larger ${\rm BSR}$ implies more severe beam squint effect.

By dividing both sides of (<ref>) to $f_{\rm{c}}$, we obtain

\begin{align}\label{eq:xi bandwidth}
\xi_k = \frac{f_k}{f_{\rm{c}}}&=1+ \left(\frac{k}{K}-\frac{K+1}{2K} \right) b
\end{align}

where $b= \frac{B}{f_{\rm c}}$. From (<ref>) and (<ref>), we can write

\begin{align}
{\rm BSR}  &=\frac{N\Delta b }{2K}  \sum\limits_{k=1}^K \left|\frac{k}{K}-\frac{K+1}{2K}\right|\int^{1}_{-1} | \theta|d\theta \nonumber \\
&=\frac{N\Delta b }{2} \cdot \frac{1}{K}\sum\limits_{k=1}^K \left|\frac{k}{K}-\frac{K+1}{2K}\right| \notag \\
&\overset{(a)}{\approx} \frac{N\Delta b }{2} \int^{1}_{0} | x-\frac{1}{2}|dx
=\frac{N\Delta b}{8},
\end{align}

where $(a)$ follows the facts that $K \gg 1$ and $\int^{1}_{0} f(x)dx=\displaystyle \lim_{K \to \infty} \frac{1}{K} \sum\limits_{k=1}^K f(\frac{k}{K})$ for a continuous real-valued function $f(\cdot)$ defined on the closed interval $[0, 1]$.

Accuracy of the closed-form BSR expression ($N=128, b=1/14, \Delta=1/2$).

The relationship between the BSR and system parameters .

Note that compared to the BSR defined in <cit.>,  that in (<ref>) is given in a more compact form and takes the normalized antennas spacing $\Delta$ into consideration. In Fig. <ref>, we compare the BSR obtained by the definition (<ref>) and the closed-form expression (<ref>) with the number of subcarriers $K$. It can be observed that as $K$ increases, the BSR in (<ref>) tends to that obtained in (<ref>), which verifies the approximation in Proposition <ref>.

[${\rm BSR}= 0.1,\Delta=1/2,N=22$ .]

[${\rm BSR}= 1,\Delta=1/2,N=224$.]

[${\rm BSR}= 0.5,\Delta=1/4,N=224$. ]

Illustration of beam squint effect in multi-carrier systems using PS-based ULA with $f_{\rm c}=140~{\rm GHz},B=10~{\rm GHz}, K=128$. Each figure presents the normalized array gain of three frequencies, namely the carrier frequency, the maximum frequency, and the minimum frequency, versus $\Bar{\theta}= \sin \theta$.  Figs. (a), (b), and (c) indicate three different BSR cases, i.e., $({\rm BSR}= 0.1,\Delta=1/2,N=22)$, $({\rm BSR}= 1,\Delta=1/2,N=112)$, and $({\rm BSR}= 0.5,\Delta=1/4,N=224)$, respectively.

Moreover, (<ref>) unveils the beam squint effect linearly increases with the number of antennas $N$, fractional bandwidth $b$, and normalized antenna spacing $\Delta$. Fig. <ref> shows how those three parameters intertwine to affect the beam squint. Two groups of lines, i.e. the solid lines and dotted lines are shown in the figure indicating two typical BSR scenarios, including the severe beam squint effect ${\rm BSR}=1$ and the mild one ${\rm BSR}=0.1$. As an example, the normalized array gain for ${\rm BSR}=0.1$ and ${\rm BSR}=1$ are shown in Figs <ref> and <ref>, respectively, with $f_{\rm c}=140~{\rm GHz},B=10~{\rm GHz}, K=128$. Some insights can be drawn based on the results in  Fig. <ref>:

* For fixed antenna spacing, the number of antennas and the system fractional bandwidth can be jointed exploited to tradeoff the beam squint effect. Specifically, when the beam squint effect is same, reducing the number of antennas allow to deploy wider fractional bandwidth, and vice versa.

* Adopting spatial oversampling, i.e., reducing the antenna spacing, can effectively mitigate the beam squint effect. The reason behind this is that the diminish of antenna spacing results in a wider beamwidth, increasing the array gain of the marginal frequencies in the bandwidth. It can be seen from Figs. <ref> and <ref> that the beamwidth is doubled when antenna spacing decreases half times, alleviating the beam squint effect.

Furthermore, based on the Proposition (<ref>), we draw following corollary.

There exists a small positive scalar $\varrho $ with which the beam squint effect can be negligible if $N \leq \frac{8\varrho}{b\Delta}$ is satisfied.

Illustration of beam squint effect in multi-carrier systems using PS-based ULA with $f_{\rm c}=140~{\rm GHz},B=10~{\rm GHz}, K=128,\Delta=1/2$.

Corollary <ref> can be useful to evaluate the beam squint effect of wideband mMIMO systems once $\varrho $ is determined. Fig. <ref> shows the normalized array gain of three BSR cases versus the subcarrier index with $f_{\rm c}=140~{\rm GHz},B=10~{\rm GHz}, K=128,\Delta=1/2$. We can observe that the array gain loss for ${\rm BSR}=0.1$ is marginal, which can also be corroborated from Fig. <ref>. Thus, we postulate that $\varrho=0.1$ is sufficiently small. In other words, as long as the system parameters satisfy $N \leq \frac{0.8}{b\Delta}$, the beam squint effect can be neglected, as demonstrated through simulations in Section <ref>.

§.§ Expected Array Gains of SW-HBF and PS-HBF Architectures

§.§.§ PS-based Array

The array gain achieved by phase shifter-based analog beamforming vector $\tilde{\wb}=\left[e^{j\varpi _1}, \cdots,e^{j\varpi_n},\cdots,e^{ j\varpi _N} \right]^T$ for physical direction $\theta$ at frequency $f$ can be given as <cit.>

\begin{align}\label{eq:func NAG}
g\left( \tilde{\wb}, \theta,f \right)&=\frac{1}{\sqrt{N}}\left| \tilde{\wb}^H \ab_{\rm r}\left( \theta,f \right) \right| \notag
\\
&=\frac{1}{N}\left|\sum\limits_{n=1}^{N}e^{j\left[\varpi_n-2\pi(n-1)\Delta  \frac{f}{f_{\rm c}}\sin \theta \right]}\right|.
\end{align}

For a multi-carrier system, the beamformer $\tilde{\wb}$ is conventionally designed for the center carrier frequency, i.e., $f_{\rm{c}}$. To obtain the highest array gain, the beamformer $\tilde{\wb}$ should align with the AoA $\theta$ of the arriving signal <cit.>. Specifically, the $n$-th phase shift in $\tilde{\wb}$ should be set to

\begin{equation}\label{eq:aligned beamformer}
\varpi_n=2\pi(n-1)\Delta \sin \theta, n=1,\cdots, N.
\end{equation}

With this, we can rewrite (<ref>) as

\begin{align}\label{eq:PS NAG bi-variable}
g\left(\Bar{\theta} , \xi_k \right)&=\frac{1}{N}\left| \sum\limits_{n=1}^{N}e^{j\left[2\pi(n-1)\Delta \left( 1-\xi_k\right) \sin \theta   \right]} \right| \nonumber\\
&=\left|  \frac{\sin\left( N\pi \Delta \left(1- \xi_k\right)\Bar{\theta}  \right)}{N \sin\left( \pi \Delta \left(1- \xi_k\right)\Bar{\theta}  \right)}\right|,
\end{align}

Then, the expected array gain (EAG) over all possible physical directions is

\begin{equation}\label{eq:EAG for PS array}
\Es_{\rm ps}[g\left( \Bar{\theta},\xi_k \right)]= \frac{1}{2K} \sum\limits_{k=1}^K\int^{1}_{-1}  g\left( \Bar{\theta},\xi_k \right) d \Bar{\theta}.
\end{equation}

The EAG of PS-based array $\Es_{\rm ps}[g\left( \Bar{\theta},\xi_k \right)]$ can be approximated as

\begin{equation}\label{eq:approxmated EAG for PS-based array}
\Es_{\rm ps}[g\left( \Bar{\theta},\xi_k \right)]\approx \frac{2}{3\rho} \int_0^{\rho} \left|{\rm sinc}(\pi x)\right|dx + \frac{1}{3},
\end{equation}

where $\rho\triangleq\frac{N\Delta b}{2}$. The EAG given in (<ref>) is monotonically decreasing with $\rho$, and is bounded as

\begin{equation}
\frac{1}{3}\leq \Es_{\rm ps}[g\left( \Bar{\theta},\xi_k \right)]\leq 1,
\end{equation}

The lower bound is attained at $\rho \rightarrow \infty$, while the upper bound is achieved at $\rho \rightarrow 0$.

See Appendix <ref>.

Proposition <ref> reveals that the EAG of PS-based array is monotonically decreasing with the number of antennas $N$, fractional bandwidth $b$ and normalized antenna spacing $\Delta$. This conclusion coincides with Proposition <ref>, which further demonstrates that $N\Delta b$ is an important indicator to evaluate the beam squint effect. Specifically, increasing value of $N\Delta b$ causes more serious EAG loss, leading to more severe beam squint effect.

§.§.§ Switch-based Array

The normalized array gain achieved by switch-based analog beamforming bector $\wb = [w_1, \ldots, w_N], w_n \in \left\{0,1\right\}, n = 1, \ldots, N$  is given by

\begin{align}\label{eq:base matrix for SW-BSR}
&g\left( \wb, \Bar{\theta},\xi_k \right)=\frac{1}{\sqrt{N}}\left| \wb^H \ab_{\rm r}\left( \Bar{\theta},f_k \right) \right| \notag\\
&=\frac{1}{N}\left|\sum\limits_{n=1}^{N} w_n e^{-j2\pi(n-1)\Delta  \xi_k \Bar{\theta} }\right|,
\end{align}

The EAG of the switch-based array can be written as

\begin{align}\label{eq:EAG of switch bases array}
\Es_{\rm sw}[g\left( \wb, \Bar{\theta}\right)] &=  \frac{1}{2K} \sum_{k=1}^K  \int^{1}_{-1} g\left( \wb, \Bar{\theta},\xi_k \right) d\Bar{\theta} \nonumber\\
&=\frac{1}{2K} \int^{1}_{-1} \left\|\Omegab (\Bar{\theta} )\wb \right\|_1 d\Bar{\theta},
\end{align}

where $\Omegab (\Bar{\theta} )$ is a $K\times N$ matrix with $n$-th column given by $\frac{1}{N} \left[ e^{-j2\Delta (n-1)\pi \xi_1  \Bar{\theta}}, e^{-j2\Delta (n-1)\pi \xi_2  \Bar{\theta}}, \ldots, e^{-j2\Delta (n-1)\pi \xi_K  \Bar{\theta}} \right]$. With the Newton polynomial fitting the integration in (<ref>) at three points $(-1,\left\|\Omegab (-1 )\wb\right\|_1 )$, $(0,\left\|\Omegab (0 )\wb \right\|_1 )$ and $(1,\left\|\Omegab (1 )\wb\right\|_1 )$, we obtain

\begin{align}\label{eq:interpolated expected array gain sw}
\Es_{\rm sw}[g\left( \wb, \Bar{\theta} \right)]
&\approx \medmath{\frac{1}{2}\int^{1}_{-1} \left[\left( \left\|\Omegab (-1 )\wb \right\|_1 -\left\|\Omegab (0 )\wb \right\|_1 \right)\Bar{\theta}^2+\left\|\Omegab (0 )\wb \right\|_1 \right] d\Bar{\theta}}  \nonumber\\
&\hspace{-1cm} =\frac{2}{3}\left\|\Omegab (-1 )\wb \right\|_1 +\frac{1}{3N}\left\|\wb \right\|_1 \overset{(i)}{\approx} \frac{1}{3N} \left\|\wb \right\|_1 \leq \frac{1}{3},
\end{align}

where approximation $(i)$ follows

\begin{align}\label{eq:switch EAG term1}
\left\|\Omegab (-1 )\wb \right\|_1&=\frac{1}{KN} \sum\limits_{k=1}^K \left|\sum\limits_{n=1}^{N} w_n e^{j2\Delta\pi\left[1+\left(\frac{k}{K}-\frac{K+1}{2K} \right) b\right](n-1)}\right| \approx 0
\end{align}

by the Law of Large Numbers for $\left\|\wb \right\|_1 \gg 1$.
Unlike the EAG of the PS-based array, that of the switch-based array in (<ref>) and (<ref>) depends on  $\wb$ and achieves the largest array gain of $1/3$ with $\wb=\mathbf{1}_N$, i.e., when each RF chain is connected to all antennas.

[EAG versus bandwidth.]

[EAG versus BSR with $N=256$.]

EAG of the PS-based array and switch-based array, marked as PS and SW in the figure, respectively, versus key factors with $f_{\rm c}=140~{\rm GHz},K=128,\Delta=1/2$. (a) exhibits the EAG versus the system bandwidth and number of antennas. (b) shows the EAG versus the BSR where the number of antennas $N$ is fixed and bandwidth increases.

Figs. <ref> and <ref> show the EAG of the PS-based array and switch-based array versus the system parameters and BSR, respectively. In Fig. <ref>, the carrier frequency and antenna spacing are fixed while the number of antennas and system bandwidth vary. It can be observed that the EAG of the PS-based array decreases with the increases of bandwidth and number of antennas and tends to its lower bound. i.e. $1/3$, verifying the Proposition <ref>. Furthermore, Fig. <ref> demonstrates the EAG is inversely proportional to the BSR, further corroborating the Proposition <ref>. Moreover, the maximum EAG of the switch-based array is presented in Figs. <ref> and <ref> for comparison.

It is observed that different from the EAG of the PS-based array, that of a switch-based array is independent of system bandwidth, indicating the SW-HBF can be less affected by the beam squint effect than PS-HBF. As shown in (<ref>) and (<ref>), the array gain function is obtained by the inner product of the steering vector and the analog beamforming vector. For the PS-based array, the analog beamforming vector has coherent phases to the steering vector, as set in (<ref>). While the coherent combining of the entries of the steering vector leads to high array gains, it also causes severe dependency on the bandwidth. Regarding the switch-based array, the analog beamforming vector is used to select the entries of the steering vector. Combining the entries of the steering vector without non-coherent phases renders the term related to bandwidth close to zero, as shown in (<ref>). As a result, the impact of bandwidth on the array gain is lessened, which implies that the switch-based array can be more robust to the beam squint effect than the PS-based array. Thus, it is expected that SW-HBF can achieve higher SE than PS-HBF for wideband mMIMO system under beam squint effect, which will be further justified by numerical results in Section <ref>.

§ SWITCH-BASED HYBRID BEAMFORMING DESIGN

As analyzed in Section <ref>, the SW-HBF is less affected by the beam squint effect than the PS-HBF. Motivated by this, in this section, we focus on the design of efficient SW-HBF scheme for wideband mMIMO systems.

§.§ Problem Formulation

Considering the Gaussian signalling, based on signal model (<ref>), the SE of the $k$-th subcarrier, denoted by $R_k$, is expressed as<cit.>

\begin{equation}\label{eq:SE for each subcarrier}
R_k=\log_{2}\left|\Ib_{N_{\rm S}}+\frac{1}{\sigma_{\rm n}^{2}} \Wb_k^{\dagger}  \Hb_k \Fb_k \Fb^H_k \Hb^H_k \Wb_k \right|.
\end{equation}

where we recall that $\Fb_k=\Fb_{\rm RF}\Fb_{\rm BB}[k]$ and $\Wb_k=\Wb_{\rm RF}\Wb_{\rm BB}[k]$.
We aim at designing the transmit and receive HBF matrices for all subcarrriers, i.e., $\left\{ \Fb_k\right\} \triangleq \{ \Fb_1, \ldots, \Fb_K\}$ and $\left\{\Wb_k\right\} \triangleq \{\Wb_1, \ldots, \Wb_K \}$, respectively, that maximize the average SE of the system under power constraints of subcarriers and binary constraints of analog beamformers.
The problem can be fomulated as

\begin{align}
\underset{ \{ \Fb_k \}, \{ \Wb_k\}}{\max} \quad  & \frac{1}{K} \sum_{k=1}^{K} R_k\\
\textrm{s.t.} \qquad  & \| \Fb_k \|^2_F \leq P_{\rm b},\ \forall k, \\
&\Fb_{\rm RF}\in \Bc ^{N_{\rm T} \times N_{\rm RF}}, \label{eq: analog precoder binary constraint}\\
&\Wb_{\rm RF}\in \Bc ^{N_{\rm R} \times N_{\rm RF}}, \label{eq: analog combiner binary constraint}\\
&\rank(\Fb_{\rm RF})=N_{\rm RF},\label{eq:precoder rank constr}\\
&\rank(\Wb_{\rm RF})=N_{\rm RF}.\label{eq:combiner rank constr}
\end{align}

The joint optimization of the transmit and receive beamformers, i.e., $\left\{ \Fb_k \right\}$ and $\left\{ \Wb_k \right\}$, in problem (<ref>) is of significantly challenging. It is a non-convex non-smooth mix-integer programming problem due to the logarithmic objective, the full-rank constraints, and the binary constraints of the analog beamformers. To overcome the challenges, we first decouple the problem into two subproblems: transmit and receive beamforming design. We first solve the transmit beamformers $\left\{ \Fb_k \right\}$ with given receive beamformers. Then, giving $\left\{ \Fb_k \right\}$, the receive beamformers $\left\{ \Wb_k \right\}$ are obtained. The procedure are elaborated next.

§.§ Transmit HBF Design

Given combiners $\left\{ \Wb_k \right\}$, the analog and digital precoders, i.e., $\Fb_{\rm RF}$ and $\left\{ \Fb_{\rm BB}[k]\right\}$, can be obtained by solving the following problem

\begin{align}
\max\limits_{\Fb_{\rm RF}, \left\{ \Fb_{\rm BB}[k]\right\} } & \frac{1}{K} \sum_{k=1}^K R_k\\
\text { s.t. } \quad & \eqref{eq: analog precoder binary constraint}~ \text{and} ~\eqref{eq:precoder rank constr}~ \text{are satisfied}, \nonumber \\
& \| \Fb_{\rm RF}\Fb_{\rm BB}[k] \|^2_F \leq P_{\rm b},\ \forall k.  \label{eq:coupled analog and digital beamformer}
\end{align}

To tackle the highly coupled variables $\Fb_{\rm RF}$ and $\Fb_{\rm BB}[k]$ in (<ref>), we first introduce auxiliary variables $\hat{\Fb}_{\rm BB}[k]\in \Cs^{N_{\rm RF} \times N_{\rm s}},\forall k $ such that $\Fb_{\rm BB}[k]\triangleq (\Fb_{\rm RF}^H\Fb_{\rm RF})^{-\frac{1}{2}}\hat{\Fb}_{\rm BB}[k] $. Thus, we have $\| \Fb_{\rm RF}\Fb_{\rm BB}[k] \|^2_F = \| \hat{\Fb}_{\rm BB}[k] \|^2_F$, and problem (<ref>) can be recast as

\begin{align}
\max\limits_{\Fb_{\rm RF}, \left\{ \hat{\Fb}_{\rm BB}[k] \right\} } & \frac{1}{K} \sum_{k=1}^{K} R_k\\
\text { s.t. } \quad & \eqref{eq: analog precoder binary constraint}~ \text{and} ~\eqref{eq:precoder rank constr}~ \text{are satisfied}, \nonumber \\
& \| \hat{\Fb}_{\rm BB}[k] \|^2_F \leq P_{\rm b},\ \forall k,
\end{align}

where $R_k$ is re-expressed as

\begin{align}\label{eq:equivalent SE for Tx}
R_k &= \log_{2}\left|\Ib+\frac{1}{\sigma_{\rm n}^{2}} \Wb_k^{\dagger} \Hb_k \Fb_{\rm RF}\Fb_{\rm BB}[k] \Fb_{\rm BB}^H[k]\Fb_{\rm RF}^H \Hb^H_k \Wb_k  \right| \notag \\
&=\log_{2}\left|\Ib+\frac{1}{\sigma_{\rm n}^{2}} \Hb_{\rm eff}[k]\hat{\Fb}_{\rm BB}[k] \hat{\Fb}_{\rm BB}[k]^H  \Hb_{\rm eff}[k]^H  \right| ,
\end{align}

with the definition of $\hat{\Fb}_{\rm BB}[k]$, $\Hb_{\rm eff}[k]\triangleq \Pb_{\rm t}[k]\Hb_k \Fb_{\rm RF}(\Fb_{\rm RF}^H\Fb_{\rm RF})^{-\frac{1}{2}}$ and $\Pb_{\rm t}[k]\triangleq \Wb_k \Wb_k^{\dagger}$. The binary constraint (<ref>) makes problem (<ref>) a challenging non-convex mixed-integer programming. To address the challenge, we decouple the designs of the analog and digital precoders next.

§.§.§ Digital Precoding Design

We first design the digital precoder with given analog one. When $\Fb_{\rm RF}$ is fixed, the optimal solution for $\hat{\Fb}_{\rm BB}[k]$ is given as <cit.>

\begin{equation}\label{eq: sub digital precoder}
\hat{\Fb}_{\rm BB}[k]=\Vb_k\Gammab_k^{\frac{1}{2}},
\end{equation}

where $\Vb_k\in \Cs^{N_{\rm RF}\times N_{\rm S}}$ is the matrix whose columns are $N_{\rm S}$ right-singular vectors corresponding to the $N_{\rm S}$ largest singular values of $\Hb_{\rm eff}[k]$, $\Gammab_k \in \Cs^{N_{\rm S}\times N_{\rm S}}$ is the diagonal matrix whose diagonal elements are the power allocated to $N_{\rm S}$ data streams at the $k$-th subcarrier, with $\tr(\Gammab_k)=P_{\rm b}$. As a result, with the given $\Fb_{\rm RF}$, the optimal digital precoder for $k$-th subcarrier is

\begin{equation}\label{eq:digital precoder}
\Fb_{\rm BB}[k]=(\Fb_{\rm RF}^H\Fb_{\rm RF})^{-\frac{1}{2}}\Vb_k\Gammab^{\frac{1}{2}}.
\end{equation}

§.§.§ Analog Precoding Design

From (<ref>), it is observe that the objective function of problem (<ref>) is a very complicated function of $\Fb_{\rm RF}$. In the following analysis, we derive an upper bound of the objective function (<ref>), which is more trivial to solve. Specifically, for moderate and high signal-to-noise-ratio (SNR) regime, it is asymptotically optimal to adopt the equal power allocation for all streams in each subcarrier <cit.>. In this regime, we have $\Gammab_k=\gamma\Ib_{N_{\rm S}}$ with $\gamma=P_{\rm b}/N_{\rm S}$. Further, by defining $\tilde{\Ib}_{N_{\rm RF}}\triangleq \left[\begin{array}{rr}
\Ib_{N_{\rm S}} & \mathbf{0} \\
\mathbf{0} & \mathbf{0}
\end{array}\right]$ as a $N_{\rm RF}\times N_{\rm RF}$ matrix, we obtain

\begin{equation}\label{eq:unitary matrix equality}
\Vb_k\Vb_k^H=\tilde{\Vb}_k \tilde{\Ib}_{N_{\rm RF}} \tilde{\Vb}_k ^H,
\end{equation}

where $\tilde{\Vb}_k \in \Cs^{N_{\rm RF}\times N_{\rm RF}}$ is an introduced unitary matrix. With (<ref>), (<ref>) and (<ref>), we obtain

\begin{align}
&\frac{1}{K}\sum_{k=1}^{K} R_k \nonumber\\
&=  \frac{1}{K}\sum_{k=1}^{K} \log_{2}\left|\Ib+\frac{\gamma}{\sigma_{\rm n}^2} \tilde{\Vb}_k^H \Hb_{\rm eff}[k]^H \Hb_{\rm eff}[k] \tilde{\Vb}_k \tilde{\Ib}_{N_{\rm RF}}  \right| \notag \\
& \overset{(a)} \leq  \frac{1}{K}\sum_{k=1}^{K} \log_{2}\left| \Ib + \frac{\gamma}{\sigma_{\rm n}^2} \tilde{\Vb}_k^H \Hb_{\rm eff}[k]^H \Hb_{\rm eff}[k] \tilde{\Vb}_k \Ib_{N_{\rm RF}} \right| \label{eq:SE relaxation}\\
&\overset{(b)}=  \frac{1}{K} \sum_{k=1}^{K} \log_{2}\left| \Ib + \frac{\gamma}{\sigma_{\rm n}^2} \Fb_{\rm RF}^{\dagger} \Hb_k^H \Pb_{\rm t}[k] \Hb_k \Fb_{\rm RF}\right|, \label{eq:SE upper bound}
\end{align}

where inequality $(a)$ follows the fact that $\tilde{\Vb}_k^H\Hb_{\rm eff}[k]^H \Hb_{\rm eff}[k]\tilde{\Vb}_k $ is a semi-definite positive matrix, and equality $(b)$ is because $\Pb_{\rm t}[k]$ can be written as $\Pb_{\rm t}[k]=\Wb_k(\Wb_k^H\Wb_k)^{-1}\Wb_k^H$ which is a Hermitian idempotent matrix.  Note that the relaxation in (<ref>) is tight when $N_s=N_{RF}$ <cit.> under which the multiplexing gain is maximized. Now, by defining $\Hb_{\rm t}[k]\triangleq \Hb_k^H\Pb_{\rm t}[k] \Hb_k$, the problem of designing the analog precoder $\Fb_{\rm RF}$ is reformulated as

\begin{align}\label{eq:analog precoder design}
\max\limits_{\Fb_{\rm RF}} \quad & f\left( \Fb_{\rm RF}\right) \triangleq  \frac{1}{K} \sum_{k=1}^{K} \log_{2}\left| \Ib + \frac{\gamma}{\sigma_{\rm n}^2} \Fb_{\rm RF}^{\dagger} \Hb_{\rm t}[k] \Fb_{\rm RF}\right|\\
\text { s.t. } \quad & \eqref{eq: analog precoder binary constraint}~ \text{and} ~\eqref{eq:precoder rank constr}~ \text{are satisfied}. \nonumber
\end{align}

Problem (<ref>) is still challenging due to the binary and rank constraints of $\Fb_{\rm RF}$. Unfortunately, there is no existing solution to it in the literature. A relaxation of the binary variables may cause significant performance loss and cannot guarantee the rank constraint. Although the optimal solution can be obtained via an exhaustive search, the resultant algorithm requires a prohibitive complexity. For such challenging mixed-integer programming, the tabu search (TS) approach <cit.> can find a near-optimal solution with reduced complexity compared to the exhaustive search. It is widely used for challenging integer problems such as symbol detection in massive MIMO systems <cit.>.

Generally, for a maximization problem in which the objective function is defined as $f(\cdot)$  with a feasible region $\Uc$, the TS procedure starts with initializing a point $\ub_0 \in \Uc$ as the best candidate $\ub_{\rm b}$. Then, in each iteration, the neighbors of the best candidate within the feasible region are examined to find the best neighbor $\ub_{\rm nb}$ that is not in the tabu list but yields the highest objective value. Following that, the best candidate is updated to $\ub_{\rm nb}$, which is then pushed into the tabu list. By iteratively searching the neighbors of the best candidate in the feasible region, the TS algorithm keeps updating $\ub_{\rm b}$ until the stopping criteria are met. During the process, a tabu list records the visited points to avoid cycling. Finally, the TS algorithm outputs the best solution $\ub_{\rm b}^*$, which is found among the best candidates over iterations. With iterative neighbor search and maintaining the tabu list, TS procedure is guaranteed to converge to a local optimum with reasonable complexity <cit.>. Furthermore, it has been well justified in the literature that as long as the number of iterations is sufficiently large, the TS procedure ensures obtaining a near-optimal solution <cit.>. Therefore, we propose a TS-aided analog precoding design to solve the problem (<ref>).

TS-based Analog Precoding Design
TS is a searching-based approach where most complexity comes from evaluating the objective function. Particularly, the computation of $f\left( \Fb_{\rm RF}\right)$ in (<ref>) involves the Moore—Penrose inverse of $\Fb_{\rm{RF}}$, which has high complexity. This can be tackled by exploiting the QR decomposition.

Specifically, let  $\Fb_{\rm RF}=\Vb_{\rm RF}\Db_{\rm t}$ be the QR decomposition of $\Fb_{\rm RF}$, where $\Vb_{\rm RF} \in \Cs^{N_{\rm T}\times N_{\rm RF}}$ has orthogonal column vectors, i.e., $\Vb_{\rm RF}^H\Vb_{\rm RF}=\Ib_{N_{\rm RF}}$, and $\Db_{\rm t}\in \Cs^{N_{\rm RF}\times N_{\rm RF}}$ is an upper triangular matrix with real non-zero diagonal elements. Because $\Fb_{\rm RF}=\Vb_{\rm RF}\Db_{\rm t}$ with $\Db_{\rm t}$ being a invertible matrix,  the objective function in (<ref>) can be rewritten as

\begin{align}
f\left( \Fb_{\rm RF}\right)
&=\frac{1}{K} \sum_{k=1}^{K} \log_{2}\left| \Ib + \frac{\gamma}{\sigma_{\rm n}^2} \Vb_{\rm RF}^H \Hb_{\rm t}[k]\Vb_{\rm RF}\right| \label{eq:O1 a},\\
&=\medmath{\frac{1}{K} \sum_{k=1}^{K} \sum_{i=1}^{N_{\rm S}} \log_2\left(1+ \frac{\gamma}{\sigma_{\rm n}^2} \lambda_{\Vb_{\rm RF}^H \Hb_{\rm t}[k] \Vb_{\rm RF},i}\right)} ,\label{eq:O1 b}
\end{align}

where $\lambda_{\Vb_{\rm RF}^H \Hb_{\rm t}[k] \Vb_{\rm RF},i}$ is the $i$-th eigenvalue of $\Vb_{\rm RF}^H \Hb_{\rm t}[k] \Vb_{\rm RF}$[Based on channel model (<ref>), one can see that $\rank(\Hb_k)=L_{\rm p},\forall k$ where $L_{\rm p}$ is the number of effective propagation paths. We assume $N_{\rm S} \leq L_{\rm p} \leq N_{\rm RF}$ because the number of RF chains is typically greater than or equal to the maximum number of channel paths in practice to ensure the hybrid beamforming performance <cit.>. ]. While the computational complexity of (<ref>) is approximately $\Oc\left(N_{\rm T}^3 \right)$, computation in (<ref>) requires approximate $\Oc\left(N_{\rm S}N_{\rm T}^2 \right)$ FLOPs with the power method <cit.> obtaining the $N_{\rm S}$ largest eigenvalues of $\Vb_{\rm RF}^H \Hb_{\rm t}[k] \Vb_{\rm RF}$. Using (<ref>) as the objective function, we propose the TS-based analog precoding scheme, which is summarized in Algorithm <ref>.

TS-aided Analog Precoding Design

Set $i \gets 0$. Generate initial candidate $\Fb_{\rm{RF}}^{(0)}$, push it to the tabu list $(\Lc_{\rm ts} \gets \Lc_{\rm ts} \cup \Fb_{\rm{RF}}^{(0)})$, and set it as the best solution $(\Fb_{\rm{RF}}^{\star} \gets \Fb_{\rm{RF}}^{(0)})$

$i \leq N_{\rm{iter}}^{\rm max}$ and not converge

Find neighbor set $\Nc(\Fb_{\rm{RF}}^{(i)})$ of the candidate $\Fb_{\rm{RF}}^{(i)}$

Find the best neighbor $\tilde{\Fb}_{\rm RF}^{(i)\star}$ in $\Nc(\Fb_{\rm{RF}}^{(i)})$ but not in the tabu list $\Lc_{\rm ts}$ i.e., $\tilde{\Fb}_{\rm RF}^{(i)\star}=\argmax\limits_{\tilde{\Fb}_{\rm RF}^{(i)} \in \Nc(\Fb_{\rm{RF}}^{(i)}) \setminus  \Lc_{\rm ts}} f(\tilde{\Fb}_{\rm RF}^{(i)})$.

$ f(\tilde{\Fb}_{\rm RF}^{(i)\star}) > f(\Fb_{\rm{RF}}^{\star}) $
Update the best solution: $\Fb_{\rm{RF}}^{\star} \gets \tilde{\Fb}_{\rm RF}^{(i)\star}$.

Push $\tilde{\Fb}_{\rm RF}^{(i)\star}$ to the tabu list: $\Lc_{\rm ts} \gets \Lc_{\rm ts} \cup \tilde{\Fb}_{\rm RF}^{(i)\star}$. Remove the first candidate in $\Lc_{\rm ts}$ if it is full.

Set $\tilde{\Fb}_{\rm RF}^{(i)\star}$ as the candidate for the next iteration, i.e., $\Fb_{\rm{RF}}^{(i+1)} \gets \tilde{\Fb}_{\rm RF}^{(i)\star}$.

$i \gets i+1$

In Algorithm <ref>, the initial candidate for the analog precoder, denoted as $\Fb_{\rm{RF}}^{(0)}$, is initialized in step 1 and aput to the tabu list $\Lc_{\rm ts}$. It is also set as the best solution. Then, in steps 2–10, the iterative neighbor search procedure is conducted. Specifically, in iteration $i$, the neighbor set of the candidate $\Fb_{\rm{RF}}^{(i)}$ is found in step 3. Here, the neighbor set of a candidate $\Fb_{\rm RF}^{(i)}$ is denoted as $\Nc(\Fb_{\rm{RF}}^{(i)})$, which contains analog beamforming matrices that have one element difference to $\Fb_{\rm RF}^{(i)}$ while satisfying the full-rank constraints (<ref>). Then in step 4, the neighbor set is examined to find the best neighbor, denoted as $\tilde{\Fb}_{\rm RF}^{(i)\star}$, which results in the largest objective value and is not in the tabu list. Then, the solution is updated to the best neighbor if a larger objective value is found, as in steps 5–7. The best neighbor is then pushed into the tabu list in step 8. Finally, the candidate is updated for the next iteration as in step 9. This iterative procedure is terminated once the objective value converges, or when the number of iterations exceeds a predefined threshold $N_{\rm{iter}}^{\rm max}$.

The convergence of this algorithm is guaranteed because over iterations, the objective values form a non-decreasing sequence, as in steps 5–7. With sufficiently large number of iterations, the TS-based analog precoding algorithm can attain a near-optimal solution. However, one potential drawback of Algorithm <ref> is that the large number of iterations could result in excessively high complexity. In each iteration of Algorithm <ref>, $|\Nc(\Fb_{\rm RF}^{(i)})| = 2^{N_{\rm T} N_{\rm RF}}$ neighbors (except for some in the tabu list) need to be examined to find the best one. This could cause prohibitively high complexity for mmWave massive MIMO systems, where $N_{\rm T}$ and $N_{\rm{iter}}$ are typically large. To overcome those limitations, we propose an improved TS-based algorithm with a good initial solution obtained based on the projected gradient ascent (PGA) method.

PGA-TS-aided Analog Precoding Design

The proposed PGA-TS scheme consists of two major steps. We first apply the PGA method to find a relaxed continuous initial solution in the closed interval $[0,1]$. As a result, the elements of the solution are either $0$, $1$, or values between $0$ and $1$. Therefore, we can treat those non-integer elements as erroneous entries. Then, we focus on optimizing those erroneous ones in the second step. With fewer entries required to be optimized than the original ones, the search space of the TS algorithm is significantly reduced, which yields faster convergence and lower complexity.

Specifically, to obtain the good initial solution for the TS scheme, the binary and rank constraints in (<ref>) are relaxed, i.e.,

\begin{align}
\max\limits_{\Fb_{\rm RF}\in \Fc} \quad & f\left( \Fb_{\rm RF}\right)
\end{align}

where $\Fc\triangleq\{\Fb_{\rm RF}| \Fb_{\rm RF}(i,j)\in [0,1], \forall i,j \} $, $\Fb_{\rm RF}(i,j)$ denotes the element of $\Fb_{\rm RF}$ at $i$-th row and $j$-th column. Problem (<ref>) can be efficiently solved via the PGA method summarized in Algorithm <ref>. Denote the output of the PGA algorithm by $\Fb_{\rm RF}^{\rm pga}$ and the solution in $i$-th iteration by $\Fb_{\rm RF}^{(i)}$. In each iteration, $\Fb_{\rm RF}^{\rm pga}$ is set as $\Fb_{\rm RF}^{(i)}$ which is updated as $\Fb_{\rm RF}^{(i)} \leftarrow [\Fb_{\rm RF}^{(i-1)}+\mu\nabla_{\Fb_{\rm RF}^{(i-1)}}f(\Fb_{\rm RF}^{(i-1)})]_{\Fc}$. Here, $[\cdot]_{\Fc}$ denotes the projection into $\Fc$, and $\mu$ is the update step-size set as the decreasing value with iterations. $\nabla_{\Fb_{\rm RF}^{(i)}}f(\Fb_{\rm RF}^{(i)})$ denotes the gradient of $f(\Fb_{\rm RF}^{(i)})$ with respect to $\Fb_{\rm RF}$ at $i$-th iteration. With $\Fb_{\rm RF}^{\dagger}=\left(\Fb_{\rm RF}^H\Fb_{\rm RF} \right)^{-1}\Fb_{\rm RF}$, $f(\Fb_{\rm RF})$ can be written as

\begin{align}
f(\Fb_{\rm RF})&=\frac{1}{K} \sum_{k=1}^{K} \log_{2}\left| \Ib + \frac{\gamma}{\sigma_{\rm n}^2} \left(\Fb_{\rm RF}^H\Fb_{\rm RF} \right)^{-1}\Fb_{\rm RF}^H \Hb_{\rm t}[k] \Fb_{\rm RF}\right| \nonumber \\
&=\frac{1}{K} \sum_{k=1}^{K} \log_{2}\left|\Fb_{\rm RF}^H \Ab_k\Fb_{\rm RF}\right|- \log_{2}\left|\Fb_{\rm RF}^H \Fb_{\rm RF}\right|
\end{align}

where $\Ab_k\triangleq\Ib+\frac{\gamma}{\sigma_{\rm n}^2} \Hb_{\rm t}[k] $. Thus, $\nabla_{\Fb_{\rm RF}^{(i)}}f(\Fb_{\rm RF}^{(i)})$ is given as

\begin{equation}\label{eq:PGA gradient}
\setlength{\arraycolsep}{1.5pt}
\medmath{
\nabla_{\Fb_{\rm RF}}f(\Fb_{\rm RF})=\frac{2}{K} \sum_{k=1}^{K} \Ab_k\Fb_{\rm RF}\left(\Fb_{\rm RF}^H\Ab_k\Fb_{\rm RF}\right)^{-1}
-2\Fb_{\rm RF}\left(\Fb_{\rm RF}^H\Fb_{\rm RF}\right)^{-1},
\end{equation}

Note that we omit the superscript $i$ in (<ref>) for simplicity. The iterations continue until the procedure satisfies the stopping criteria $\frac{\left\|\nabla_{\Fb_{\rm RF}^{(i)}}f(\Fb_{\rm RF}^{(i)}) \right\|_F}{\left\| \nabla_{\Fb_{\rm RF}^{(0)}}f(\Fb_{\rm RF}^{(0)}) \right\|_F}\leq \epsilon$ where $\epsilon$ is a small positve scalar.  With proper diminishing step-size, the PGA algorithm is guaranteed to converge to a stationary point of the non-convex problem (<ref>) <cit.>.

PGA algorithm for problem (<ref>)

$\Fb_{\rm RF}^{\rm pga}$
Initialize $c=3,\Fb_{\rm RF}^{(0)}\in \Fc$, set $i=0$ and $\Fb_{\rm RF}^{\rm pga}=\Fb_{\rm RF}^{(i)}$
stopping criteria is satisfied
$\Fb_{\rm RF}^{(i)} \leftarrow [\Fb_{\rm RF}^{(i-1)}+\mu\nabla_{\Fb_{\rm RF}^{(i-1)}}f(\Fb_{\rm RF}^{(i-1)})]_{\Fc}$, and set $\Fb_{\rm RF}^{\rm pga}$ as $ \Fb_{\rm RF}^{(i)}$

Since the elements of $\Fb_{\rm RF}^{\rm pga}$ are continuous in $[0,1]$, their values could be $0$, $1$, or non-integer in $(0,1)$. We treat those non-integer elements as erroneous entries that must be optimized by the TS algorithm. As such, the search space of the TS algorithm can be considerably reduced. Furthermore, it is of high probability that those erroneous entries very close to $0$ or $1$ can be approximated as the nearest integer without incurring significant performance degradation. To realize the approximation, we define the following refinement function, i.e.,

\begin{align}\label{eq:definition of refinement function}
r\left(x\right)\triangleq \begin{cases}
\delta & \text{ if } \left|x-\delta \right|\leq \varepsilon,  \\
x & {\rm otherwise},
\end{cases},
\end{align}

where $\delta\in \{0,1\}$ and $\varepsilon$ is a threshold rounding an entry to a binary value with an error upper bound $\varepsilon$. Here, $\varepsilon$ is a small positive scalar. By applying the refinement function $r(\cdot)$ to the entries of $\Fb_{\rm RF}^{\rm pga}$, we obtain a refined solution $\tilde{\Fb}_{\rm RF}^{\rm pga}$, i.e., $\tilde{\Fb}_{\rm RF}^{\rm pga}(i,j)=r\left(\Fb_{\rm RF}^{\rm pga}(i,j)\right), \forall i,j$. Thus, the number of non-integer elements of $\tilde{\Fb}_{\rm RF}^{\rm pga}$ is smaller than that of $\Fb_{\rm RF}^{\rm pga}$, which further reduces the search space of the TS algorithm.

For example, to find an effective analog precoder for a system with $N_{\rm T}=8$ antennas and $N_{\rm RF}=4$ RF chains in Tx, the search space for the TS algorithm is $2^{32}$. Equation (<ref>) presents the simulation results in one experiment. Running the PGA algorithm and then applying the refinement function ($\varepsilon=0.1$) yield $\Fb_{\rm RF}^{\rm pga}$ and $\tilde{\Fb}_{\rm RF}^{\rm pga}$, respectively. There are only seven non-integer elements in $\Fb_{\rm RF}^{\rm pga}$ based on which the search space for the TS algorithm is reduced to $2^{7}$. Furthermore, with four non-integer entries left in $\tilde{\Fb}_{\rm RF}^{\rm pga}$, the search space for the TS algorithm is only $2^{4}$, which has decreased $2^{28}$ times compared to the initial case. Therefore, with the refined PGA-aided initial solution, the search space of the TS algorithm can be considerably shrunken, leading to a significant computational complexity reduction and faster convergence.

\begin{equation}\label{eq:example Fpga}
\setlength{\arraycolsep}{1.5pt}
\medmath{\Fb_{\rm RF}^{\rm pga}=  \begin{bmatrix}
&0.99	&1	&0	&0\\
&0	&0	&1	&1\\
&0.79	&1	&0	&0\\
&0.28	&0	&1	&1\\
&0.44	&1	&0.05	&0\\
&0.97	&0	&1	&1\\
&0	&1	&0.17	&0\\
&1	&0	&1	&1
\end{bmatrix},
\tilde{\Fb}_{\rm RF}^{\rm pga}=  \begin{bmatrix}
&{\color{red} 1}	&1	&0	&0\\
&0	&0	&1	&1\\
&0.79	&1	&0	&0\\
&0.28	&0	&1	&1\\
&0.44	&1	&{\color{red}0}	&0\\
&{\color{red} 1}	&0	&1	&1\\
&0	&1	&0.17	&0\\
&1	&0	&1	&1
\end{bmatrix},}
\end{equation}

The modified TS algorithm with the refined PGA-aided initial solution, termed the PGA-aided TS (PGA-TS) algorithm, is summarized in Algorithm <ref>. Generally, the difference between the TS and PGA-TS algorithm lies in the initialization and the search space. The PGA-TS algorithm first obtains a reliable initial solution $\tilde{\Fb}_{\rm RF}^{\rm pga}$ in the continuous interval $[0,1]$ by applying the PGA algorithm and the refinement. To facilitate the recovery of the analog precoder, the indexes of non-integer entries of $\tilde{\Fb}_{\rm RF}^{\rm pga}$ can be stored into the set $\Sc_I$. Then, the PGA-TS algorithm performs the TS procedure in the reduced search space $\Bc^{|\Sc_I|}$ which can be significantly smaller than $\Bc^{N_{\rm T}N_{\rm RF}}$. Therefore, the PGA-TS algorithm can find the solution to the analog precoder with reduced complexity and faster convergence. Let $\qb^{(i)} \in \Bc^{|\Sc_I|}$ be the candidate of the TS procedure in the $i$-th iteration. Based on $\Sc_I$, the entries of $\qb^{(i)}$ can be inserted back the corresponding non-integer position in $\tilde{\Fb}_{\rm RF}^{\rm pga}$ to recover the analog precoder, denoted by $\Fb_{\rm RF}^{(i)}$. Other steps of the TS procedure performs similar to that of the TS algorithm.

PGA-TS for designing analog precoder

$\Fb_{\rm RF}^*$
Run the PGA algorithm to obtain $\Fb_{\rm RF}^{\rm pga}$, and
apply the refinement function $m(\cdot)$ to obtain $\tilde{\Fb}_{\rm RF}^{\rm pga}$.
Record the index of those non-integer elements of $\tilde{\Fb}_{\rm RF}^{\rm pga}$ into $\Sc_I$
Set $i \gets 0$, and generate the initial candidate $\qb^{(0)} \in \Bc^{|\Sc_I|}$, then recover the analog precoder $\Fb_{\rm RF}^{(0)}$ based on $\qb^{(0)}$ and $\Sc_I$
Set $\qb^{(0)}$ as the best candidate $\qb^*$ and $\Fb_{\rm RF}^{(0)}$ accordingly as the best solution $\Fb_{\rm RF}^{*}$, i.e., $\qb^* \gets \qb^{(0)}, \Fb_{\rm RF}^{*} \gets \Fb_{\rm RF}^{(0)} $, and then push $\qb^{(0)}$ to the tabu list $\Lc_{\rm ts} \gets \Lc_{\rm ts} \cup \qb^{(0)}$
$i \leq N_{\rm iter}^{\rm max}\; \&$  not converge
Find neighbor set $\Nc(\qb^{(i)})$ of the candidate $\qb^{(i)}$
Recover the analog precoder of each neighbor in $\Nc(\qb^{(i)})$ with the aid of $\Sc_I$, and cache the analog precoder into $\Fc_{\qb}^{(i)}$ if it has full rank
Find the best neighbor $\tilde{\qb}^{(i)\star}$ in $\Nc(\qb^{(i)})$ but not in the tabu list $\Lc_{\rm ts}$ . The according analog precoder of the best neighbor $\tilde{\qb}^{(i)\star}$, denoted as $\tilde{\Fb}_{\rm RF}^{(i)\star}$ achieves the largest objective value $f(\tilde{\Fb}_{\rm RF}^{(i)\star})$ among $\Fc_{\qb}^{(i)}$

$ f(\tilde{\Fb}_{\rm RF}^{(i)\star}) > f(\Fb_{\rm{RF}}^{\star}) $
Update the best solution: $\Fb_{\rm{RF}}^{\star} \gets \tilde{\Fb}_{\rm RF}^{(i)\star}$

Push $\tilde{\qb}^{(i)\star}$ to the tabu list: $\Lc_{\rm ts} \gets \Lc_{\rm ts} \cup \tilde{\qb}^{(i)\star}$. Remove the first candidate in $\Lc_{\rm ts}$ if it is full

Set $\tilde{\qb}^{(i)\star}$ as the candidate for the next iteration, i.e., $\qb^{(i+1)} \gets \tilde{\qb}^{(i)\star}$

$i \gets i+1$

Dimension of feasible space for TS procedure versus the number of antennas ($\varepsilon=0.1,N_{\rm T}=N_{\rm R},N_{\rm S}=N_{\rm RF}=4,K=64,f_{\rm c}=140~{\rm GHz},B=10~{\rm GHz}$).

We note that compared to the TS algorithm, the PGA-TS algorithm has significantly reduced complexity while guaranteeing to find a reliable solution to the analog precoder. Even though the efficiency of the PGA-TS algorithm depends on the $|\Sc_I|$ which is not controllable, the efficiency of the PGA-TS algorithm is manifested by the extensive simulations in Section <ref>. Moreover, we evaluate the dimension reduction of the feasible space of the TS procedure resulted by applying PGA algorithm and the refinement function with Monte Carlo method, as shown in Fig. <ref>. The red solid and dashed lines in the figure represents the original dimension of the feasible space of the TS procedure and the half of it, i.e., $N_{\rm T}N_{\rm RF}$ and $N_{\rm T}N_{\rm RF}/2$, respectively. It can be observed that $|\Sc_I|$ empirically satisfies $|\Sc_I| < N_{\rm T}N_{\rm RF}/2$, verifying the efficacy of the PGA-TS algorithm in shrinking the search space. Furthermore, the $|\Sc_I|$-dimension vector $\qb$ saves memory for the tabu list $\Lc_{\rm ts} $ compared to $N_{\rm T}N_{\rm RF}$-dimension $\Fb_{\rm RF}$.

Although the TS procedure of Algorithm <ref> is implemented in a reduced search space, it involves calculating the objective value $|\Sc_I|$ times in each iteration, which can lead to a high computational load. Therefore, we propose replacing the neighbor search with random sampling in each iteration to enhance computational efficiency further. Based on the tabu list and random sampling with the refined PGA-aided initial solution, we obtain the PGA-TbRs algorithm summarized in Algorithm <ref>.

PGA-TbRs Algorithm
Specifically, with the refined PGA-aided solution $\tilde{\Fb}_{\rm RF}^{\rm pga}$, the non-integer entries of $\tilde{\Fb}_{\rm RF}^{\rm pga}$ is first extracted to form a basis candidate $\qb_b$. In the $i$-th iteration, based on $\qb_b$, the statistical rounding method is employed to generate a sample $\qb_s^{(i)}$ which is pushed into the tabu list if not appeared in the previous iterations. The best candidate $\qb_s^{*}$ is updated if $\qb_s^{(i)}$ achieves a higher objective value. As such, the objective function is only calculated once in each iteration leading to significant complexity reduction.

As is shown in Algorithm <ref>, the first step is the same as that in Algorithm <ref>. The non-integer elements of $\tilde{\Fb}_{\rm RF}^{\rm pga}$ are then extracted to form the basis candidate $\qb_b$ in step 2. The best solution $\Fb_{\rm RF}^*$ is initialized by rounding entries of $\tilde{\Fb}_{\rm RF}^{\rm pga}$ to the nearest integer, i.e. $0$ or $1$. In steps 3-14, repeated sampling based on $\qb_b$ is conducted to find a better solution. Specifically, in iteration $i$, a new sample $\qb_s^{(i)} \in \Bc^{|\Sc_I|}$ is generated by the stochastic rounding method based on $\qb_b$ in step 5. Then, its corresponding analog precoder $\Fb_{\rm RF}^{(i)}$ is recovered with the aid of $\Sc_I$ in step 6. As is required in problem (<ref>), steps 5-6 are repeated until a full-rank analog precoder $\Fb_{\rm RF}^{(i)}$ is found or reaching the maximum sampling times $N_{\rm smp}^{\rm max}$. In steps 8-10, the sampled solution $\Fb_{\rm RF}^{(i)}$ is then updated to the best solution $\Fb_{\rm RF}^*$ if it achieves a higher objective value. Then, the sample $\qb_s^{(i)}$ is pushed into the tabu list to avoid being sampled again, as in steps 8-11. This repeated sampling procedure is terminated once the objective value converges or when the number of iterations exceeds the predefined threshold $N_{\rm iter}^{\rm max}$.

We note that the repeated sampling and storing samples into the tabu list in iterations increases the likelihood of finding a satisfactory local optimum. Therefore, the PGA-TbRs algorithm is potential to achieve comparable performance to the PGA-TS algorithm as long as the tabu list is long enough. Furthermore, for the small-scale problem, the PGA-TaRs algorithm could be highly efficient due to the relatively smaller feasible space, as is verified in Section <ref>.

PGA-TbRs Algorithm

$\Fb_{\rm RF}^*$
Run the PGA algorithm to obtain $\Fb_{\rm RF}^{\rm pga}$, and
apply the refinement function $r(\cdot)$ to obtain $\tilde{\Fb}_{\rm RF}^{\rm pga}$.
Record the index of those non-integer elements of $\tilde{\Fb}_{\rm RF}^{\rm pga}$ into $\Sc_I$
Extract the non-integer elements of $\tilde{\Fb}_{\rm RF}^{\rm pga}$ to form the basis candidate $\qb_b$; Initialize $\Fb_{\rm RF}^*$ by rounding entries of $\tilde{\Fb}_{\rm RF}^{\rm pga}$ to the nearest binary value; Set $i \gets 0$ and $ \Lc_{\rm ts} =\varnothing$
$i \leq N_{\rm iter}^{\rm max} \; \&$  not converge
$j \gets 0$
$\Fb_{\rm RF}^{(i)}$ is full-rank or $j\geq N_{\rm smp}^{\rm max}$
Based on $\qb_b$, generate a sample $\qb_s^{(i)} \in \Bc^{|\Sc_I|}$ that is not in the tabu list $\Lc_{\rm ts}$   Recover its corresponding analog precoder $\Fb_{\rm RF}^{(i)}$ with the aid of $\Sc_I$
$j \gets j+1$

$f(\Fb_{\rm RF}^{(i)}) > f(\Fb_{\rm RF}^*) $
$\Fb_{\rm RF}^* \gets \Fb_{\rm RF}^{(i)}$

$\Lc_{\rm ts} \gets \Lc_{\rm ts} \cup \qb_s^{(i)}$
Remove the first candidate in $\Lc_{\rm ts}$ if it is full

$i \gets i+1$

§.§ Receive HBF Design

For a given analog combiner $\Wb_{\rm RF}$, the optimal digital combiner of each subcarrier is the MMSE solution <cit.>, i.e.,

\begin{equation}\label{eq:optimal baseband combiner}
\Wb_{\rm BB}[k]=\left(\Jb_k\Jb^H_k+\sigma_{\rm n}^2\Wb_{\rm RF}^H\Wb_{\rm RF}\right)^{-1}\Jb_k,
\end{equation}

where $\Jb_k\triangleq \Wb_{\rm RF}^H \Hb_k\Fb_k$. Then, the problem of analog combiner design can be formulated as

\begin{align}\label{eq:analog combiner problem}
\max\limits_{\Wb_{\rm RF}} & \quad \frac{1}{K}\sum_{k=1}^{K}\log_2\left|\Ib+\frac{1}{\sigma_{\rm n}^2}\Wb_{\rm RF}^{\dagger}\Tb_k\Wb_{\rm RF}\right| \\
\text { s.t. }  &\quad \eqref{eq: analog combiner binary constraint}~ \text{and} ~\eqref{eq:combiner rank constr}~ \text{are satisfied}, \nonumber
\end{align}

where $\Tb_k\triangleq \Hb_k\Fb_k\Fb^H_k\Hb^H_k$. It is observed that problem (<ref>) and (<ref>) are similar. Therefore, the algorithms proposed in Section <ref> are applicable to solve problem (<ref>).

§.§ Overall SW-HBF Transceiver Design and Complexity Analysis

Our proposed SW-HBF scheme is summarized in Algorithm <ref>, in which, $\Fb_{\rm RF},\{\Fb_{\rm BB}\}[k]$ and $\Wb_{\rm RF},\{\Wb_{\rm BB}\}[k]$ are solved sequentially. Here, while the analog beamforming matrices $\Fb_{\rm RF}$ and $\Wb_{\rm RF}$ are obtained with Algorithm <ref>, <ref> or <ref>, the digital ones $\{\Fb_{\rm BB}\}[k]$ and $\{\Wb_{\rm BB}\}[k]$ admits closed-form solutions, as shown in (<ref>) and (<ref>). Algorithm <ref> requires the initialization of the combiner for each subcarrier $\Wb_k$, which can be set to the optimal digital combiners, i.e.,  $\Wb_k=\Ub_k$ whose columns are the right singular vectors of $\Hb_k$. Correspondingly, $\Pb_{\rm t}[k]$ is an identity matrix, i.e., $\Pb_{\rm t}[k]=\Ib_{N_{\rm R}}$. In this sense, the precoder is found by maximizing the mutual information <cit.>.

Proposed SW-HBF scheme

$\Hb[k],\Pb_{\rm t}[k]= \Wb_k \Wb_k^{\dagger},\forall k$,$\gamma=\frac{P_{\rm b}}{N_{\rm S}}, \sigma^2_{\rm n}$
$\Fb_{\rm RF},\Wb_{\rm RF} $, $\Fb_{\rm BB}[k],\Wb_{\rm BB}[k], \forall k$
$\Hb_{\rm t}[k]=\Hb_k^H \Pb_{\rm t}[k]\Hb_k,\forall k$
Obtain the analog precoder $\Fb_{\rm RF}$ using Algorithm <ref>, <ref> or <ref> with $\{\Hb_{\rm t}[k]\}_{k=1}^K$
Obtain the digital precoder $\Fb_{\rm BB}[k]$ for each subcarrier according to (<ref>)
$\Fb_k = \Fb_{\rm RF}\Fb_{\rm BB}[k], \Tb_k=\Hb_k\Fb_k\Fb^H_k\Hb^H_k, \forall k$
Obtain the analog combiner $\Wb_{\rm RF}$ using Algorithm <ref>, <ref> or <ref> with $\{\Tb_k\}_{k=1}^K$
Obtain the digital combiner $\Wb_{\rm BB}[k]$ for each subcarrier according to (<ref>)

We first evaluate the complexity of obtaining the analog beamformer with Algorithm <ref>, <ref> and <ref>. The complexity to compute the objective function $ f\left( \cdot \right) $ is approximate

\begin{equation}\label{eq:cmplx obj}
\Oc\left( KN_{\rm RF}\left(2N_{\rm T}^2+4N_{\rm T} N_{\rm RF}+N_{\rm RF}^2\right) \right) \approx \Oc\left( 2KN_{\rm RF}N_{\rm T}^2\right),
\end{equation}

which are mainly caused by the QR decomposition, eigenvalue decomposition (EVD) and matrix multiplications. Equation (<ref>) is due to the fact that the number of RF chains in mmWave (THz) transceiver in practice is much less than the number of antennas.
Since Algorithm <ref> requires to compute the objective function about $N_{\rm T}N_{\rm RF}$ times in each iteration, its complexity is approximate $\Oc \left( 2N_{\rm iter}^{\rm max} KN_{\rm RF}^2N_{\rm T}^3 \right)$ FLOPs. In Algorithm <ref>, the complexity of the PGA algorithm is approximate

\begin{align}
&\Oc\left(N_{\rm pga}KN_{\rm RF}\left(4N_{\rm T}^2+4N_{\rm T} N_{\rm RF}+N_{\rm RF}^2\right) \right) \nonumber \\ &\approx  \Oc\left(4N_{\rm pga}KN_{\rm RF}N_{\rm T}^2 \right)
\end{align}

due to the matrix inversion and $N_{\rm RF} \ll N_{\rm T}$. Here, $N_{\rm pga}$ denotes the number of iterations when the PGA algorithm converges. Considering that Algorithm <ref> requires to perform $|\Sc_I|$ times computation of the objective function in each iteration, we use the empirically upper bound $N_{\rm T}N_{\rm RF}/2$ as $|\Sc_I|$ to evaluate the complexity of Algorithm <ref>, which is given by

\begin{align}
&\Oc\left(4N_{\rm pga}KN_{\rm RF}N_{\rm T}^2 \right)+N_{\rm iter}^{\rm max} \cdot N_{\rm T}N_{\rm RF}/2 \cdot \Oc\left( 2KN_{\rm RF}N_{\rm T}^2\right)\nonumber \\
&=\Oc\left(4N_{\rm pga}KN_{\rm RF}N_{\rm T}^2+ N_{\rm iter}^{\rm max}KN_{\rm RF}^2N_{\rm T}^3 \right)
\end{align}

As Algorithm <ref> uses repeated sampling instead of the neighbor search procedure, the objective function is computed only once. Thus, its complexity is approximate

\begin{align}
&\Oc\left(4N_{\rm pga}KN_{\rm RF}N_{\rm T}^2 \right)+N_{\rm iter}^{\rm max} \cdot \Oc\left( 2KN_{\rm RF}N_{\rm T}^2\right)\nonumber \\
&=\Oc\left((4N_{\rm pga} +2N_{\rm iter}^{\rm max}) KN_{\rm RF}N_{\rm T}^2 \right)
\end{align}

By replacing $N_{\rm T}$ with $N_{\rm R}$, the complexity of obtaining analog combiner $\Wb_{\rm RF}$ can be derived similarly.
The complexity of the other steps of Algorithm <ref> is given by $\Oc\left(2KN_{\rm T}N_{\rm R}\left(N_{\rm T} +N_{\rm R}\right) \right)$, which is mainly caused by matrix multiplications

Finally, the complexities of the Algorithm <ref>, Algorithm <ref> and Algorithm <ref> are summarized in TABLE <ref> where $N$ could be $N_{\rm T}$ or $N_{\rm R}$ for obtaining the analog precoder or combiner, respectively.

Computational complexity of proposed algorithms.

Algorithm     Computational complexity

TS     $\Oc \left( 2N_{\rm iter}^{\rm max} KN_{\rm RF}^2N_{\rm T}^3 \right)$

PGA-TS      $\Oc\left(4N_{\rm pga}KN_{\rm RF}N_{\rm T}^2+ N_{\rm iter}^{\rm max}KN_{\rm RF}^2N_{\rm T}^3 \right) $

PGA-TbRs       $\Oc\left((4N_{\rm pga} +2N_{\rm iter}^{\rm max}) KN_{\rm RF}N_{\rm T}^2 \right)  $

Other steps of Algorithm <ref>     $\Oc\left(2KN_{\rm T}N_{\rm R}\left(N_{\rm T} +N_{\rm R}\right) \right)$

§ SIMULATION RESULTS

In this section, we present the numerical results to evaluate the performance and computational complexity of the proposed SW-HBF schemes. In the simulations, unless otherwise stated, we set $L_p=4, d_s=\frac{\lambda_c}{2},f_c=140~{\rm GHz}$, and $N_s=N_{RF}=4$. The AoA/AoDs are uniformly distributed over $\left[-\frac{\pi}{2},\frac{\pi}{2}\right]$, and the pulse shaping filter is modeled by the raised cosine function<cit.>

\begin{align*}
\begin{cases}
\frac{\pi}{4} \operatorname{sinc}\left(\frac{1}{2 \beta}\right), &\text{if~} t=\pm \frac{T_{\mathrm{s}}}{2 \beta} \\
\operatorname{sinc}\left(\frac{t}{T_{\mathrm{s}}}\right) \frac{\cos \left(\frac{\pi \beta t}{T_{\mathrm{s}}}\right)}{1-\left(\frac{2 \beta t}{T_{\mathrm{s}}}\right)^{2}}, &\text {otherwise},
\end{cases}
\nthis \label{eq:pulse shaping}
\end{align*}

with $T_{\rm s}$ being the sampling period and the roll-off factor $\beta=1$. The path delay is uniformly distributed in $[0,(D-1)T_{\rm s}]$ where $D$ is the cyclic prefix length, given by $D=K/4$ according to 802.11ad specification. The SNR is defined as SNR$\triangleq \frac{P_{\rm b}}{K\sigma_{\rm n}^2}$. All reported results are averaged over $10^3$ channel realizations.

§.§ Convergence Analysis

Fig. <ref> shows the convergence of the proposed Algorithm <ref>, <ref>, <ref>, and <ref>, respectively. The Algorithm <ref>, <ref> and <ref> are terminated if the best objective value keeps the same over $N_{\rm end}$ consecutive iterations. We set $N_{\rm end}=10$ and set the maximum number of iterations as $N_{\rm iter}^{\rm max}=200$. It can be observed in Fig. <ref> that the proposed algorithms maintain good convergence properties. Specifically, it takes $50$ iterations for the PGA algorithm to converge. As expected, the PGA-TS algorithm converges significantly faster than the TS algorithm. It takes $40$ iterations for the former to converge, while the latter converges after $90$ iterations. Furthermore, we can observe that the PGA-TbRs algorithm requires only $15$ iterations to converge, which saves many iterations than the PGA-TS algorithm and the TS algorithm, revealing its high computational efficiency.

convergence of proposed algorithms with $N_{\rm T}=N_{\rm R}=64,K=128,B=10~{\rm GHz}$, and ${\rm SNR}=10~{\rm dB}$.

§.§ Performance of Proposed SW-HBF Schemes

In this section, we show the SE and EE of the proposed SW-HBF schemes. The EE is defined as <cit.>

\begin{equation}
{\rm EE}=\frac{{\rm SE}}{P_{\rm total}},
\end{equation}

where $P_{\rm total}$ is the total power consumption of the transceiver.  The total power consumption of DBF, PS-HBF, and SW-HBF schemes are given as

\begin{align}\label{eq:architecture power model}
P^{\rm DBF}_{\rm total}&=(N_{\rm T}+N_{\rm R})(P_{\rm LNA}+P_{\rm RF}+2P_{\rm ADC}),\\
P^{\rm PS-HBF}_{\rm total}& =(N_{\rm T}+N_{\rm R})(P_{\rm LNA}+P_{\rm SP}+N_{\rm RF}P_{\rm PS}) \nonumber\\
&\qquad\qquad +2N_{\rm RF}(P_{\rm RF}+P_{\rm C}+2P_{\rm ADC}), \\
P^{\rm SW-HBF}_{\rm total} &=(N_{\rm T}+N_{\rm R})(P_{\rm LNA}+P_{\rm SP}+N_{\rm RF}P_{\rm SW}) \nonumber\\
&\qquad\qquad +2N_{\rm RF}(P_{\rm RF}+P_{\rm C}+2P_{\rm ADC}),
\end{align}

respectively. Here, $P_{\rm RF}$ represents the power consumption of an RF chain, which can be given as <cit.>

\begin{equation}\label{eq:RF power model}
P_{\rm RF}=P_{\rm M}+P_{\rm LO}+P_{\rm LPF}+P_{\rm BBamp} ,
\end{equation}

where $P_{\rm M}$, $P_{\rm LO}$, $P_{\rm LPF}$, and $P_{\rm BBamp}$ are the power consumption of the mixer, the local oscillator, the low pass filter, and the baseband amplifier, respectively. The assumptions on the component power consumptions are given in the Table <ref> <cit.>.

Moreover, we include the performance of the state-of-the-art PS-HBF algorithm in <cit.> with low-resolution PSs for comparison. The discrete phase values of each PS are within the predefined phase set, i.e., $\Pc=\left\{ e^{j\frac{2\pi n}{2^{b_{\rm ps}}}}\left|n=0,1,\cdots,2^{b_{\rm ps}}-1 \right.\right\}$ where $b_{\rm ps}$ controls the PS resolution. For $b_{\rm ps}\geq 4$, the performance loss induced by the phase quantization is negligible <cit.>. Furthermore, the power consumption of a 4-bit, 2-bit and 1-bit PS are set to $P_{\rm PS}^{\rm 4-bit}=40~{\rm mW}$, $P_{\rm PS}^{\rm 2-bit}=20~{\rm mW}$ and $P_{\rm PS}^{\rm 1-bit}=10~{\rm mW}$, respectively <cit.>.

Power consumption of each device

Device      Notation       Value

Low Noise Amplifier (LNA)     $P_{\rm LNA}$    39mW

Splitter             $P_{\rm SP}$    19.5mW

Combiner            $P_{\rm C}$        19.5mW

Switch              $P_{\rm SW}$     5mW

Mixer                $P_{\rm M}$         19mW

Local oscillator      $P_{\rm LO}$       5mW

Low pass filter       $P_{\rm LPF}$      14mW

Base-band amplifier     $P_{\rm BBamp}$     5mW

ADC                  $P_{\rm ADC}$       240mW

The SE and EE versus the BSR are shown in Figs. <ref> and <ref>, respectively, with $N_{\rm T}=N_{\rm R}=256,K=64,{\rm SNR}=10~{\rm dB}$. With the fixed number of antennas, carrier frequency, and antenna spacing, the BSR is linear with the system bandwidth according to (<ref>). To intuitively show the impact of the beam squint, we display the SE and EE versus the BSR instead of the system bandwidth. Furthermore, simulations of the optimal DBF and the SW-HBF scheme in which the analog beamforming matrices are randomly generated are also provided for comparison. The following observations are noted:

* It can be observed that when the ${\rm BSR}=0.1$, the PS-HBF scheme with 4-bit PSs can achieve 96.7% of the optimal SE obtained by the DBF. Such results validate the postulation in Section <ref> that $\varrho=0.1$ is sufficiently small to overlook the beam squint effect. Further, based on the Corollary <ref>, the beam squint effect is negligible if the system parameters satisfy $N \leq \frac{0.8}{b\Delta}$.

* Both the performance of the SW-HBF and PS-HBF schemes degrades with the BSR. However, it can be observed from Fig. <ref> that compared to the SE of the PS-HBF technique, that of the proposed SW-HBF schemes decreases with the BSR more slowly, verifying that the SW-HBF architectures are more robust to beam squint effect.

* The proposed SW-HBF schemes can achieve higher EE than the PS-HBF scheme with low-resolution PSs. The proposed PGA-TS algorithm performs better than the TS and PGA-TbRs algorithms. Furthermore, we can observe that the SW-HBF scheme using the PGA-TS algorithm outperforms the PS-HBF scheme with 4-bit PSs when ${\rm BSR}\geq 0.6$ and outperform the PS-HBF scheme with 1-bit PSs when ${\rm BSR}\geq 0.3$. Regarding the TS and PGA-TbRs algorithms, they can achieve higher SE than the PS-HBF scheme with 4-bit PSs when ${\rm BSR}\geq 1.3$ and the PS-HBF scheme with 1-bit PSs when ${\rm BSR}\geq 0.5$. To summarize, SW-HBF techniques can outperform the PS-HBF schemes in terms of both the SE and EE under the beam squint effect. For example, Table <ref> shows the performance of the SW-HBF schemes and PS-HBF techniques at ${\rm BSR}=0.7$. It can be observed that the SW-HBF with the PGA-TS algorithm outperforms the  4-bit PS-HBF scheme with 1bps/Hz SE improvement and 280% EE enhancement. Compared to the PS-HBF scheme with 1-bit PSs, the SW-HBF with PGA-TbRs algorithm can achieve 3bps/Hz SE improvement and 138% EE enhancement.

[Average SE versus BSR.]

[Average EE versus BSR.]

Average SE and EE versus BSR with $N_{\rm T}=N_{\rm R}=256,K=64
$, and ${\rm SNR}=10~{\rm dB}$.

Performance comparison for ${\rm BSR}=0.7$.

2*Metric   3c|SW-HBF    3c|PS-HBF

TS     PGA-TS   PGA-TbRs     1 bit    2 bits     4 bits

SE (bps/Hz)    23.39    27.82     24.28    21.25    25.25    26.08

EE (Gbits/J)    0.53   0.62     0.54    0.39    0.34    0.22

Fig. <ref>  shows the SE and EE versus the SNR with $N_{\rm T}=N_{\rm R}=256,K=64$, and $B=10~{\rm GHz}$, and ${\rm BSR}=1.14$. It can be observed that the SW-HBF schemes achieve higher SE than the PS-HBF techniques in the high SNR regime. Specifically, the SW-HBF scheme with the PGA-TS algorithm achieves 0.79bps/Hz SE improvement compared to the 4-bit PS-HBF scheme for ${\rm  SNR}= 5~ {\rm dB}$. Furthermore, the SE enhancement becomes larger when SNR increases. The same behavior is seen for the TS and PGA-TbRs algorithms for ${\rm  SNR}\geq 13~ {\rm dB}$. Such results are because the TS, PGA-TS, and PGA-TbRs algorithms obtain the analog beamformer by directly optimizing the objective function $f(\cdot)$,i.e., the average SE for $N_{\rm S}=N_{\rm RF}$, while the PS-HBF scheme aims at maximizing the upper bound of the average SE. Moreover, it can be observed that the PGA-TbRs algorithm achieves comparable SE to the TS algorithm, manifesting its efficacy in designing the analog beamformers. In addition, compared to Fig. <ref> where the SW-HBF scheme with the PGA-TS algorithm outperforms the 4-bit PS-HBF scheme for ${\rm BSR}\geq 0.6$ in a ${\rm SNR}=10~{\rm dB}$ case,  it is for ${\rm BSR}=1.14$ in a ${\rm  SNR}= 5~ {\rm dB}$ scenario here that SW-HBF schemes achieve higher SE than the PS-HBF techniques. Such results reveal that the SW-HBF scheme is potential to perform better than the PS-HBF technique in lower SNR scenarios when the beam squint effect is more severe.

Average SE versus SNR with $N_{\rm T}=N_{\rm R}=256,K=64$, and $B=10~{\rm GHz}$.

Fig. <ref> compares the performance of proposed SW-HBF schemes with proposed algorithms to that using the optimal exhaustive search (ES) method with $N_{\rm T}=N_{\rm R}=8,N_{\rm S}=N_{\rm RF}=2, K=128$ and ${\rm BSR}=0.04$. It can be seen that the proposed algorithms achieve comparable performance as the SW-HBF scheme with the ES method, which verifies their effectiveness in designing analog beamformers. However, the TS algorithm performs better than the PGA-TS and PGA-TbRs algorithms. This is because the small search space increases the chance of finding a good local optimum by simple random initialization. In contrast, instead of improving the SE, the PGA-aided initialization limits the set of potential solutions, resulting in a marginal performance loss of the PGA-TS and PGA-TbRs algorithms. In other words, the PGA-aided initialization only benefits the performance in systems with sizable arrays.

Achievable SE of SW-HBF schemes versus SNR with $N_{\rm T}=N_{\rm R}=8,N_{\rm S}=N_{\rm RF}=2$, and $K=128$.

[Complexity versus number of antennas ($N_{\rm T}=N_{\rm R},N_{\rm S}=N_{\rm RF}=4,K=128$).]

[Complexity versus number of RF chains ($N_{\rm T}=N_{\rm R}=256,N_{\rm S}=N_{\rm RF},K=128$).]

[Complexity versus number of subcarriers ($N_{\rm T}=N_{\rm R}=256,N_{\rm S}=N_{\rm RF}=4$).]

Comparison of computational complexity.

§.§ Complexity Analysis

Fig. <ref> demonstrates the complexities of the proposed SW-HBF schemes in comparison with the optimal ES method and the PS-HBF algorithm. As shown in Fig. <ref>, the proposed algorithms have much lower complexity than the optimal ES method. Furthermore, the complexity of proposed algorithms slowly increases with the number of antennas. As anticipated, the PGA-TbRs algorithm has the lowest complexity among the proposed algorithms while the complexity of the TS algorithm is the highest one, as seen in Figs. <ref>–<ref>. This is because the search space for PGA-TS algorithms is considerably reduced compared to the TS algorithm. In the PGA-TbRs algorithm, the repeated sampling is utilized rather than the neighbor search method, resulting in much lower complexity than the PGA-TS algorithm. Moreover, we can observe from Fig. <ref> that the complexity of the PS-HBF algorithm slightly increases with the number of RF chain. This is because the complexity of obtaining the analog beamformer in the PS-HBF algorithm is approximate $\Oc(N_{\rm RF}N^2)$ with $N\in \{N_{\rm T},N_{\rm R}\}$ which is much smaller than the complexity of other steps in the PS-HBF algorithm, i.e., $\Oc\left(2KN_{\rm T}N_{\rm R}\left(N_{\rm T} +N_{\rm R}\right) \right)$. Furthermore, because of the sparse nature of the channel in mmWave and THz mMIMO systems, the multiplexing gain is limited by the number of effective paths, i.e., $L_{\rm p}$ in (<ref>), which is usually small. To achieve the maximum multiplexing gain, the required number of RF chain can be set to $N_{\rm RF}=L_{\rm p}$. As such, the complexity gap between the proposed SW-HBF algorithms and PS-HBF algorithms is limited. Finally, it can be observed from Fig. <ref> that the complexities of the proposed SW-HBF algorithms and the PS-HBF algorithm slowly increase with the number of subcarriers.

§ CONCLUSIONS

This paper investigates the SW-HBF technique in wideband mMIMO systems considering the beam squint effect. We derive a closed-form expression that quantifies the beam squint effect for PS-HBF schemes. The results show that the beam squint effect is scaled with the number of antennas, fractional bandwidth, and the normalized antenna separation. Further, we show that the SW-HBF architectures are less affected by the beam squint effect than the PS-HBF schemes by analyzing the EAG of the PS-based array and switch-based array. Moreover, we have developed several efficient SW-HBF algorithms by maximizing the average SE. The analog beamformers are efficiently found by proposed search algorithms, while the digital beamformers are solved via closed-form expression. The numerical results show that the proposed SW-HBF schemes can attain near-optimal performance and outperform the PS-HBF schemes in terms of both SE and EE when the beam squint effect is severe. The results reveal the superiority of the SW-HBF schemes over the PS-HBF schemes, especially in wideband mMIMO systems. Future works in this area may include channel estimation, robust SW-HBF methods, spatial oversampling, and low-resolution ADCs/DACs.

§ PROOF OF PROPOSITION <REF>

As it is difficult to calculate the integration in (<ref>), we utilize the Newton polynomial to fit it by three points $(-1,g \left( -1,\xi_k\right) ),(0,g\left( 0,\xi_k\right)),(1,g \left( 1,\xi_k\right))$, which leads

\begin{align}\label{eq:Newton polynomial interpolation}
&\int^{1}_{-1}  g\left( \Bar{\theta},\xi_k \right) d\Bar{\theta} \nonumber \\
&\approx \int^{1}_{-1} \left[\left(g \left( -1,\xi_k\right) -g\left( 0,\xi_k\right)\right)\Bar{\theta}^2+g\left( 0,\xi_k\right)\right] d\Bar{\theta}  \nonumber\\
&= \frac{4}{3}g\left( -1,\xi_k\right)  + \frac{2}{3}g\left( 0,\xi_k\right).
\end{align}

Since the normalized sinc function achieves the maximum value 1 when the variable is 0, i.e., $g\left( 0,\xi_k\right)=1$, the EAG of PS-based array is expressed as

\begin{align}\label{eq:expection array gain PS interpolated}
\Es_{\rm ps}[g\left( \Bar{\theta},\xi_k \right)]&=\frac{2}{3K} \sum\limits_{k=1}^K g\left( -1,\xi_k\right) +\frac{1}{3}.
\end{align}

Substituting (<ref>) into (<ref>), the EAG of PS-based array can be written as

\begin{align}
\Es_{\rm ps}[g\left( \Bar{\theta},\xi_k \right)]&=\frac{2}{3K}\sum\limits_{k=1}^K\left|  \frac{\sin\left( N \Delta b \pi c_k   \right)}{N \sin\left(  \Delta b \pi c_k   \right)}\right|+ \frac{1}{3} \nonumber \\
&\overset{(i)}{\approx} \frac{2}{3K}\sum\limits_{k=1}^K\left|  \frac{\sin\left( 2\rho \pi c_k  \right)}{2\rho \pi c_k }\right|+ \frac{1}{3}, \label{eq:EAG med expression}\\
&\overset{(ii)}{\approx}\frac{2}{3} \int_0^1 \left| \frac{\sin\left( 2\rho \pi (x-\frac{1}{2})   \right)}{2 \rho \pi (x-\frac{1}{2})}\right|dx+ \frac{1}{3} \nonumber \\
&=\frac{2}{3\rho} \int_0^{\rho} \left|{\rm sinc}(\pi x)\right|dx + \frac{1}{3} ,
\end{align}

where $\rho\triangleq\frac{N\Delta b}{2}$, $c_k\triangleq\left|\frac{k}{K}-\frac{K+1}{2K}\right|$; $(i)$ and $(ii)$ are due to $\Delta b \pi c_k \ll 1$ and $K \gg 1$, respectively. Defining $f(\rho)\triangleq\frac{1}{\rho} \int_0^{\rho} \left|{\rm sinc}(\pi x)\right|dx$, we claim that $f(\rho)$ is monotonically decreasing with $\rho$. The first derivative of $f(\rho)$ is

\begin{align}
f'(\rho)&=\frac{1}{\rho}\left( \left|{\rm sinc}(\pi \rho) \right|- \frac{1}{\rho} \int_0^{\rho} \left|{\rm sinc}(\pi x)\right|dx\right) \nonumber\\
&\overset{(i)}{=}\frac{1}{\rho}\left( \left|{\rm sinc}(\pi \rho) \right|- \left|{\rm sinc}(\varepsilon  \pi) \right| \right),
\end{align}

where $(i)$ is due to the utilization of the Mean Value Theorem with $\varepsilon \in (0,\rho)$.

* For $\rho \in (0,\pi)$, the ${\rm sinc}(\pi x)$ is non-negative and monotonically decreasing, thus $f'(\rho)<0$.

* For $\rho \in [\pi,\infty)$, applying $\left|{\rm sinc}(\pi \rho) \right|\leq \frac{1}{\pi \rho}$, we obtain

\begin{equation}
f'(\rho) \leq \frac{1}{\rho}\left( \frac{1}{\pi \rho}- \frac{1}{\pi \varepsilon  } \right) < 0
\end{equation}

In summary, the fist derivative of $f(\rho)$ is negative in $(0,\infty)$, which indicates $f(\rho)$ is monotonically decreasing with $\rho$ in $(0,\infty)$. From (<ref>), we can obtain

\begin{equation}
\frac{1}{3}\leq \Es_{\rm ps}[g\left( \Bar{\theta},\xi_k \right)]\leq 1.
\end{equation}

The lower bound is attained at $\rho \rightarrow \infty$, while the upper bound is achieved at $\rho \rightarrow 0$.