The explosive growth of mobile data, as well as the emergence of innovative technologies such as virtual reality and holographic communication, require future wireless communication systems with extremely high data rates, seamless coverage, and ultra-high reliability. Fortunately, Reconfigurable Intelligent Surface (RIS) can enable future wireless communication systems to meet these requirements [1–3]. Therefore, RIS has been widely studied.

The authors in [4] focus on an RIS-aided multicell multiple-input single-output (MISO) system and maximize the minimum weighted signal-to-interference-plus-noise ratio (SINR) through jointly optimizing the transmit beamforming vectors at the base stations (BS) and the reflective beamforming vector of the RIS. The authors in [5] investigate a RIS backscatter-based uplink coordinated transmission strategy and aim to maximize the weighted sum rate (WSR). In [6], the authors suggest RIS be deployed at the boundaries of multiple cells to mitigate inter-cell interference and assist cell-edge users in communicating with BS, aiming to maximize the WSR by jointly optimizing the phase shifts of RIS and the precoding matrix at the BS. In [7], the authors consider the phase shifts design of RIS and joint user scheduling in RIS-aided multicell downlink systems. The authors in [8] study an intelligent omni-surface-aided multicell MIMO communication system and design a distributed hybrid beamforming scheme to maximize the sum rate. In [9], besides considering the joint beamforming problems for a RIS-aided multicell MISO system, authors further consider the SINR balancing beamforming design to enhance the fairness of users by maximizing the minimum SINR among all users.

However, the transceiver hardware impairments (THWIs) in the scenario of RIS-aided multicell massive MIMO system have not been considered. Therefore, we investigate an RIS-aided multicell mMIMO system with THWIs. More importantly, the phase optimization of RIS is based on the more practical and feasible two-timescale design scheme. Specifically, the phase shifts of the RIS are designed based on the statistical channel state information (CSI). In contrast, while the beamforming at the BS is designed based on the instantaneous CSI, which can significantly reduce the frequency of reconfiguring RIS [10]. We first derive and analyze the expression of uplink achievable data rate. Then, we use a genetic algorithm (GA) to optimize the phase shifts of RIS to maximize the sum and minimum user rates. Finally, we verify the accuracy of derivation through Monte Carlo (MC) simulation. The simulation shows that in the RIS-aided multicell mMIMO system with THWIs, appropriately increasing the number of RIS’s reflecting elements can compensate for the performance loss caused by inter-cell interference and THWIs, and also reduce the demand for the number of BS antennas, which can significantly reduce the hardware costs at the BS.

As shown in Figure 1, we consider an RIS-aided uplink multicell mMIMO communication system. The number of cells is L, and the RIS comprises N reflecting elements. The base station (BS) in each cell is equipped with M receiving antennas, both the BS and RIS adopt the uniform square planar array. In each cell, the direct link between users and BS is blocked by obstacles, K single-antenna users can only communicate with BS through RIS, and both the transmitters and receivers have hardware impairments.

