# Beamforming Design and Association Scheme for Multi-RIS Multi-User mmWave Systems Through Graph Neural Networks

Mengbing Liu, Chongwen Huang, Ahmed Alhammadi, Marco Di Renzo,
Mérouane Debbah, Chau Yuen

<figure><img src="image_1.png" /></figure>

## To cite this version:

Mengbing Liu, Chongwen Huang, Ahmed Alhammadi, Marco Di Renzo, Mérouane Debbah, et al. Beamforming Design and Association Scheme for Multi-RIS Multi-User mmWave Systems Through Graph Neural Networks. IEEE Transactions on Wireless Communications, 2025, 24 (9), pp.7940-7954. 10.1109/TWC.2025.3563529. hal-05378146

HAL Id: hal-05378146

https://hal.science/hal-05378146v1

Submitted on 22 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Beamforming Design and Association Scheme for Multi-RIS Multi-User mmWave Systems Through Graph Neural Networks

Mengbing Liu, *Graduate Student Member, IEEE*, Chongwen Huang, *Member, IEEE*, Ahmed Alhammadi, *Member, IEEE*, Marco Di Renzo, *Fellow, IEEE*, Mérouane Debbah, *Fellow, IEEE*, and Chau Yuen, *Fellow, IEEE*

**Abstract**—Reconfigurable intelligent surface (RIS) is emerging as a promising technology for next-generation wireless communication networks, offering a variety of merits such as the ability to tailor the communication environment. Moreover, deploying multiple RISs helps mitigate severe signal blocking between the base station (BS) and users, providing a practical and efficient solution to enhance the service coverage. However, fully reaping the potential of a multi-RIS aided communication system requires solving a non-convex optimization problem. This challenge motivates the adoption of learning-based methods for determining the optimal policy. In this paper, we introduce a novel heterogeneous graph neural network (GNN) to effectively leverage the graph topology of a wireless communication environment. Specifically, we design an association scheme that selects a suitable RIS for each user. Then, we maximize the weighted sum rate (WSR) of all the users by iteratively optimizing the RIS association scheme, and beamforming designs until the considered heterogeneous GNN converges. Based on the proposed approach, each user is associated with the best RIS, which is shown to significantly improve the system capacity in multi-RIS multi-user millimeter wave (mmWave) communications. Specifically, simulation results demonstrate that the proposed heterogeneous GNN closely approaches the performance of the high-complexity alternating optimization (AO) algorithm in the considered multi-RIS aided communication system, and it outperforms other benchmark schemes. Moreover, the performance improvement achieved through the RIS association scheme is shown to be of the order of 30%.

**Index Terms**—Beamforming, heterogeneous graph neural network, reconfigurable intelligent surfaces.

Mengbing Liu and Chau Yuen are with the School of Electrical and Electronics Engineering, Nanyang Technological University, Singapore 639798 (E-mails: mengbing001@e.ntu.edu.sg; chau.yuen@ntu.edu.sg).

Chongwen Huang is with College of Information Science and Electronic Engineering, Zhejiang University, Hangzhou 310027, China, and Zhejiang Provincial Key Laboratory of Multi-Modal Communication Networks and Intelligent Information Processing, Hangzhou 310027, China (E-mail: chongwenhuang@zju.edu.cn).

Ahmed Alhammadi is with the Technology Innovation Institute, Masdar City, Abu Dhabi, United Arab Emirates (E-mail: Ahmed.Alhammadi@tii.ae).

Marco Di Renzo is with Université Paris-Saclay, CNRS, CentraleSupélec, Laboratoire des Signaux et Systèmes, 3 Rue Joliot-Curie, 91192 Gif-sur-Yvette, France. (E-mail: marco.di-renzo@universite-paris-saclay.fr), and with King's College London, Centre for Telecommunications Research - Department of Engineering, WC2R 2LS London, United Kingdom (E-mail: marco.di-renzo@kcl.ac.uk).

Mérouane Debbah is with KU 6G Research Center, Department of Computer and Information Engineering, Khalifa University, Abu Dhabi 127788, UAE and also with CentraleSupelec, University Paris-Saclay, 91192 Gif-sur-Yvette, France (E-mail: merouane.debbah@ku.ac.ae).

Parts of this paper have been presented in part at the IEEE International Conference on Communication (ICC), 2023 [1].

## I. INTRODUCTION

Millimeter wave (mmWave) communications play an indispensable role in providing high-bandwidth and high-throughput services for multimedia applications in current and next-generation wireless systems [2]. Although mmWave communications allow for higher data rates due to the wider bandwidth, signals transmitted in high-frequency bands experience serious propagation attenuation and penetration losses [3]. As a remedy, reconfigurable intelligent surface (RIS) has emerged as a promising technology to enhance the link quality of service and restore communication in dead zones. Compared to alternatives like ultra-dense networks and massive active antenna arrays, RIS is an effective and cost-efficient paradigm for assisting mmWave networks through a large number of nearly passive reflective elements [4]–[7].

The integration of RIS into various communication systems may improve the network capacity and reduce the energy consumption. However, jointly designing the beamforming of base stations (BSs) and RISs remains a critical challenge, as the corresponding optimization problem for distributed RIS-aided systems is non-convex and computationally intensive. Consequently, significant efforts have been directed toward developing efficient algorithms to address this challenge in RIS-aided systems [8]–[14]. Recently, the beamforming design in a multi-RIS aided system has been investigated to compensate for the severe propagation attenuation in high-frequency bands [15]–[19]. These problem-specific algorithms, however, are laborious and require extensive knowledge, particularly in multi-RIS aided systems.

Inspired by the success of deep learning in various application domains, e.g., computer vision, deep neural network (DNN) models have been employed in wireless communication systems in recent years [20]–[25]. For instance, an unsupervised deep learning approach was proposed in [22] to obtain the optimal phase configuration in an RIS-aided single-user wireless communication system. In [23], a deep quantization neural network was developed for optimizing the beamforming matrices at the BS and RIS, by assuming discrete phase shifts and imperfect channel state information (CSI).

However, such DNN models may not be well-suited for applications in wireless communications, as they fail to capture the underlying network topology, often leading to dramatic performance degradation when the network size increases. To
enhance the scalability and generalization of deep learning for wireless applications, it is crucial to incorporate the network structure into the neural network architecture [26].

From this perspective, graph neural networks (GNNs) have emerged as a suitable model to exploit the inherent graph topology of wireless communication networks [27]–[29]. A GNN is a type of neural network designed to operate on graphs and other structured data, allowing for the utilization of abundant information with structural relations among the network nodes [30]. For example, recent research works apply GNNs for resource allocation in device-to-device (D2D) communication networks [31], [32]. In these works, the service links and interference relations are modeled as nodes and edges, respectively. Based on these graph embedding features, the beamforming design for maximizing the weighted sum rate (WSR) can be solved efficiently. Recently, GNN architectures have been applied to the design of RIS-aided systems [33]–[36]. In [33], a GNN was employed to learn both BS and RIS beamforming directly from received pilots, thereby maximizing the WSR without explicit channel estimation. Building on this GNN architecture, [34] proposed a three-stage algorithm with two GNNs to jointly determine beamforming matrices and user scheduling when the number of users greatly exceeds the number of BS antennas. Recently, in [35], a GNN-based framework jointly optimizes user association, BS beamforming, and RIS phase adjustment without explicit CSI, boosting sum rate and enhancing load balancing. Meanwhile, [36] proposes an enhanced-GNN (E-GNN) featuring improved angular feature extraction, graph attention, and transfer learning to overcome low signal-to-noise ratio (SNR) instability and multi-path challenges, achieving near-optimal sum-rate performance.

In this paper, we leverage the advantages of GNNs for beamforming optimization in a multi-RIS multi-user mmWave communication system. Besides beamforming optimization, we consider a more general problem in which each user selects the most suitable RIS to improve the link quality of service and energy efficiency. In single RIS systems, the BS must spread its power via one RIS to serve all the users, which is not efficient. By assigning each user to a specific RIS, instead, the system can adjust the RIS phase configurations to concentrate the signal on each user, reducing power waste. Furthermore, the RIS association scheme improves system flexibility and scalability while reducing the computational overhead.

For this reason, an RIS association design is introduced into the beamforming design problem to select the most suitable RIS for each user so as to maximize the WSR, i.e., aiming to provide the best performance while reducing the number of RISs. Though there exist some research works analyzing the BS-RIS-user association problem [37]–[39], they focus on the association between multiple BSs and users in RIS-aided and full-connected multi-RIS aided communication scenarios. Different from previous works, we primarily consider the association scheme between multiple RISs and multiple users to improve the system’s performance while improving the flexibility and efficiency of the system by selecting the most appropriate RIS.

Meanwhile, existing GNN structures with single-type nodes

have only been applied to the beamforming design in RIS-aided communication systems, and do not incorporate the association scheme between multiple RISs and multiple users [31], [32]. Considering the strong coupling between the beamforming and association designs for multiple RISs, we propose a new graph modeling approach to better map the beamforming design and RIS association scheme in the proposed system, i.e., a heterogeneous GNN, which contains multi-type nodes with specific node features tailored to the system requirements [40]. Overall, the main contributions of this work are the following:

* We first present a heterogeneous GNN structure and map it onto an optimization problem tailored to multi-RIS aided systems. The BS, RISs, and users are represented as three different types of nodes and we embed different optimization variants as features of the corresponding types of nodes, i.e., the beamforming matrix in the BS node, the phase configuration matrix in each RIS node, and the RIS association matrix in the user node, respectively.

* Two different loss functions are adopted to address the challenge of learning continuous (the beamforming matrix and phase configuration matrices) and discrete variables (the RIS association matrix) simultaneously. One is the WSR of the multi-RIS aided communication system, and maximizing it is the target of the considered optimization problem. The second one is the categorical cross-entry (CE) loss, which helps the learning process of the RIS association scheme. In order to balance the significance between two loss functions during training, a penalty term is utilized for the loss function.

* To thoroughly verify the effectiveness of the proposed structure, a series of commonly used benchmarks are illustrated. The optimized results are very close to those obtained by using the traditional alternating optimization (AO) algorithm. In addition, it is shown that the proposed heterogeneous GNN scheme performs 10 times better than a DNN-based scheme.

* We conduct a comparative experiment, where the proposed RIS association scheme is determined in advance according to the distance. Then, the heterogeneous GNN network is used to learn the optimal beamforming and phase configuration matrices. Extensive experiments demonstrate that the WSR obtained by jointly optimizing the three variables is higher than that of using the predetermined scheme. Finally, compared to the case study when the users are associated with the single RIS, the performance improvement achieved by the proposed association scheme is around 30%.

The remainder of this paper is organized as follows: Section II describes the multi-RIS aided system model and the channel model for mmWave bands. Also, the problem formulation for maximizing the WSR is introduced. Section III presents the message-passing paradigm and normalization layers for the proposed heterogeneous GNN scheme. Section IV illustrates the complete algorithmic frameworks with different predetermined RIS association strategies. Section V reports the
simulation results of the proposed GNNs and the comparison with relevant benchmarks verifying the effectiveness of the proposed approach. Finally, Section VI concludes the article.