The phase shifts matrix of RIS is $\mathbf{\Theta }=\mathrm{diag}\left\{e^{j{\theta }_1},...,e^{j{\theta }_n},...,e^{j{\theta }_N}\right\}$, where ${\theta }_n\in [0,2\pi )$, $n=1,2,...,N$, ${\theta }_n$ represents the phase shift of element n. In the l-th cell, the channel from the users to the RIS is denoted by ${\mathbf{H}}_{l,ur}=\left[{\mathbf{h}}_{l,1},...,{\mathbf{h}}_{l,k},...,{\mathbf{h}}_{l,K}\right]\in {\mathbb{C}}^{N\times K}$, $l=1,...,L$, where ${\mathbf{h}}_{l,k}$ represents the channel from the k-th user to the RIS. The channel from the RIS to the BS is denoted by ${\mathbf{H}}_{l,rb}\in {\mathbb{C}}^{M\times N}$. The expressions of ${\mathbf{H}}_{l,rb}$ and ${\mathbf{h}}_{l,k}$ are where ${\nu }_l$ and ${\mu }_{l,k}$ are large scale path loss factors, ${\rho }_l$ and ${ϵ}_{l,k}$ are Rician factors, ${\tilde{\mathbf{H}}}_{l,rb}$ and ${\tilde{\mathbf{h}}}_{l,k}$ are the none-line-of-sight (NLoS) components of the channel whose elements are independent and identical distribution (i.i.d.) random variables following $\mathcal{C}\mathcal{N}\left(0,1\right)$. ${\overline{\mathbf{H}}}_{l,rb}$ and ${\overline{\mathbf{h}}}_{l,k}$ are LoS components of channel, and the expressions are where ${\psi }_{l,rb}^a$, ${\psi }_{l,rb}^e$ respectively represent the arrival azimuth and elevation angles from RIS to the BS in the l-th cell. ${\varphi }_{l,rb}^a$, ${\varphi }_{l,rb}^e$ respectively represent the departure azimuth and elevation angles from the RIS to the BS in the l-th cell. ${\psi }_{l,kr}^a$, ${\psi }_{l,kr}^e$ respectively represent the arrival azimuth and elevation angles from user k to the RIS in the l-th cell. In addition, the expressions of ${\mathbf{a}}_M\left(\begin{array}{c}{\psi }_{l,rb}^a,{\psi }_{l,rb}^e\end{array}\right)$ is where $0\le p,q\le \sqrt{M}-1$, d represents the elements spacing, $\lambda$ represents wavelength.

Based on the above definitions, the cascaded channel in the l-th cell can be denoted by ${\mathbf{G}}_l={\mathbf{H}}_{l,rb}\mathbf{\Theta }{\mathbf{H}}_{l,ur}=\left[{\mathbf{g}}_{l,1},...{\mathbf{g}}_{l,k},...{\mathbf{g}}_{l,K}\right]$, where ${\mathbf{g}}_{l,k}={\mathbf{H}}_{l,rb}\mathbf{\Theta }{\mathbf{h}}_{l,k}$.

In the l-th cell, the transmit distortion of users is denoted by ${\mathbf{z}}_{l,t}={\left[z_{l,1},...z_{l,k},...,z_{l,K}\right]}^T$, where $z_{l,k}$ represents the transmit distortion of the user k. The receive distortion of BS is denoted by ${\mathbf{z}}_{l,r}$.

Therefore, the signal received by the BS of the l-th cell is where ${\mathbf{P}}_l=\mathrm{diag}\left(\sqrt{p_{l,1}},...,\sqrt{p_{l,K}}\right)$ is the transmit power of users. ${\mathbf{x}}_l=\left[x_{l,1},...,x_{l,K}{]}^T\right.$ is the signal vector of users, and $\mathbb{E}\left\{|x_{l,i}{|}^2\right\}=1$. ${\mathbf{n}}_l\sim \mathcal{C}\mathcal{N}\left(0,{\sigma }^2{\mathbf{I}}_N\right)$ represents the additive white Gaussian noise (AWGN). The transceiver distortions are described in terms of conditional distributions with respect to the channel realizations as $z_{l,i}\sim \mathcal{C}\mathcal{N}\left(0,k_{l,u}{p}_{l,i}\right)$ and ${\mathbf{z}}_{l,r}\sim \mathcal{C}\mathcal{N}\left(0,k_{l,b}{\displaystyle \sum _{j=1}^L}{\displaystyle \sum _{i=1}^K}{p}_{j,i}{\mathbf{I}}_M{\mathbf{g}}_{j,i}{\mathbf{g}}_{j,i}^H\right)$, where $k_{l,u}$ and $k_{l,b}$ represent the severity of the residual impairments of the transceiver.

Then, we employ the MRC technology, and the received signal at the BS can be expressed as

Therefore, in the l-th cell, the k-th user’s signal received by the BS can be expressed as

The uplink data rate of user k can be denoted by $R_{l,k}=\mathbb{E}\left\{{log}_2(1+{\mathrm{SINR}}_{l,k})\right\}$, where ${\mathrm{SINR}}_{l,k}$ is the signal-to-interference-plus-noise ratio (SINR) of the user k in l-th cell. The expression of ${\mathrm{SINR}}_{l,k}$ is given as