**Notation:** Uppercase (lowercase) boldface letters denote matrix (vector) and the field of complex (real) numbers is denoted by $\mathbb{C}(\mathbb{R})$. $(\cdot)^T$ and $(\cdot)^H$ indicate the operations of transpose, and Hermitian transpose. $[\cdot]_{m,n}$ denotes the elements at the $m$-th row and the $n$-th column of a matrix, and diag($\mathbf{x}$) is a diagonal matrix whose diagonal elements are formed with the elements of $\mathbf{x}$. Besides, $\mathbf{A} \otimes \mathbf{B}$ denotes the Kronecker product of $\mathbf{A}$ and $\mathbf{B}$. For a complex number, $\mathrm{Re}\{\cdot\}$ and $\mathrm{Im}\{\cdot\}$ denote the real part and imaginary part. $|\cdot|$ and $\|\cdot\|_{\mathcal{F}}$ represent the absolute value and Frobenius norm, respectively. $v \sim \mathcal{CN}(0, \sigma^2)$ means that $v$ follows the complex Gaussian distribution with zero-mean and variance $\sigma^2$.

## II. MULTI-RIS MULTI-USER MMWAVE COMMUNICATION NETWORKS

In this section, we introduce a multi-RIS aided mmWave communication system model and the RIS association scheme. Also, the channel model for mmWave systems and the WSR maximization problem with respect to the beamforming design and the RIS association scheme are presented.

### A. System Model

We consider a multi-RIS aided system with an $N_t$-antenna BS, a set of single-antenna users $\mathcal{K} = \{1, 2, \dots, K\}$ and a set of RISs $\mathcal{R} = \{1, 2, \dots, R\}$. Each RIS is equipped with $M = M_x \times M_y$ reflecting elements, where $M_x$ and $M_y$ represent the number of horizontal and vertical elements, respectively. Assume that the BS transmits $K$ data streams to the users through beamforming optimization. Thus, the transmitted signal $\mathbf{x}$ from the BS can be expressed as

$$ \mathbf{x} \triangleq \mathbf{F}\mathbf{s} = \sum_{k=1}^{K} \mathbf{f}_k s_k, \quad (1) $$

where $\mathbf{s} = [s_1, \dots, s_K]^T \in \mathbb{C}^{K \times 1}$ indicates the symbol vector and $\mathbb{E}(\mathbf{s}\mathbf{s}^H) = \mathbf{I}_K$. $\mathbf{F} = [\mathbf{f}_1, \dots, \mathbf{f}_K] \in \mathbb{C}^{N_t \times K}$ denotes the beamforming matrix at the BS and the power constraint is denoted by $\|\mathbf{F}\|_{\mathcal{F}}^2 \le P_{\max}$, where $P_{\max}$ is the maximum transmit power at the BS.

The direct BS-user links are assumed to be blocked by obstacles. The BS can only communicate with the users via the reflections from the RISs. Therefore, it is crucial to investigate the joint design of RIS association, BS beamforming, and the phase configurations of the RISs. In addition, we assume that the CSI of all the links is perfectly known to the BS<fnref n="1" />, and all the channels experience quasi-static flat fading.

In this paper, each user is associated with exactly one RIS, which is selected from the available RISs to optimize the communication performance by maximizing the WSR of the multi-RIS aided communication system. To easily express the association between the RISs and users, a matrix

<figure><img src="image_3.png" /><figcaption>Fig. 1. Illustration of the multi-RIS aided system model.</figcaption></figure>

$\mathbf{U} = [\mathbf{u}_1, \dots, \mathbf{u}_K]^T \in \mathbb{R}^{K \times R}$ consisting of binary variables $u_{k,i}$ is defined as

$$ u_{k,i} = \begin{cases} 1, & \text{user}_k \text{ is assisted by RIS}_i \\ 0, & \text{otherwise.} \end{cases} \quad (2) $$

Based on the RIS association scheme, the equivalent channel $\mathbf{h}_k$ from the BS to user$_k$, including the cascaded links, can be expressed as

$$ \mathbf{h}_k = \sum_{i \in \mathcal{R}} u_{k,i} \mathbf{h}_{ik} \mathbf{\Phi}_i \mathbf{G}_i, \quad (3) $$

where $\mathbf{G}_i \in \mathbb{C}^{M \times N_t}$ and $\mathbf{h}_{ik} \in \mathbb{C}^{1 \times M}$ denote the BS-RIS$_i$ and RIS$_i$-user$_k$ channel, respectively. $\mathbf{\Phi}_i = \mathrm{diag}(\theta_{i,1}, \dots, \theta_{i,M}) \in \mathbb{C}^{M \times M}$ with $\theta_{i,m} = e^{j\phi_{i,m}}$ denotes the phase shift matrix of RIS$_i$.

In the considered model, the RIS association scheme optimizes the signal reflections by matching users with specific RISs, focusing on enhancing the signal quality for the associated users. Considering the distributed arrangement of RISs and the consequent signal attenuation due to propagation distances, the signals reflected toward non-associated users are markedly attenuated compared to those directed at associated users. Therefore, some interfering signals, for ease of analysis, are not considered due to their weak impact. Thus, the signal $\mathbf{y}_k$ received by the user$_k$ can be rewritten as

$$ \begin{aligned} \mathbf{y}_k &= \mathbf{h}_k \mathbf{x} + \mathbf{n}_k \\ &= \left( \sum_{i \in \mathcal{R}} u_{k,i} \mathbf{h}_{ik} \mathbf{\Phi}_i \mathbf{G}_i \right) \mathbf{f}_k s_k \\ &\quad + \sum_{j \in \mathcal{K}, j \neq k}^{K} \left( \sum_{i \in \mathcal{R}} u_{k,i} \mathbf{h}_{ik} \mathbf{\Phi}_i \mathbf{G}_i \right) \mathbf{f}_j s_j + \mathbf{n}_k. \end{aligned} \quad (4) $$

According to (4), the signal-to-interference-plus-noise-ratio (SINR) at user$_k$ is given as

$$ \mathrm{SINR}_k = \frac{\left| \left( \sum_{r \in \mathcal{R}} u_{r,k} \mathbf{h}_{rk} \mathbf{\Phi}_r \mathbf{G}_r \right) \mathbf{f}_k \right|^2}{\sum_{j \in \mathcal{K}, j \neq k}^{K} \left| \left( \sum_{r \in \mathcal{R}} u_{r,k} \mathbf{h}_{rk} \mathbf{\Phi}_r \mathbf{G}_r \right) \mathbf{f}_j \right|^2 + \sigma^2}, \quad (5) $$

where $\sigma^2$ denotes the effective noise variance.

<footnote n="1">¹Several research papers can be found in the literature with focus on channel estimation for RIS-aided channels [41]–[47]. Therefore, we do not focus on these aspects in the present paper.</footnote>
## B. Channel Model

In this paper, we consider the typical Saleh-Valenzuela mmWave channel model [48] with a limited number of scattering paths to characterize the sparse nature of mmWave channels. We assume that the BS is equipped with a uniform linear array (ULA) whose array response vector is given as

$$ \mathbf{a}_{\text{ULA}}(\phi) = \frac{1}{\sqrt{N}} \left[ e^{-j \frac{2\pi d}{\lambda} \phi_i} \right]_{i \in \mathcal{I}(N)}, \quad (6) $$

$$ \mathcal{I}(N) = \{n - (N - 1)/2, n = 0, 1, \dots, N - 1\}, \quad (7) $$

where $\lambda$ and $d$ represent the signal wavelength and the antenna spacing, respectively. Besides, the number of antennas is $N$ and $\phi$ denotes the angle of arrival (AoA) or the angle of departure (AoD).

Each RIS is modeled as a uniform planar array (UPA) with $M = M_y \times M_x$ elements. The corresponding array response vectors are given by

$$ \mathbf{a}_{\text{UPA}}(\phi^1, \phi^2) = \mathbf{a}_{\text{ULA}}^{\text{x}}(\phi^1) \otimes \mathbf{a}_{\text{ULA}}^{\text{y}}(\phi^2), \quad (8) $$

where $\phi^1$ and $\phi^2$ denote the azimuth and elevation angles, respectively, and $\mathbf{a}_{\text{ULA}}^{\text{x}}(\phi^1)$ and $\mathbf{a}_{\text{ULA}}^{\text{y}}(\phi^2)$ are defined as $\mathbf{a}_{\text{ULA}}(\phi)$. The array element spacings of the considered ULA and UPA are assumed to be $\lambda/2$. Due to the severe path loss, the transmit power of two or more reflections can be ignored [49]. Hence, according to the considered array response vectors, the effective channel $\mathbf{h}_{ik}$ for the RIS$_i$-user$_k$ and $\mathbf{G}_i$ for the BS-RIS$_i$ can be represented as

$$ \mathbf{h}_{ik} = \sum_{l_h=1}^{L_h} \beta_{l_h} \mathbf{a}_r^H \left( \phi_{l_h}^{r,1}, \phi_{l_h}^{r,2} \right), \quad (9) $$

$$ \mathbf{G}_i = \sum_{l_g=1}^{L_g} \beta_{l_g} \mathbf{a}_r \left( \phi_{l_g}^{r,1}, \phi_{l_g}^{r,2} \right) \mathbf{a}_b^H \left( \phi_{l_g}^b \right), \quad (10) $$

where $L_g$ ($L_h$) denotes the number of paths, including one line-of-sight (LoS) path and $L_g(L_h) - 1$ non-line-of-sight (NLoS) paths for the cascaded links $\mathbf{G}_i$ ($\mathbf{h}_{ik}$). Similarly, $\beta_{l_g}$ ($\beta_{l_h}$) denotes the complex path gain of the $l_g$ ($l_h$) path of the corresponding channels. Note that mmWave channels normally consist of only a small number of dominant multi-path components [50], [51], while the scattering at sub-6 GHz is generally rich. Here, $\phi_{l_g}^b$ denotes the AoD at the BS while $\phi_{l_h}^{r,2}$ and $\phi_{l_g}^{r,2}$ denote the elevation angles of the cascaded links $\mathbf{h}_{ik}$ and $\mathbf{G}_i$, respectively. In addition, $\phi_{l_h}^{r,1}$ and $\phi_{l_g}^{r,1}$ represent the azimuth of the cascaded links $\mathbf{h}_{ik}$ and $\mathbf{G}_i$, respectively.

## C. Problem Formulation

Our main objective is to maximize the WSR by designing the beamforming matrix $\mathbf{F}$, the RIS association matrix $\mathbf{U}$ as well as the beamforming matrices $\mathbf{\Phi}_i, \forall i \in \mathcal{R}$. First, the WSR of the multi-RIS aided system is given by

$$ W_R = \sum_{k=1}^{K} w_k W_{R_k} \quad (11) $$

$$ = \sum_{k=1}^{K} w_k \log_2(1 + \text{SINR}_k), $$

where $w_k > 0$ denotes the weight for user$_k$ with $\sum_{k=1}^K w_k = 1$. Subsequently, the WSR maximization problem can be formulated as

$$ \text{P1:} \quad \max_{\mathbf{F}, \mathbf{U}, \mathbf{\Phi}_i, \forall i \in \mathcal{R}} \sum_{k=1}^{K} w_k \log_2(1 + \text{SINR}_k) \quad (12a) $$

$$ \text{s.t.} \quad \|\mathbf{F}\|_{\mathcal{F}}^2 \le P_{\max} \quad (12b) $$

$$ |\mathbf{\Phi}_{i,m}| = 1, \forall m = 1, 2, \dots, M \quad (12c) $$

$$ \sum_{i \in \mathcal{R}} u_{k,i} = 1, \forall k \in \mathcal{K} \quad (12d) $$

$$ u_{k,i} \in \{0, 1\}, \forall k \in \mathcal{K}, \forall i \in \mathcal{R}. \quad (12e) $$

where (12a) denotes the optimization objective for maximizing the WSR, and achieving this goal depends on the beamforming matrix $\mathbf{F}$, phase configuration matrices $\mathbf{\Phi}_i, \forall i \in \mathcal{R}$ and the RIS association matrix $\mathbf{U}$. The constraint (12b) imposes that the maximum transmit power of the BS is $P_{\max}$, while the constraint (12c) ensures that the element in each RIS satisfies the unit-modulus requirement. In addition, the constraints (12d) and (12e) are related to the RIS association scheme, where (12d) requires each user to choose the most appropriate RIS as the cascaded link, and (12e) accounts for the binary nature of the element $u_{k,i}$ in $\mathbf{U}$.

Solving problem (12) is computationally challenging because it is a non-convex problem involving both discrete variables (the RIS association matrix) and continuous variables (the reflection coefficients at the RISs and the beamforming matrix at the BS). In addition, it is difficult to decouple the optimization variables as they are closely interrelated. To tackle this problem, a GNN model is utilized as a powerful function approximator [28] to learn the transmission scheme from the CSI data. The structure of the proposed heterogeneous GNN will be introduced in detail in the next section.

## III. WEIGHTED SUM RATE MAXIMIZATION VIA HETEROGENEOUS GRAPHS NEURAL NETWORKS

In order to solve the rate maximization problem P1 via a heterogeneous GNN, we need to map the problem onto a GNN. Then, a message-passing paradigm is utilized to learn the optimal policy based on the channel state information.

### A. Graph Modeling of Multi-RIS Aided Networks

To solve P1, we model the considered network as a graph to solve P1. A basic GNN can be formally characterized by a tuple $\mathcal{G} = (\mathcal{V}, \mathcal{E})$, where $\mathcal{V}$ and $\mathcal{E}$ are sets containing nodes and edges, respectively. A node $i \in \mathcal{V}$ represents an entity, and an edge $(i, j) \in \mathcal{E}$ defines a directed relation from node $i$ to node $j$. Denote the neighboring set of node $i$ as $\mathcal{N}_i = \{j \in \mathcal{V} | (j, i) \in \mathcal{E}\}$. The attributes of node $i$ and edge $(i, j)$ are characterized by $\mathbf{v}_i$ and $\mathbf{e}_{ij}$, respectively.

Specifically, the optimization problem P1 encompasses a variety of variables each bearing distinct characteristics, including $\mathbf{F}$ related to the beamforming design of the BS, $\mathbf{\Phi}_i$ related to the phase configuration of the RIS$_i$, and $\mathbf{U}$ related to the association scheme of the users. To effectively manage this diversity, we introduce a heterogeneous GNN model
<figure><img src="image_4.png" /><figcaption>Fig. 2. The graph representation in the considered heterogeneous GNN, where the black solid edge indicates that information can be exchanged bidirectionally, and the blue dotted line indicates the directional output to the BS node.</figcaption></figure>

[40], which is adept at encapsulating a wealth of information through structural relations among many types of nodes and the unstructured content inherent to each node, as illustrated in Fig. 2. A heterogeneous graph $\mathcal{G} = (\mathcal{V}, \mathcal{E}, \mathcal{O}_V, \mathcal{R}_\mathcal{E})$, it is characterized by two additional sets, which are the predefined type of nodes $\mathcal{O}_V$ and type of nodes $\mathcal{R}_\mathcal{E}$ with $|\mathcal{O}_V| + |\mathcal{R}_\mathcal{E}| > 2$ [30]. The set $\mathcal{O}_V$ includes the BS node $\mathbf{v}_b$, RIS$_i$ node $\mathbf{v}_{r_i}$ and user$_k$ node $\mathbf{v}_{u_k}$.

According to the considered wireless communication scenario, the attributes of user$_k$ node include the weight $w_k$, noise variance $\sigma^2$ as well as the RIS association matrix $\mathbf{u}_k$. Similarly, the attributes of the BS node and the RIS$_r$ node include the beamforming matrix $\mathbf{F}$ and the phase configurations $\mathbf{\Phi}_i$, respectively. In addition, the edges in the network are only used to represent the connection between two nodes without attributes. Therefore, the heterogeneous GNN utilizes only one type of edge, leading to the omission of the types of edges set $\mathcal{R}_\mathcal{E}$ in this study. Given the heterogeneous GNN $\mathcal{G} = (\mathcal{V}, \mathcal{E}, \mathcal{O}_V)$, the objective of this paper is to find a policy $p_\Theta(\cdot)$ with parameters $\Theta$, mapping the graph $\mathcal{G}$ to the optimal $\mathbf{F}, \mathbf{\Phi}_i, \mathbf{U}$ in multi-RIS aided communication systems.

## B. The Initialization Layer of the Heterogeneous GNN

The initial layer of the heterogeneous GNN takes the CSI data as the input features for the whole heterogeneous GNN. In order to better learn the phase configuration of each RIS $\mathbf{\Phi}_i$, we extract the phase shift matrix $\theta_i$ and the equivalent channel $\mathbf{h}_k$ can be expressed as

$$ 
\begin{aligned}
\mathbf{h}_k &= \sum_{i \in \mathcal{R}} u_{k,i} \theta_i \mathrm{diag}(\mathbf{h}_{ik}) \mathbf{G}_i \quad (13) \\
&= \sum_{i \in \mathcal{R}} u_{k,i} \theta_i \mathbf{H}_{ik},
\end{aligned}
$$

where $\mathbf{H}_{ik} \in \mathbb{C}^{M \times N_t}$ denotes the cascaded channel and $\theta_i = [\theta_{i,1}, \dots, \theta_{i,M}] \in \mathbb{C}^{1 \times M}$ represents the phase coefficient of RIS$_i$. Thus, the CSI $\mathbf{H}_{ik}$ from all the RISs for each user $k$ is concatenated as $\mathbf{H}_{ik} \in \mathbb{C}^{RM \times N_t}$, which is considered as the input feature $\mathbf{v}_{u_k}^{[0]}$ for each user node, i.e., $\mathbf{v}_{u_k}^{[0]} = [\mathrm{Re}(\mathbf{H}_{ik}), \mathrm{Im}(\mathbf{H}_{ik})]$. Given the input feature representation $\mathbf{v}_{u_k}^{[0]}$, we use an encoder block $f_{\Theta}^{[1]}(\cdot)$, consisting of one layer of fully connected neural networks and a Leaky ReLU activation

function [52], to produce $\mathbf{v}_{u_k}^{[1]}$ for the user$_k$ node and $\mathbf{v}_{r_i}^{[1]}$ for the RIS$_i$ node.

In a multi-RIS multi-user communication network, each user experiences interference from the other users. In order to analyze the interference, the channels of the interfering links are aggregated to update the features through the encoder block $f_{\Theta_u}^{[1]}(\cdot)$, which is given as

$$ \mathbf{v}_{u_k}^{[1]} = f_{\Theta_u}^{[1]} ([\mathbf{v}_{u_1}^{[0]}, \dots, \mathbf{v}_{u_k}^{[0]}, \dots, \mathbf{v}_{u_K}^{[0]}]), \quad k' \in \mathcal{K}, k' \neq k. \quad (14) $$

For each RIS node, all the users receive the signals reflected through the RIS, so all the channels of the users contain an equal amount of information about the RIS. The update for the RIS node by the encoder block $f_{\Theta_r}^{[1]}(\cdot)$ is defined as

$$ \mathbf{v}_{r_i}^{[1]} = f_{\Theta_r}^{[1]} ([\mathbf{v}_{u_1}^{[0]}, \dots, \mathbf{v}_{u_k}^{[0]}, \dots, \mathbf{v}_{u_K}^{[0]}]), \quad k \in \mathcal{K}. \quad (15) $$

## C. Message-Passing Paradigm

The graph network (GN) block is adopted as the basic computation unit over graphs [53], which includes the update functions $\phi$ and aggregation functions $\rho$, as illustrated in Fig. 3. The message-passing for the RIS nodes and user nodes are marked in green lines and purple lines, respectively. The update of the representation vector in the (T-1)-th layer is based on combining its previous representation and the aggregation of the representations from its neighboring nodes.

1) Update for RIS Nodes: The updated feature $\mathbf{v}_{r_i}^{[t]}$ for the RIS$_i$ node at step $t$ ($2 \le t \le T-1$) is obtained by combining the features from the last state $\mathbf{v}_{r_i}^{[t-1]}$ and all the user nodes $\mathbf{v}_{u_k}^{[t-1]}$, $\forall k \in \mathcal{K}$. Then, the message $\xi_{rr_i}$ from its own last state and the message $\xi_{ur_i}$ from the user nodes are given as

$$ \xi_{rr_i} = \phi_{r_i \to r_i}^{(\mathrm{up})}(\xi_{r_i}, \mathbf{v}_{r_i}^{[t-1]}), \quad (16) $$

$$ \xi_{ur_i} = \phi_{u \to r_i}^{(\mathrm{up})}(\rho_{u \to r_i}(\xi_{u_1}, \dots, \xi_{u_K}), \mathbf{v}_{r_i}^{[t-1]}), \quad (17) $$

with

$$ \xi_{\mathrm{src}} = \phi_{\mathrm{src} \to r_i}(\mathbf{v}_{\mathrm{src}}^{[t-1]}), \mathrm{src} \in \{r_i, u_1, \dots, u_K\}, \quad (18) $$

where $\phi_{\mathrm{src} \to r_i}(\cdot)$ and $\phi_{\mathrm{src} \to r_i}^{(\mathrm{up})}(\cdot)$ are two-layer multi-layer perceptrons (MLPs) for message combining and updating, respectively, with the superscript (up) denoting “update”. Besides, $\rho_{u \to r_i}$ is the mean aggregation function for all the different user nodes [54], i.e., $\rho_{u \to r_i}(\xi_{u_1}, \dots, \xi_{u_K}) = 1/K \sum_{k=1}^K \xi_{u_k}$.

Then, in the considered heterogeneous GNN, the updated RIS node $\mathbf{v}_{r_i}^{[t]}$ aggregates the feature information from different-type neighbors, which is given as

$$ \mathbf{v}_{r_i}^{[t]} = \rho_{r_i}(\xi_{rr_i}, \xi_{ur_i}) + \mathbf{v}_{r_i}^{[t-1]}, \quad (19) $$

where $\rho_{r_i}(\cdot)$ denotes the operation of mean aggregation, which aggregates the features from the RIS node and user nodes. Besides, a residual operation is added in this step to avoid the degradation and gradient explosion of the considered network [55]. Thus, the node features $\mathbf{v}_{r_i}^{[t-1]}$ at the last iteration can be retained through this residual operation.
<figure><img src="image_5.png" /><figcaption>Fig. 3. The message-passing paradigm in the heterogeneous GNN network. The purple line on the left represents the update process for the features of the user nodes and the green line on the right represents the features for the update process of the RIS nodes. In addition, the black line at the bottom denotes the output layers for the features of the final nodes.</figcaption></figure>