The expressions of ${\mathbb{E}}_{\mathrm{signal}}^{l,k}$, ${\mathbb{E}}_{\mathrm{intra}}^{l,k}$, ${\mathbb{E}}_{\mathrm{inter}}^{j,i}$, ${\mathbb{E}}_{\mathrm{noise}}^{l,k}$ and ${\mathbb{E}}_{\mathrm{thwis}}^{l,k}$ are respectively given by (11)-(15), where $a_{l,k}=\frac{{\nu }_l{\mu }_{l,k}}{\left({\rho }_l+1\right)\left({ϵ}_{l,k}+1\right)}$, ${\mathbf{\Phi }}_{l,k}\left(\mathbf{\Theta }\right)={\mathbf{a}}_N^H\left({\varphi }_{l,rb}^a,{\varphi }_{l,rb}^e\right)\mathbf{\Theta }{\stackrel{-}{\mathbf{h}}}_{l,k}$, ${\mathbf{g}}_{l,k_m}$ is the m-th element of ${\mathbf{g}}_{l,k}$,

By observing (10), We find that the phase shift term ${\mathbf{\Phi }}_{l,k}\left(\mathbf{\Theta }\right)$ is always related to the Rician factors ${\rho }_l$, ${ϵ}_{l,k}$, ${ϵ}_{l,i}$, $\forall l,k,i$. We can conclude that if the channel in this system is fully Rayleigh channel, namely ${\rho }_l={ϵ}_{l,k}={ϵ}_{l,i}=0$, $\forall l,k,i$, the performance of this system will be independent of the phase shifts and only depends on the number of RIS elements N and BS antennas M. Hence, it is important to reveal the relation between the system performance and the number of RIS elements N and BS antennas M.

It can be seen that the noise term ${\mathbb{E}}_{\mathrm{noise}}^{l,k}$ is not included in equation (19) and the rate will tend to be a constant. This shows that in the RIS-aided uplink multicell mMIMO communication system with THWIs when $N\to \infty$ and the phase shifts of RIS are random, we can eliminate the effect of noise.

As same as Corollary 1, when the number of BS antennas M tends to be infinite, the rate of this system will tend to a constant too. Besides, the effect of noise can be ignored.

By observing (10), we found that the rate only depends on statistical CSI including the Rician factors, LoS components of the channel and path loss coefficients. Therefore, we can optimize the phase shifts of RIS based on statistical CSI.

Specifically, we respectively formulate two optimization problems the sum rate maximization problem to guarantee overall system performance and the minimum user rate maximization problem to guarantee fairness between users as (27), (28). where $R_{l,k}$ is given in (10).

Due to the complexity of the rate expression, the existing gradient descent method and convex optimization method cannot solve this problem. The genetic algorithm (GA) does not need to calculate the first-order derivative of the objective function as the conventional gradient descent method does. Instead, we only need to calculate the objective function in each iteration, so we use GA. The detailed steps of the algorithm are as follows

Initialize population: Generate S individuals, and each individual has N chromosomes randomly generated from $[0,2\pi )$, corresponding to the phase shifts matrix of RIS.

Calculate fitness: In the current population, we first calculate the fitness of each individual through the objective function of the optimization problem (32) or (33). Then, we sort them in descending order according to fitness.

Select elite: We select the top $S_1$ individuals as elites and pass them to the next generation.

Mutation: We select the last $S_3$ individuals as parents and mutate them with a probability of 0.1 to generate $S_3$ offspring. The mutation algorithm is shown as follows.

Crossover: We use stochastic universal sampling to generate $2S_2$ parents from the remaining $S_2=S-S_1-S_3$ individuals. Then, we use the two-points crossover method to generate $S_2$ offspring from $2S_2$ parents. The two-points crossover method is shown as follows.

Finally, we combine $S_1$ elites and $S_2+S_3$ offspring to generate the next generation population. When the change value of average fitness is less than the preset value or the number of iterations reaches the preset maximum number of iterations, the algorithm stops. In the current population, the chromosome of the individual with the highest fitness corresponds to the optimal phase shifts of RIS.

Algorithm 1: Mutation Algorithm

1. for $s=1:S_3$ do;

2. for $n=1:N$ do;

3. $r=rand\left(1\right)$;

4. if $r<0.1$ then ;