2) Update for User Nodes: The updated feature $\mathbf{v}_{u_k}^{[t]}$ is obtained by combining the features from all the user nodes and all the RIS nodes. In the same manner, the message $\xi_{ru_k}$ combined from the RIS nodes and the message $\xi_{uu_k}$ combined from the user nodes are given as

$$ \xi_{ru_k} = \phi_{r \to u_k}^{(\text{up})} (\rho_{r \to u_k} (\xi_{r_1}, \dots, \xi_{r_R}), \mathbf{v}_{u_k}^{[t-1]}), \quad (20) $$

$$ \xi_{uu_k} = \rho_{u_k}^{(u)}(\phi_{u \to u_k}^{(\text{up})}(\xi_u), \phi_{u_k \to u_k}^{(\text{up})}(\xi_{u_k}, \mathbf{v}_{u_k}^{[t-1]})), \quad (21) $$

where both $\phi_{src \to u_k}(\cdot)$ and $\phi_{src \to u_k}^{(\text{up})}(\cdot)$ are two-layer MLPs for the message combination and update process of $\mathbf{v}_{u_k}^{[t]}$, respectively. Moreover, $\xi_u = \rho_{u \to u_k}(\xi_{u_1}, \dots, \xi_{u_K})$ denotes the message aggregated from all the user nodes. The aggregation function $\rho_{u \to u_k}(\cdot)$ is a max aggregation function, i.e., $\rho_{u \to u_k}(\xi_{u_1}, \dots, \xi_{u_K}) = \max(\xi_{u_1}, \dots, \xi_{u_K})$, where the max-pooling is adopted as it performs well and effectively captures the strongest interference signal among the users [33]. Besides, $\xi_{src} = \phi_{src \to r_i}(\mathbf{v}_{src}^{[t-1]})$, $src \in \{r_1, \dots, r_R, u_1, \dots, u_K\}$ denotes the message from each neighbor node and $\rho_{r \to u_k}$ is the mean aggregation function.

Then, in heterogeneous GNNs, the updated user node $\mathbf{v}_{u_k}^{[t]}$ aggregates the feature information from different neighbors, which is given as

$$ \mathbf{v}_{u_k}^{[t]} = \rho_{u_k}(\xi_{ru_k}, \xi_{uu_k}) + \mathbf{v}_{u_k}^{[t-1]}, \quad (22) $$

where $\rho_{u_k}(\cdot)$ denotes the operation of mean aggregation function.

3) Output Layers for Node Features: To obtain the final feature of each node $\mathbf{v}_b^{[T]}$, $\mathbf{v}_{r_i}^{[T]}$ and $\mathbf{v}_{u_k}^{[T]}$, a decoder layer is used to output all the optimized features, which is given as

$$ \mathbf{v}_{\text{dst}}^{[T]} = \rho_{\text{dst}}^{(\text{out})}(\xi_{\text{udst}}, \xi_{\text{rdst}}), \\ \text{dst} \in \{b, r_1, \dots, r_R, u_1, \dots, u_K\}, \quad (23) $$

where $\xi_{udst} = \phi_{u \to \text{dst}}^{(\text{out})}(\mathbf{v}_{u_1}^{[T-1]}, \dots, \mathbf{v}_{u_K}^{[T-1]})$ and $\xi_{rdst} = \phi_{r \to \text{dst}}^{(\text{out})}(\mathbf{v}_{r_1}^{[T-1]}, \dots, \mathbf{v}_{r_R}^{[T-1]})$ denotes the destination node features extracted from the user nodes and RIS nodes, respectively. As shown by the black line in Fig. 3, the final state $\mathbf{v}_{\text{dst}}^{[T]}$ is obtained by a mean aggregation function $\rho_{\text{dst}}^{(\text{out})}(\cdot)$. According to the defined types of nodes already mentioned, the BS node outputs the beamforming matrix $\mathbf{F}$ as $\mathbf{v}_b^{[T]} = [\text{Re}\{\mathbf{F}\}, \text{Im}\{\mathbf{F}\}]$. Similarly, the RIS$_r$ node and the user$_k$ node output the phase configuration matrix $\bar{\mathbf{\Phi}}_i$ and the association scheme $\bar{\mathbf{u}}_k$, respectively. Specifically, $\mathbf{v}_{r_i}^{[T]} = [\text{Re}\{\mathbf{\Phi}_i\}, \text{Im}\{\mathbf{\Phi}_i\}]$ and $\mathbf{v}_{u_k}^{[T]} = \bar{\mathbf{u}}_k$. After obtaining the final state of each node, the node features are fed into the normalization layers to meet the constraints in P1 (12b) ~ (12e). Then, the WSR-related loss functions are used to update the network parameters to obtain the desired variables in P1.

Suppose that T iterations are executed during the training process, and we define each update as HGN$_t$, $\forall t \in T$. An encoder-process-decoder architecture is utilized to facilitate the training [53], including one encoder block HGN$_1$ to initialize the node features through (14) and (15), one decoder block HGN$_T$ to output the node features through (23), and one core block (a concatenation of blocks from HGN$_2$, ..., HGN$_t$, ..., HGN$_{T-1}$) to implement the message-
passing paradigm through (16) to (22). After that, we can obtain the final result for the optimal **F**, **U** and **Φ**<sub>i</sub>, ∀i ∈ ℝ.

### D. Model Architecture of the Heterogeneous GNN

In addition to the message-passing paradigm, the proposed heterogeneous GNN needs some additional layers to complete the whole learning process, and the functions of these layers are listed as follows.

1) **The Normalization Layers:** Typical neural networks operate on real numbers, hence it is not possible to learn the mapping of the channel information to the values of **F** and **Φ**<sub>i</sub>, which are all complex values in wireless communication systems. Instead, we learn the vector consisting of the real and imaginary parts of the beamforming matrix **F** and the phase configuration matrices **Φ**<sub>i</sub>, ∀i ∈ ℝ. Therefore, before calculating the WSR loss function, additional layers are necessary to merge the real-valued results into complex-valued results.

* To satisfy the constraint on the transmit power in (12b), a normalization layer is employed at the output feature of the BS node $\mathbf{v}_{\text{b}}^{[T]} = [\text{Re}\{\bar{\mathbf{F}}\}, \text{Im}\{\bar{\mathbf{F}}\}]$, and the function of this layer can be expressed as

$$ \mathbf{F} = \sqrt{P_{\max}} \frac{\text{Re}\{\bar{\mathbf{F}}\} + j \text{Im}\{\bar{\mathbf{F}}\}}{\sqrt{\|\text{Re}\{\bar{\mathbf{F}}\}\|_{\mathcal{F}}^2 + \|\text{Im}\{\bar{\mathbf{F}}\}\|_{\mathcal{F}}^2}}. \quad (24) $$

* Similarly, the phase configuration of each RIS needs to fulfill the unit-modulus constraint (12c), and a normalization layer is utilized to tackle the real and imaginary components of each reflection coefficient $\bar{\theta}_{i,m}$ in the output feature of the RIS<sub>i</sub> node $\mathbf{v}_{\text{r}_i}^{[T]} = [\text{Re}\{\bar{\mathbf{\Phi}}_i\}, \text{Im}\{\bar{\mathbf{\Phi}}_i\}]$, ∀i ∈ ℝ. Thus, the function can be represented as

$$ \theta_{i,m} = \frac{\text{Re}\{\bar{\theta}_{i,m}\} + j \text{Im}\{\bar{\theta}_{i,m}\}}{\sqrt{|\text{Re}\{\bar{\theta}_{i,m}\}|^2 + |\text{Im}\{\bar{\theta}_{i,m}\}|^2}}. \quad (25) $$

2) **The RIS Association Scheme:** Each user chooses a suitable RIS as the reflection link to enhance the link quality, as dictated by constraint (12d). Thus, we utilize a SoftMax layer on the output feature of each user node $\mathbf{v}_{\bar{u}_k}^{[T]} = \bar{u}_k$ to implement the association scheme. For each element $\bar{u}_{k,i}$ in $\bar{u}_k$, the SoftMax score can be calculated by

$$ c_{u_{k,i}} = \text{SoftMax}(\bar{u}_{k,i}) = \frac{e^{\bar{u}_{k,i}}}{\sum_{i=1}^{R} e^{\bar{u}_{k,i}}}. \quad (26) $$

After obtaining the SoftMax score $c_{u_{k,i}}$, we sort the scores of each user, i.e., each column in **U**. Leveraging the Softmax function for the RIS association scheme inherently allows for the possibility of associating multiple RISs with each user. We select the positions of the largest values and set them to 1, indicating that the user<sub>k</sub> selects the most appropriate RIS to assist the communication links. Therefore, the selection scheme is obtained

$$ u_{k,i} = \begin{cases} 1, & \text{if } c_{u_{k,i}} = \max_{i \in \mathcal{R}}(c_{u_{k,i}}), \forall k \in \mathcal{K} \\ 0, & \text{otherwise.} \end{cases} \quad (27) $$

For ease of understanding, a network deployment with 4 RISs and 5 users is considered as an example to explain the

<figure><img src="image_6.png" /><figcaption>Fig. 4. An example of the RIS association scheme in the multi-RIS multi-user system, considering 4 RISs and 5 users.</figcaption></figure>

structure of the matrix **U** in Fig. 4. It is apparent that each user selects only one RIS as the cascaded communication link, whereas each RIS can serve multiple users or none of them. For instance, $\mathbf{u}_1 = [1, 0, 0, 0]^T$ means that user<sub>1</sub> chooses RIS<sub>1</sub> as the cascaded link. Meanwhile, we notice that the third row of each user node $\mathbf{u}_k, \forall k \in \mathcal{K}$ is 0, which means that no user chooses RIS<sub>3</sub> as a cascaded link. Accordingly, RIS<sub>3</sub> is set to the 'OFF' state, which implies that it does not reflect any signals to any users.

## IV. THE TRAINING PROCESS OF THE HETEROGENEOUS GNN IN MULTI-RIS MULTI-USER SYSTEMS

Based on the message-passing paradigm and the normalization layers introduced in the previous section, we focus on maximizing the WSR in multi-RIS aided communication systems. Aiming to address the challenge of learning continuous and discrete variables in the maximization problem, a pre-training process is adopted to obtain the penalty term. Therefore, the overall training processes and the pre-training process are elaborated in the following section.

### A. The Training Process of the Heterogeneous GNN

As defined in problem P1, in addition to the beamforming matrix and phase configuration at the RISs, the heterogeneous GNN is required to output the RIS association scheme **U**. However, the association matrix is composed of binary elements, which are discrete variables. Thus, the process to determine the optimal value of **U**, **F**, and **Φ**<sub>i</sub>, ∀i ∈ ℝ is very difficult to learn. In order to better learn the RIS association matrix, it is important to have a suitable optimization criterion. Therefore, two different loss functions are adopted to solving the maximization problem P1. One is related to the WSR of the multi-RIS aided communication systems for maximizing the target of the considered problem, which is defined as