5. The n-th chromosome of individual mutates to $2\pi *r$;

6. end if

7. end for

8. end for

In this section, we provide simulation results to validate the derived expressions and analyze various system parameters. The parameter settings were referred to [6] and [12]. Unless otherwise stated, we set it as follows: $N=49,L=2,M=49,k_{l,u}=k_{l,b}=0.01,K=4,p_{l,k}=30$dBm, ${\sigma }^2=-104$dBm, ${\rho }_l=1,{ϵ}_{l,k}=10,\forall l,k$. We assume that the distance from RIS to BS is $d_{rb}=800m,\forall l$, users are distributed on a circle with RIS as the center and a radius of $d_{ur}=30m$. Therefore, ${\nu }_l={10}^{-3}{d}_{rb}^{-2.5}$, ${\mu }_{l,k}={10}^{-3}{d}_{ur}^{-2},\forall l,k$. In addition, the element spacing of $d=\lambda /2$ is set to avoid spatial correlation. The various angles at BS and RIS are randomly generated from $[0,2\pi ]$ and fixed after initial generation.

As shown in Figure 2, we plot the achievable rate obtained from optimization problems (32) and (33) and carry out MC simulation for verification. The MC simulation results are consistent with the derived results, which shows the accuracy of the derived results.

In Figure 3, we verified that the optimization problem (33) can ensure users’ fairness. The relationship between the two vertical coordinates is four times, and the two curves coincide with each other, while the number of users $K=4$, means that the minimum user rate reaches the average user rate, i.e., all users have the same rate.

Figure 4 depicts the achievable rate in different scenarios. As shown in Figure 4, scheme ’$N=49$, Sum rate by max-sum’ is superior to scheme ’$N=100$, Sum rate by random phase’, and as N increases, the performance gap between the two will become larger. This is consistent with our Corollary 2 that the rate will tend to be a constant when BS antennas $M\to \infty$.

Besides, it indicates that it is necessary to optimize the phase shifts of RIS.

In Figure 5, we plot the achievable rate versus the Rician factor of the RIS-BS channel. As ${\rho }_l$ increases, the sum and minimum user rates decrease. This is because, with the increase of ${\rho }_l$, the channels between users have a stronger correlation, which reduces the spatial multiplexing gain and also increases the interference between users.

In Figure 6, we plot the achievable rate with different THWIs coefficients. As the transmit power increases, the sum rate and minimum user rate gradually increase and tend towards a limit value. In addition, with the increase of power, the performance gap between the two THWIs scenarios also widens because the THWIs are tightly coupled with the transmitted signal and become stronger as the THWIs coefficients $k_{l,u},k_{l,b}$ increase.

Figure 7 depicts the achievable rate under different schemes. Through comparison, we find that the rate performance of scheme ’$M=256$, $N=16$’ is the same as that of scheme ’$M=N=64$’. This means that we can significantly reduce the demand for the number of BS antennas by appropriately increasing the number of RIS’s reflecting elements, effectively reducing the hardware cost at BS.

To ensure the accuracy of scheme comparison, we further draw Figure 8. As shown in Figure 8, the scheme ’$M=N=64$’ is superior to the scheme ’$M=256$, $N=16$’ in different THWIs and power scenarios, which shows the accuracy of the scheme comparison.

In Figure 9, we plot the sum rate for three scenarios. In scenario $k_{l,u}=k_{l,b}=0.1,L=4$, Sum rate by max-sum’, the impact of inter-cell interference and THWIs can be compensated by increasing the number of components. This is meaningful because increasing N is more cost-effective and practical compared to replacing transceivers and reducing the number of cells.

This paper has focused on an RIS-aided uplink multicell mMIMO communication system with THWIs. We have derived and analyzed the achievable rate expression, and the phase shifts of RIS have been optimized based on the two-time scale design scheme. Through simulation, the correctness of the derived expression has been verified, and the following conclusion has been obtained. In RIS-aided uplink multicell mMIMO communication system with THWIs, we can compensate for the performance loss caused by inter-cell interference and THWIs by appropriately increasing the number of RIS’s reflecting elements, and also significantly reduce the demand for the number of BS receiving antennas, which will effectively reduce the hardware cost at the BS.