$$ \begin{aligned} \text{Loss}_1 &= -\sum_{k=1}^{K} w_k \log_2(1 + \text{SINR}_k) \\ &= -\sum_{k=1}^{K} w_k \log_2\left(1 + \frac{|\mathbf{h}_k \mathbf{f}_k|^2}{\sum_{j \in \mathcal{K}, j \neq k} |\mathbf{h}_k \mathbf{f}_j|^2 + \sigma^2}} \right). \end{aligned} \quad (28) $$
<figure><img src="image_7.png" /><figcaption>Fig. 5. The overall structure of the proposed heterogeneous GNN. The left part is the message-passing paradigm, which has been illustrated in Section III. The orange rectangle identifies the pre-training process for the maximization problem with the predetermined RIS association scheme.</figcaption></figure>

**Algorithm 1** The overall process of the heterogeneous GNN

**Input:**

The heterogeneous graph $\mathcal{G} = (\mathcal{V}, \mathcal{E}, \mathcal{O}_\mathcal{V})$;

The training and testing samples $\mathcal{S} = \{\mathbf{H}_{ik}\}, \forall i \in \mathcal{R}, k \in \mathcal{K}$.

**Output:**

The final node embedding feature $\mathbf{F}, \mathbf{\Phi}_i, \forall i \in \mathcal{R}$, and $\mathbf{U}$;

1: **for** epoch = 1, ..., $N_e$ **do**

2: Initialize the node features via the encoder block HGN0 from samples for all the RIS nodes and user nodes using (14) and (15);

3: **for** t = 1, 2, ..., T - 1 **do**

4: **for** all the user nodes and RIS nodes **do**

5: Combine the message from the different source node features for all the RIS nodes and user nodes using (16)-(18) and (20)-(21);

6: Aggregate the message from different-type neighbors for all the RIS nodes and user nodes using (19) and (22);

7: **end for**

8: Update the node feature $\mathbf{v}_p^{[t+1]}$, $\forall p \in \mathcal{O}_\mathcal{V}$ using (23);

9: **end for**

10: Calculate the loss function using (30):

$$ \text{Loss} = - \sum_{k=1}^{K} \left( w_k \log_2(1 + \text{SINR}_k) + \eta \mathbf{u}_k^{(1)} \log \mathbf{c}_{\text{uk}} \right); $$

11: Perform backpropagation and update the parameters of the heterogeneous GNN;

12: **end for**

13: **return** $\mathbf{F}, \mathbf{\Phi}_i, \forall i \in \mathcal{R}$, and $\mathbf{U}$.

Then, to achieve more efficient matching based on the nearest-distance principle, we employ the CE loss function to guide the learning process, defined as

$$ \text{Loss}_2(\mathbf{u}_k^{(1)}, \mathbf{c}_{u_k}) = - \sum_{k=1}^{K} \mathbf{u}_k^{(1)} \log \mathbf{c}_{u_k}, \quad (29) $$

where $\mathbf{u}_k^{(1)}$ and $\mathbf{c}_{u_k}$ denote the distance-dependent label and the output for the user $k$, respectively. The process of obtaining the distance-dependent label $\mathbf{u}_k^{(1)}$ requires the calculation of the distances between the users and RISs using their location

information from the dataset. For each user $k$, we identify the nearest RIS, denoted by the index $\hat{r}$. In the label vector $\mathbf{u}_k^{(1)}$, we set the element corresponding to $\mathbf{u}_{k, \hat{r}}$ to 1, which indicates the association with the nearest RIS, while all the other elements are set to 0, which implies no association.

Based on the output variables and the CSI information, the loss function is defined as a weighted sum of Loss1 and Loss2, which can be expressed as

$$ \text{Loss} = \text{Loss}_1 + \eta \text{Loss}_2 \\ = - \sum_{k=1}^{K} \left( w_k \log_2(1 + \text{SINR}_k) + \eta \mathbf{u}_k^{(1)} \log \mathbf{c}_{\text{uk}} \right). \quad (30) $$

The first term in the loss function, Loss1, optimizes the overall system's weighted sum rate, while the second term, Loss2, is a distance-dependent loss that encourages the selection of the nearest RIS. The coefficient $\eta$ is a penalty term, and its role is to strike a balance between the two components of the loss function, which can be calculated by

$$ \eta = \frac{\text{WSR}_p}{\text{WSR}_{p_0}}, \quad (31) $$

where WSRp is the result in the pre-training process by using the loss function in (28) and WSRp0 is the result under the distance-dependent label $\mathbf{u}_k^{(1)}$ when $P_{\text{max}} = 30$ dBm. For ease of understanding, the overall process of the heterogeneous GNN for solving P1 is shown as Algorithm 1 and the pre-training process is given in the next subsection.

## B. The Pre-training Process for the Penalty Term

To obtain WSRp, we introduce a pre-training process where the learning is based on the distance-dependent label $\mathbf{u}_k^{(1)}$, assuming all users connect to their nearest RIS already. Thus, the message-passing paradigm remains the same as illustrated in Section III-C, but without returning $\mathbf{U}$, as the orange rectangle in Fig. 5. Based on the predetermined $\mathbf{u}_k^{(1)}$, the cascaded channel $\mathbf{h}_k^{(1)}$ can be rewritten as

$$ \mathbf{h}_k^{(1)} = \sum_{i \in \mathcal{R}} u_{k,i}^{(1)} \theta_i^{(1)} \mathbf{H}_{ik}, \quad (32) $$

where $\mathbf{H}_{ik} = \text{diag}(\mathbf{h}_{ik})\mathbf{G}_i$ is the input of the heterogeneous GNN. Therefore, the WSR maximization problem, given the
**Algorithm 2** The pre-training process

**Input:**

A predetermined scheme $\mathbf{U}^{(1)} = [\mathbf{u}_1^{(1)}, \dots, \mathbf{u}_K^{(1)}]$;

The heterogeneous graph $\mathcal{G} = (\mathcal{V}, \mathcal{E}, \mathcal{O}_\mathcal{V})$;

The training and testing samples $\mathcal{S} = \{\mathbf{H}_{ik}\}, \forall i \in \mathcal{R}, k \in \mathcal{K}$.

**Output:**

The final node embedding feature $\mathbf{F}^{(1)}$ and $\theta_i^{(1)}, \forall i \in \mathcal{R}$.

1: **for** epoch = 1, ..., $N_e$ **do**

2: Initialize the node features via the encoder block $\text{HGN}_0$ as $\mathbf{v}_{\text{src}}^{[1]}$ from the training sample set for all the RIS nodes and user nodes using (14) and (15);

3: **for** $t = 2, \dots, T-1$ **do**

4: **for** all the user nodes and RIS nodes **do**

5: Combine the message from the different source node features for all the RIS nodes and user nodes using (16)-(18) and (20)-(21);

6: Aggregate the message from different-type neighbors for all the RIS nodes and user nodes using (19) and (22);

7: **end for**

8: Update the features of the BS node $\mathbf{v}_b^{[T]}$ and RIS nodes $\mathbf{v}_{r_i}^{[T]}, \forall i \in \mathcal{R}$ using (23);

9: **end for**

10: Calculate the WSR loss using (28):
$$ \text{Loss} = - \sum_{k=1}^{K} w_k \log_2(1 + \text{SINR}_k^{(n)}) $$

11: Perform backpropagation and update parameters in the heterogeneous GNN;

12: **end for**

13: **return** $\mathbf{F}^{(1)}$ and $\theta_i^{(1)}, \forall i \in \mathcal{R}$.

predetermined association scheme $\mathbf{u}_k^{(1)}$, aims to find the optimal $\mathbf{F}^{(1)}$ and $\theta_i^{(1)}$, which can be expressed as

$$ \text{P2: } \max_{\mathbf{F}^{(1)}, \theta_i^{(1)}, \forall i \in \mathcal{R}} \sum_{k=1}^{K} w_k \log_2(1 + \text{SINR}_k^{(1)}) \quad (33a) $$

$$ \text{s.t.} \quad \left\| \mathbf{F}^{(1)} \right\|_{\mathcal{F}}^{2} \leq P_{\max} \quad (33b) $$

$$ \left| \theta_{i,m}^{(1)} \right| = 1, \forall m = 1, 2, \dots, M. \quad (33c) $$

For ease of understanding, the overall pre-training process of the heterogeneous GNN for solving the maximization problem P2 is shown as Algorithm 2, i.e. the training process in the orange rectangle in Fig. 5.

## V. NUMERICAL RESULTS

### A. Simulation Setup

In this section, numerical simulations are provided to verify the effectiveness of the proposed heterogeneous GNN in Table I. The BS is assumed to be located at the origin, and two RISs are positioned at coordinates $[30m, 25m]$ and $[30m, -25m]$. Additionally, two users are randomly located within a rectangular region defined by the range $[40m: 50m, -25m: 25m]$. According to [48], the number of paths $L$ is 3, including one LoS path and $L-1$ NLoS paths. Accordingly, the channel gain for the LoS path is $\beta_{LoS} \sim$

<table><caption>TABLE I
SIMULATION AND TRAINING HYPERPARAMETERS</caption><thead><tr><td>Parameter</td><td>Value</td><td>Parameter</td><td>Value</td></tr></thead><tbody><tr><td>Reflecting elements (M)</td><td>16</td><td>Transmit antenna elements (N<sub>t</sub>)</td><td>8</td></tr><tr><td>Number of users (K)</td><td>2</td><td>Number of RISs (R)</td><td>2</td></tr><tr><td>Noise power (σ<sub>n</sub><sup>2</sup>)</td><td>-85 dBm</td><td>Initial learning rate</td><td>5 × 10<sup>-4</sup></td></tr><tr><td>Weighting decay</td><td>2.5 × 10<sup>-5</sup></td><td>Epochs</td><td>400</td></tr><tr><td>Batch size</td><td>128</td><td>Hidden size</td><td>512</td></tr></tbody></table>

$\mathcal{CN}(0, 10^{-0.1\text{PL}(r)})$ where $\text{PL}(r) = \varrho_a + 10\varrho_b\log_{10}(r) + \xi$ with $\xi \sim \mathcal{N}(0, \sigma_\xi^2)$. In addition, the channel gain for the NLoS path is $\beta_{NLoS} = 0.01\beta_{LoS}$ and the channel realizations are obtained by setting $\varrho_a = 61.4$, $\varrho_b = 2$, and $\sigma_\xi = 5.8$ dB [9]. Without loss of generality, the weights $w_k$ are set to be equal. The proposed heterogeneous GNN is trained and tested separately for each SNR point to evaluate its performance robustly across varying conditions. A total of 110,000 samples, split into 100,000 for training and 10,000 for validation, are utilized. As for the training process, the Adam optimizer [56] is adopted and an NVIDIA RTX 2080Ti GPU is utilized.

### B. Baseline Case Studies

We introduce three case studies to show the significance of selecting a suitable RIS for the users, i.e., considering the scenario when the RIS association scheme is predefined and the heterogeneous GNN just needs to output the beamforming matrix and phase configurations of the RISs. The cases are listed as follows

* **Case 1:** All the users are connected to the same RIS, i.e., the system reduces to a single RIS-aided system. Therefore, the association matrix $\mathbf{U}$ is reduced to the vector $\mathbf{u}^{(c_1)}$ and does not need to be optimized. In this case, the message-passing paradigm is illustrated in Section III-B, but without returning $\mathbf{u}^{(c_1)}$. Therefore, the cascaded channel $\mathbf{h}_k^{(c_1)}$ in the single RIS-aided communication system can be rewritten as

$$ \begin{aligned} \mathbf{h}_k^{(c_1)} &= u_k^{(c_1)} \mathbf{h}_{1k} \mathbf{\Phi}_1^{(c_1)} \mathbf{G}_1 \\ &= u_k^{(c_1)} \theta_1^{(c_1)} \mathbf{H}_{1k}, \end{aligned} \quad (34) $$

where $\mathbf{H}_{1k} = \text{diag}(\mathbf{h}_{1k})\mathbf{G}_1$ is the input of the heterogeneous GNN. In the single RIS-aided scenario, $\mathbf{u}^{(c_1)}$ is an all-ones vector, so the channel can be written as $\mathbf{h}_{1k} = \theta_1^{(c_1)}\mathbf{H}_{1k}$ with $\theta_1^{(c_1)} = [\theta_1^{(c_1)}, \dots, \theta_{1,M}^{(c_1)}]$. Therefore, the WSR maximization problem aims to learn the beamforming matrix $\mathbf{F}^{(c_1)}$ and the phase vector $\theta^{(c_1)}$, which is reformulated as

$$ \text{P3: } \max_{\mathbf{F}^{(c_1)}, \theta_1^{(c_1)}} \sum_{k=1}^{K} w_k \log_2(1 + \text{SINR}_k^{(c_1)}) \quad (35a) $$

$$ \text{s.t.} \quad \left\| \mathbf{F}^{(c_1)} \right\|_{\mathcal{F}}^{2} \leq P_{\max} \quad (35b) $$

$$ \left| \theta_{1,m}^{(c_1)} \right| = 1, \forall m = 1, 2, \dots, M. \quad (35c) $$

Then, according to the same message passing paradigm and the normalization layers in the orange rectangle in
Fig. 5, the optimal $\mathbf{F}^{(c_1)}$ and $\theta_1^{(c_1)}$ with the WSR loss function in (28).

* **Case 2:** All the users are connected to the nearest RIS, then the heterogeneous GNN outputs $\mathbf{F}^{(c_2)}$ and $\theta_i^{(c_2)}, \forall i \in \mathcal{R}$. Thus, the training process is the same as the pre-training process in the last section and is not illustrated again here.

* **Case 3:** All the users are connected to the farthest RIS, which is considered as the worst case. Then the heterogeneous GNN outputs $\mathbf{F}^{(c_3)}$ and $\theta_i^{(c_3)}, \forall i \in \mathcal{R}$. Similar to Case 2, the phase configuration matrices at multiple RISs need to be optimized, and the cascaded link is expressed as $\mathbf{h}_k^{(c_3)} = \sum_{i \in \mathcal{R}} v_{k,i}^{(c_3)} \theta_i^{(c_3)} \mathbf{H}_{ik}$. Also, the loss function is the same as in (28) by replacing the channel information with $\mathbf{h}_k^{(c_3)}$.

### C. Baseline Algorithms

For performance comparison, we consider two typical benchmarks and one deep learning benchmark, including the weighted minimum mean-square error (WMMSE) algorithm with random phases, the AO algorithm for the joint optimization of the BS beamforming matrix at the BS and phase configurations at the RISs, and a DNN algorithm that also employs an encoder-processing-decoder framework to intelligently optimize the BS beamforming matrix at the BS, the phase configurations at the RISs and the RIS association scheme. These methods are illustrated as follows.

* **Benchmark 1 (Random Phase):** The WMMSE [57] is utilized to optimize the beamforming matrix with a random phase configuration at each RIS, which has been studied extensively in the literature. Then, the iterative updating rule is given as

$$ \chi_k = \left( \sum_{i=1}^{K} |\mathbf{h}_k^{\mathrm{H}} \mathbf{f}_i|^2 + \sigma_0^2 \right)^{-1} \mathbf{h}_k^{\mathrm{H}} \mathbf{f}_k, \quad (36) $$

$$ \kappa_k = (1 - \chi_k^* \mathbf{h}_k^{\mathrm{H}} \mathbf{f}_k)^{-1}, \quad (37) $$

$$ \mathbf{f}_k = \omega_k \chi_k \kappa_k \left( \lambda \mathbf{I}_M + \sum_{i=1}^{K} \omega_i |\chi_i|^2 \kappa_i \mathbf{h}_i \mathbf{h}_i^{\mathrm{H}} \right)^{-1} \mathbf{h}_k, \quad (38) $$

where $\lambda \ge 0$ is the optimal dual variable for the transmit power constraint. The complexity of the WMMSE algorithm is $O(I_\lambda I_f K N_t^3)$, where $I_\lambda$ and $I_f$ are the numbers of iterations for searching $\lambda$ and the three-step updating loop, respectively.

* **Benchmark 2 (AO):** The AO algorithm is utilized as a performance benchmark, which jointly optimizes the beamforming matrix at the BS and phase configurations at the RISs [58]. Based on the AO method, the joint optimization problem is decomposed into two subproblems: one is the conventional beamforming design problem for $\mathbf{F}$ at the BS by using the WMMSE algorithm, and the other is the phase optimization problem for $\mathbf{\Phi}$ according to given optimized beamforming vectors by using the Riemannian conjugate gradient algorithm (RCG) [59]. Thus, the total complexity of the alternating optimization

<figure><img src="image_8.png" /><figcaption>Fig. 6. Illustration of the intelligent DNN benchmark.</figcaption></figure>

approach is $O(I_O(I_\lambda I_f K N_t^3 + I_R K^2 M^2))$, where $I_O$ and $I_R$ denote the number of iterations of the outer loop, and the number of iterations of the inner RCG algorithm, respectively.

* **Benchmark 3 (Intelligent DNN):** To validate the advantages of the message-passing scheme for the considered heterogeneous GNN, we compare it with an intelligent DNN benchmark that also employs an encoder-processing-decoder framework. Similar to the proposed heterogeneous GNN, the intelligent DNN begins by concatenating all the CSI as inputs. An encoder block with one MLP and a Leaky ReLU layer extracts the features, which are then refined through a processing block that is constituted by four MLPs and Leaky ReLU layers. A decoder block outputs the desired features for problem P1, as shown in Fig. 6. For a fair comparison, the size of the MLP layers of the intelligent DNNs is the same as that of the GNNs.

### D. Performance Versus $P_{\max}$

Fig. 7 presents the WSR as a function of the maximum transmit power $P_{\max}$ of the BS in a multi-RIS aided system. Specifically, the performance of three different approaches is compared: the AO algorithm and the heterogeneous GNN algorithm (optimizing the beamforming matrix $\mathbf{F}$ and phase shifts $\theta_i, \forall i \in \mathcal{R}$), and WMMSE with random phase under three distinct baseline case studies, as shown in Fig. 7(a)-(c). From the comparison of these subfigures, it is evident that the WSR increases with rising $P_{\max}$ across all cases. However, different RIS association schemes lead to significant performance variations depending on the proximity of the RIS to the associated users. As far as Case 2 in Fig. 7(b) is concerned, the WSR consistently outperforms the other case studies. Conversely, when users are connected to RISs located far away as per Case 3, the performance degrades due to the more severe path loss. Specifically, the RIS association scheme considered in Case 2 offers a notable performance improvement over Case 1 (18.2%) and Case 3 (63.4%), highlighting the importance of selecting an optimal RIS association strategy.
<figure><img src="image_9.png" /><figcaption>Fig. 7. The WSR versus $P_{max}$ with $N_t = 8$ and $M = 16$.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 8. The WSR versus $P_{max}$ with $N_t = 8$ and $M = 16$.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 9. The WSR versus $M$ with $P_{max} = 20$ dBm and $N_t = 8$.</figcaption></figure>

Fig. 8 presents a performance comparison across different network architectures, including the combination of three case studies with the heterogeneous GNN, the intelligent heterogeneous GNN, and the intelligent DNN. Here, 'intelligent' indicates that the association scheme U is learned by the neural network simultaneously with the beamforming matrix F and phase configurations $\theta_i, \forall i \in \mathcal{R}$, as opposed to the two-step methods utilized by the baseline case studies. As shown in Fig. 8, the intelligent heterogeneous GNN outperforms the other architectures, indicating that this intelligent learning approach yields significant improvements due to the coupled optimization of all variables. Specifically, the intelligent heterogeneous GNN demonstrates a 10.7% improvement compared with Case 2 + heterogeneous GNN, and the intelligent heterogeneous GNN is shown to promote the association of each user to its nearest RIS in light of the multiplicative path loss effects discussed in [60], as guided by the CE loss function in (29) during training. Additionally, it is evident that the intelligent DNN performs the worst, primarily due to its inability to effectively account for the underlying network structure.

## E. Performance Versus $M$

<table><caption>TABLE II
COMPUTATION TIME COMPARISON FOR DIFFERENT ALGORITHMS WITH $P_{max} = 20$ dBm AND $M = 16$ [SEC]</caption><thead><tr><th></th><th>Algorithm</th><th>$N_t = 4$</th><th>$N_t = 8$</th><th>$N_t = 16$</th></tr></thead><tbody><tr><td rowspan="2">CPU</td><td>Case 2 + AO</td><td>1.2381</td><td>1.4633</td><td>2.0497</td></tr><tr><td>Intelligent Heterogeneous GNN</td><td>4.984 E-2</td><td>5.248 E-2</td><td>5.351 E-2</td></tr><tr><td rowspan="2">GPU</td><td>Case 2 + Heterogeneous GNN</td><td>2.297 E-2</td><td>2.295 E-2</td><td>2.389 E-2</td></tr><tr><td>Intelligent Heterogeneous GNN</td><td>2.244 E-2</td><td>2.250 E-2</td><td>2.374 E-2</td></tr></tbody></table>

Fig. 9 shows the WSR versus the number of reflecting elements at each RIS $M$. We see that the WSR improves monotonically as $M$ increases. Besides this, the proposed heterogeneous GNN structure can tackle the optimization problem P1 for a different number of reflecting elements ($M = 9, 16, 25$), showcasing the generality and effectiveness of the proposed architecture. Then, comparing the results of random phase under Case 1, we see that $M$ has little influence
<figure><img src="image_12.png" /><figcaption>Fig. 10. The WSR versus $N_t$ with $P_{\text{max}} = 20$ dBm and $M = 16$.</figcaption></figure>

on the performance because the random phase does not optimize the phase configurations of the RISs. In general, the RIS association scheme results in a performance gain of about 32.3% for all values of $M$.

### F. Performance Versus $N_t$

Fig. 10 illustrates the WSR versus the total number of antennas $N_t$. It is evident that the WSR increases with the number of antennas and the proposed intelligent heterogeneous GNN consistently illustrates strong performance across different values of $N_t$. Moreover, the results show that selecting an appropriate RIS for each user significantly enhances performance. The performance improvement ranges from approximately 20% to 35% depending on the number of antennas, highlighting the effectiveness of the intelligent heterogeneous GNN in selecting the suitable association scheme.

### G. Complexity and Convergence of the Heterogeneous GNN

The complexity of different algorithms and the convergence of the proposed intelligent heterogeneous GNN are illustrated in Table II and Fig. 11. In Table II, the running time of different methods for different values of $N_t$ are reported.

We executed the AO algorithm considering Case 2 and the intelligent heterogeneous GNN on an 11th Gen Intel Core i7-1165G7 @ 2.80GHz CPU, as well as on an NVIDIA RTX 2080Ti GPU. Running times were evaluated across 128 problem instances, corresponding to the batch size utilized during training. The intelligent heterogeneous GNN run on the GPU utilizes parallel computing, unlike the AO algorithm, which is sequential by nature. As $N_t$ increases, the running time of the AO algorithm increases significantly. Notably, the computation time of the GPU for the intelligent heterogeneous GNN, including the additional loss function for user-RIS association learning, is comparable to the heterogeneous GNN under Case 2, indicating minimal impact on time complexity. On the CPU, the running time of the GNN is about twice that of the GPU but still only 2.5% of the time required by the AO algorithm. Additionally, the intelligent heterogeneous

<figure><img src="image_13.png" /><figcaption>Fig. 11. The testing WSR versus training epochs with $N_t = 8$, $M = 16$ and $P_{\text{max}} = 30$ dBm.</figcaption></figure>

GNN comprises roughly 5.8 million parameters, a moderate count that underscores the network's lower complexity and supports the running time analysis as a metric for validation. In conclusion, the comparison of running times underscores the advantage of the intelligent heterogeneous GNN, which achieves performance comparable to the AO algorithm under Case 2 while significantly reducing the computation time.

Then we analyze the test results. Specifically, we plot the testing WSR evaluated on different numbers of training samples against the training epochs, to evaluate the convergence of the proposed intelligent heterogeneous GNN in Fig. 11. The intelligent heterogeneous GNN is trained for 400 epochs under different numbers of training samples. We see that the WSR initially increases sharply with the training epoch but it then converges to a stable value after 50 training epochs. Therefore, a few training epochs are sufficient to achieve satisfactory performance. If $1 \times 10^5$ training samples are considered, for instance, greater than 93% of the WSR is achieved with only 50 training epochs.

### H. Performance Versus the Training Parameters

Next, we present the WSR performance as a function of the number of training samples and analyze the impact of different model sizes for the proposed intelligent heterogeneous GNNs and intelligent DNNs. In contrast to the baseline case studies that utilize a two-step method—initially determining the association scheme $\mathbf{U}$ and subsequently calculating the beamforming matrix $\mathbf{F}$ and phase configurations $\theta_i, \forall i \in \mathcal{R}$, both intelligent neural networks output $\mathbf{U}, \mathbf{F}$, and $\theta_i, \forall i \in \mathcal{R}$ simultaneously. As shown in Fig. 12, the performance of the intelligent GNN improves with the number of training samples and the size of the layers. Notably, the difference in performance between the intelligent heterogeneous GNNs with hidden sizes of 512 and 1024 is minor, indicating that a hidden size of 512 is sufficient for the considered optimization problem. The intelligent heterogeneous GNNs significantly outperform the intelligent DNNs. Assuming $2 \times 10^5$ training
<figure><img src="image_14.png" /><figcaption>Fig. 12. The WSR versus the training samples with $N_t = 8$, $M = 16$ and $P_{max} = 30$ dBm.</figcaption></figure>

samples, for example, the intelligent GNNs perform 10 times better than the DNNs. Conversely, increasing the hidden sizes or training samples does not enhance the performance of DNNs. This is because the RIS association scheme strengthens the impact of structural information on system optimization, and DNNs are incapable of capturing the underlying network structure. Unlike the DNNs, which require large amounts of training data, the proposed intelligent GNNs provide better results with fewer samples. Also, their performance can be improved as the number of training samples increases.

## VI. CONCLUSION

In this work, we proposed a heterogeneous GNN structure to solve the NP-hard rate optimization problem in multi-RIS aided mmWave communication systems. The RIS association scheme was shown to significantly improve the system performance while increasing the optimization complexity. By mapping the structure of the wireless environment into a heterogeneous GNN, the optimal beamforming matrices and RIS phase configurations were obtained to maximize the WSR in the considered multi-RIS aided communication system. Extensive simulation results were performed to validate the performance of the proposed GNN. Specifically, it was obtained that the proposed heterogeneous GNN outperformed conventional DNNs by a factor of 10, as the latter is not capable of capturing the underlying structure of the wireless environment. Moreover, simulation results showed that the RIS association scheme can improve the system performance by about 30%. Potential extensions of this study could involve accommodating a variable number of users, assessing deployment costs, and evaluating interference from multiple non-associated RISs. Additionally, exploring the use of soft labels, such as normalized inverse distance distributions instead of traditional binary one-hot coding, may offer further enhancements in performance.

## REFERENCES

[1] M. Liu, C. Huang, M. Di Renzo, M. Debbah, and C. Yuen, "Cooperative beamforming and RISs association for multi-RISs aided multi-users mmwave MIMO systems through graph neural networks," in *Proc. IEEE Int. Conf. Commun. (ICC)*, Rome, Italy, May 2023, pp. 4286-4291.

[2] M. Xiao, S. Mumtaz, Y. Huang, L. Dai, Y. Li, M. Matthaiou, G. K. Karagiannidis, E. Björnson, K. Yang, I. Chih-Lin et al., "Millimeter wave communications for future mobile networks," *IEEE J. Sel. Areas Commun.*, vol. 35, no. 9, pp. 1909-1935, Sep. 2017.

[3] H. Tataria, M. Shafi, A. F. Molisch, M. Dohler, H. Sjöland, and F. Tufvesson, "6G wireless systems: Vision, requirements, challenges, insights, and opportunities," *Proc. IEEE*, vol. 109, no. 7, pp. 1166-1199, Jul. 2021.

[4] C. Huang, A. Zappone, G. C. Alexandropoulos, M. Debbah, and C. Yuen, "Reconfigurable intelligent surfaces for energy efficiency in wireless communication," *IEEE Trans. Wirel. Commun.*, vol. 18, no. 8, pp. 4157-4170, Aug. 2019.

[5] M. Di Renzo, A. Zappone, M. Debbah, M.-S. Alouini, C. Yuen, J. de Rosny, and S. Tretyakov, "Smart radio environments empowered by reconfigurable intelligent surfaces: How it works, state of research, and the road ahead," *IEEE J. Sel. Areas in Commun.*, vol. 38, no. 11, pp. 2450-2525, Nov. 2020.

[6] C. Huang, Z. Yang, G. C. Alexandropoulos, K. Xiong, L. Wei, C. Yuen, Z. Zhang, and M. Debbah, "Multi-hop RIS-empowered terahertz communications: A DRL-based hybrid beamforming design," *IEEE J. Sel. Areas Commun.*, vol. 39, no. 6, pp. 1663-1677, Jun. 2021.

[7] B. Yang, X. Cao, C. Huang, C. Yuen, L. Qian, and M. Di Renzo, "Intelligent spectrum learning for wireless networks with reconfigurable intelligent surfaces," *IEEE Trans. Veh. Technol.*, vol. 70, no. 4, pp. 3920-3925, Apr. 2021.

[8] B. Ning, Z. Chen, W. Chen, and J. Fang, "Beamforming optimization for intelligent reflecting surface assisted MIMO: A sum-path-gain maximization approach," *IEEE Wirel. Commun. Lett.*, vol. 9, no. 7, pp. 1105-1109, Jul. 2020.

[9] Y. Cao, T. Lv, and W. Ni, "Intelligent reflecting surface aided multi-user mmWave communications for coverage enhancement," in *IEEE 31st Annu. Int. Symp. Pers., Indoor Mobile Radio Commun. (PIMRC)*. London, U.K., Aug. 2020, pp. 1-6.

[10] G. Chen, Q. Wu, W. Chen, D. W. K. Ng, and L. Hanzo, "IRS-aided wireless powered MEC systems: TDMA or NOMA for computation offloading?" *IEEE Trans. Wirel. Commun.*, vol. 22, no. 2, pp. 1201-1218, Feb. 2023.

[11] T. Gong, P. Gavrilidis, R. Ji, C. Huang, G. C. Alexandropoulos, L. Wei, Z. Zhang, M. Debbah, H. V. Poor, and C. Yuen, "Holographic MIMO communications: Theoretical foundations, enabling technologies, and future directions," *IEEE Commun. Surveys Tuts.*, 1st Quart. 2024.

[12] J. An, C. Xu, D. W. K. Ng, G. C. Alexandropoulos, C. Huang, C. Yuen, and L. Hanzo, "Stacked intelligent metasurfaces for efficient holographic MIMO communications in 6G," *IEEE J. Sel. Areas Commun.*, Aug. 2023.

[13] E. Shi, J. Zhang, H. Du, B. Ai, C. Yuen, D. Niyato, K. B. Letaief, and X. Shen, "Ris-aided cell-free massive mimo systems for 6g: Fundamentals, system design, and applications," *IEEE Proc.*, vol. 112, no. 4, pp. 331-364, Apr. 2024.

[14] E. Basar, G. C. Alexandropoulos, Y. Liu, Q. Wu, S. Jin, C. Yuen, O. A. Dobre, and R. Schober, "Reconfigurable intelligent surfaces for 6G: Emerging hardware architectures, applications, and open challenges," *IEEE Veh. Technol. Mag.*, vol. 19, no. 3, pp. 2-22, Sep. 2024.

[15] X. Cao, B. Yang, H. Zhang, C. Huang, C. Yuen, and Z. Han, "Reconfigurable-intelligent-surface-assisted MAC for wireless networks: Protocol design, analysis, and optimization," *IEEE Internet Things J.*, vol. 8, no. 18, pp. 14171-14186, Sep. 2021.

[16] T. V. Nguyen, D. N. Nguyen, M. Di Renzo, and R. Zhang, "Leveraging secondary reflections and mitigating interference in multi-IRS/RIS aided wireless networks," *IEEE Trans. Wirel. Commun.*, vol. 22, no. 1, pp. 502-517, Jan. 2023.

[17] A. Abrardo, D. Dardari, and M. Di Renzo, "Intelligent reflecting surfaces: Sum-rate optimization based on statistical position information," *IEEE Trans. Commun.*, vol. 69, no. 10, pp. 7121-7136, Oct. 2021.

[18] A. Abrardo, D. Dardari, M. Di Renzo, and X. Qian, "MIMO interference channels assisted by reconfigurable intelligent surfaces: Mutual coupling aware sum-rate optimization based on a mutual impedance channel model," *IEEE Wirel. Commun. Lett.*, vol. 10, no. 12, pp. 2624-2628, Oct. 2021.
[19] N. S. Perović, L.-N. Tran, M. Di Renzo, and M. F. Flanagan, "On the maximum achievable sum-rate of the RIS-aided MIMO broadcast channel," *IEEE Trans. Signal Process.*, vol. 70, pp. 6316–6331, Jan. 2023.

[20] X. Wei, D. Wu, L. Zhou, and M. Guizani, "Cross-modal communication technology: A survey," *Fundam Res.*, Sep. 2023. [Online]. Available: https://www.sciencedirect.com/science/article/pii/S2667325823002297

[21] Y. Huang, "Challenges and opportunities of sub-6 GHz integrated sensing and communications for 5G-advanced and beyond," *Chin. J. Electron.*, vol. 33, no. 2, pp. 323–325, Mar. 2024.

[22] J. Gao, C. Zhong, X. Chen, H. Lin, and Z. Zhang, "Unsupervised learning for passive beamforming," *IEEE Commun. Lett.*, vol. 24, no. 5, pp. 1052–1056, May 2020.

[23] W. Xu, L. Gan, and C. Huang, "A robust deep learning-based beam-forming design for RIS-assisted multiuser MISO communications with practical constraints," *IEEE Trans. Cognit. Commun. Networking*, Jun. 2021.

[24] G. C. Alexandropoulos, K. Stylianopoulos, C. Huang, C. Yuen, M. Ben-nis, and M. Debbah, "Pervasive machine learning for smart radio environments enabled by reconfigurable intelligent surfaces," *Proc. IEEE*, vol. 110, no. 9, pp. 1494–1525, Sep. 2022.

[25] X. Cao, B. Yang, C. Huang, C. Yuen, M. Di Renzo, Z. Han, D. Niyato, H. V. Poor, and L. Hanzo, "AI-assisted MAC for reconfigurable intelligent-surface-aided wireless networks: challenges and opportunities," *IEEE Commun. Mag.*, vol. 59, no. 6, pp. 21–27, Jun. 2021.

[26] Y. Bengio, A. Lodi, and A. Prouvost, "Machine learning for combinatorial optimization: a methodological tour d'horizon," *Eur. J. Oper. Res.*, vol. 290, no. 2, pp. 405–421, Apr. 2021.

[27] S. He, S. Xiong, Y. Ou, J. Zhang, J. Wang, Y. Huang, and Y. Zhang, "An overview on the application of graph neural networks in wireless networks," *IEEE Open J. Commun. Soc.*, vol. 2, pp. 2547–2565, Nov. 2021.

[28] Y. Shen, J. Zhang, S. H. Song, and K. B. Letaief, "Graph neural networks for wireless communications: From theory to practice," *IEEE Trans. WIREL. COMM.*, vol. 22, no. 5, pp. 3554–3569, May 2023.

[29] S. He, S. Xiong, W. Zhang, Y. Yang, J. Ren, and Y. Huang, "GBLinks: GNN-based beam selection and link activation for ultra-dense D2D mmWave networks," *IEEE Trans. Commun.*, vol. 70, no. 5, pp. 3451–3466, May 2022.

[30] M. Schlichtkrull, T. N. Kipf, P. Bloem, R. van den Berg, I. Titov, and M. Welling, "Modeling relational data with graph convolutional networks," in *The Semantic Web*, A. Gangemi, R. Navigli, M.-E. Vidal, P. Hitzler, R. Troncy, L. Hollink, A. Tordai, and M. Alam, Eds. Cham: Springer International Publishing, Jun. 2018, pp. 593–607.

[31] Y. Shen, Y. Shi, J. Zhang, and K. B. Letaief, "Graph neural networks for scalable radio resource management: Architecture design and theoretical analysis," *IEEE J. Sel. Areas in Commun.*, vol. 39, no. 1, pp. 101–115, Jan. 2021.

[32] X. Zhang, H. Zhao, J. Xiong, X. Liu, L. Zhou, and J. Wei, "Scalable power control/beamforming in heterogeneous wireless networks with graph neural networks," in *Proc. IEEE Global Commun. Conf. (GLOBECOM)*, Madrid, Spain, Dec. 2021, pp. 01–06.

[33] T. Jiang, H. V. Cheng, and W. Yu, "Learning to reflect and to beamform for intelligent reflecting surface with implicit channel estimation," *IEEE J. Sel. Areas Commun.*, vol. 39, no. 7, pp. 1931–1945, Jul. 2021.

[34] Z. Zhang, T. Jiang, and W. Yu, "Learning based user scheduling in reconfigurable intelligent surface assisted multiuser downlink," *IEEE J. Sel. Top. Sign. Proces.*, Aug. 2022.

[35] K.-L. Chan, F.-T. Chien, and R. Y. Chang, "Beamforming and load-balanced user association in RIS-aided systems via graph neural networks," in *IEEE Int. Conf. Commun.*, Denver, Colorado, USA, Jun. 2024, pp. 4997–5002.

[36] T.-J. Yeh, W.-C. Tsai, C.-W. Chen, and A.-Y. Wu, "Enhanced-GNN with angular CSI for beamforming design in IRS-assisted mmwave communication systems," *IEEE Commun. Lett.*, Apr. 2024.

[37] Y. Zhu, M. Li, Y. Liu, Q. Liu, Z. Chang, and Y. Hu, "DRL-based joint beamforming and BS-RIS-UE association design for RIS-assisted mmWave networks," in *IEEE Wireless Commun. and Netw. Conf (WCNC)*, Austin, TX, USA, Apr. 2022, pp. 345–350.

[38] W. He, D. He, X. Ma, X. Chen, Y. Fang, and W. Zhang, "Joint user association, resource allocation, and beamforming in RIS-assisted multi-server MEC systems," *IEEE Trans. WIREL. COMM.*, vol. 23, no. 4, pp. 2917–2932, Apr. 2024.

[39] Y. Huang, W. Mei, and R. Zhang, "Empowering base stations with co-site intelligent reflecting surfaces: User association, channel estimation and reflection optimization," *IEEE Trans. Commun.*, vol. 70, no. 7, pp. 4940–4955, Jul. 2022.

[40] C. Zhang, D. Song, C. Huang, A. Swami, and N. V. Chawla, "Heterogeneous graph neural network," in *Proc. 25th ACM SIGKDD Int. Conf. Knowl. Discov. Data Mining*, Jul. 2019, pp. 793–803.

[41] X. Ma, Z. Gao, F. Gao, and M. Di Renzo, "Model-driven deep learning based channel estimation and feedback for millimeter-wave massive hybrid MIMO systems," *IEEE J. Sel. Areas in Commun.*, vol. 39, no. 8, pp. 2388–2406, Aug. 2021.

[42] J. He, H. Wymeersch, M. Di Renzo, and M. Juntti, "Learning to estimate RIS-aided mmwave channels," *IEEE WIREL. COMMUN. Lett.*, vol. 11, no. 4, pp. 841–845, Apr. 2022.

[43] M. Liu, X. Li, B. Ning, C. Huang, S. Sun, and C. Yuen, "Deep learning-based channel estimation for double-RIS aided massive MIMO system," *IEEE WIREL. COMMUN. Lett.*, vol. 12, no. 1, pp. 70–74, Jan. 2023.

[44] L. Wei, C. Huang, G. C. Alexandropoulos, C. Yuen, Z. Zhang, and M. Debbah, "Channel estimation for RIS-empowered multi-user MISO wireless communications," *IEEE Trans. Commun.*, vol. 69, no. 6, pp. 4144–4157, Jun. 2021.

[45] H. Guo and V. K. Lau, "Uplink cascaded channel estimation for intelligent reflecting surface assisted multiuser MISO systems," *IEEE Trans. Signal Process.*, vol. 70, pp. 3964–3977, Jul. 2022.

[46] C. Pan, G. Zhou, K. Zhi, S. Hong, T. Wu, Y. Pan, H. Ren, M. Di Renzo, A. Lee Swindlehurst, R. Zhang, and A. Y. Zhang, "An overview of signal processing techniques for RIS/IRS-aided wireless systems," *IEEE J. Sel. Top. Sign. Proces.*, vol. 16, no. 5, pp. 883–917, Aug. 2022.

[47] K. Zhi, C. Pan, H. Ren, K. Wang, M. Elkashlan, M. Di Renzo, R. Schober, H. V. Poor, J. Wang, and L. Hanzo, "Two-timescale design for reconfigurable intelligent surface-aided massive MIMO systems with imperfect CSI," *IEEE Trans. Inf. Theory*, vol. 69, no. 5, pp. 3001–3033, May 2023.

[48] M. R. Akdeniz, Y. Liu, M. K. Samimi, S. Sun, S. Rangan, T. S. Rappaport, and E. Erkip, "Millimeter wave channel modeling and cellular capacity evaluation," *IEEE J. Sel. Areas in Commun.*, vol. 32, no. 6, pp. 1164–1179, Jun. 2014.

[49] L. Wei, C. Huang, Q. Guo, Z. Yang, Z. Zhang, G. C. Alexandropoulos, M. Debbah, and C. Yuen, "Joint channel estimation and signal recovery for RIS-empowered multiuser communications," *IEEE Trans. Commun.*, vol. 70, no. 7, pp. 4640–4655, Jul. 2022.

[50] A. Alkhateeb, O. El Ayach, G. Leus, and R. W. Heath, "Channel estimation and hybrid precoding for millimeter wave cellular systems," *IEEE J. Sel. Top. Signal Proces.*, vol. 8, no. 5, pp. 831–846, Oct. 2014.

[51] P. Wang, J. Fang, X. Yuan, Z. Chen, and H. Li, "Intelligent reflecting surface-assisted millimeter wave communications: Joint active and passive precoding design," *IEEE Trans. Veh. Technol.*, vol. 69, no. 12, pp. 14960–14973, Dec. 2020.

[52] A. L. Maas, A. Y. Hannun, A. Y. Ng et al., "Rectifier nonlinearities improve neural network acoustic models," in *Proc. Int. Conf. Mach. Learn. Workshop Deep Learn. Audio, Speech Lang. Process. (ICML WDLASL)*. Atlanta, Georgia, USA, Jun. 2013, pp. 1–6.

[53] P. W. Battaglia, J. B. Hamrick, V. Bapst, A. Sanchez-Gonzalez, V. Zam-baldi, M. Malinowski, A. Tacchetti, D. Raposo, A. Santoro, R. Faulkner et al., "Relational inductive biases, deep learning, and graph networks," *arXiv preprint arXiv:1806.01261*, Jun. 2018.

[54] M. Fey and J. E. Lenssen, "Fast graph representation learning with pytorch geometric," *arXiv preprint arXiv:1903.02428*, Apr. 2019.

[55] K. He, X. Zhang, S. Ren, and J. Sun, "Deep residual learning for image recognition," in *Proc. IEEE Conf. Comput. Vis. Pattern Recognit. (CVPR)*, Jun. 2016, pp. 770–778.

[56] D. P. Kingma and J. Ba, "Adam: A method for stochastic optimization," in *Proc. 3rd Int. Conf. Learn. Represent. (ICLR)*, San Diego, CA, USA, May 2015.

[57] Q. Shi, M. Razaviyayn, Z.-Q. Luo, and C. He, "An iteratively weighted MMSE approach to distributed sum-utility maximization for a MIMO interfering broadcast channel," *IEEE Trans. Signal Process.*, vol. 59, no. 9, pp. 4331–4340, May 2011.

[58] H. Guo, Y.-C. Liang, J. Chen, and E. G. Larsson, "Weighted sum-rate maximization for reconfigurable intelligent surface aided wireless networks," *IEEE Trans. WIREL. Commun.*, vol. 19, no. 5, pp. 3064–3076, May 2020.

[59] N. Boumal, B. Mishra, P.-A. Absil, and R. Sepulchre, "Manopt, a matlab toolbox for optimization on manifolds," *J. Mach. Learn. Res.*, vol. 15, no. 1, pp. 1455–1459, Jan. 2014.

[60] E. Basar and H. V. Poor, "Present and future of reconfigurable intelligent surface-empowered communications [perspectives]," *IEEE Signal Process. Mag.*, vol. 38, no. 6, pp. 146–152, Nov. 2021.