Multi-User MISO with Stacked Intelligent Metasurfaces: A DRL-Based Sum-Rate Optimization Approach

As shown in Fig. 1, we consider a downlink MISO system comprised of a BS, a SIM, and multiple UEs, in which the SIM is integrated with the BS to improve system performance. The SIM consists of multiple programmable metasurfaces, each of which is composed of a large number of meta-atoms that are linked to an intelligent controller, e.g., an FPGA [17]. By appropriately tuning the transmission coefficients of the SIM via this controller, the SIM is capable of directly manipulating the electromagnetic (EM) behaviors of the propagating waves [34]. As a result, SIM can tailor a spatial waveform shape at the output metasurface layer. Note that traditional communication systems without SIM require many active RF chains to generate the desired EM waveforms. In contrast, employing SIM can significantly reduce the number of active RF chains [18] and enable low-precision DACs¹¹1In practical communication systems, low-precision DACs can be utilized by employing modulation schemes with lower DAC requirements, such as binary phase-shift keying (BPSK). This approach accordingly reduces overall system power consumption and cost.. This is achieved by directly steering the EM waves as they pass through the low-cost and energy-efficient metasurfaces.

Let $M$ represent the number of UEs with the corresponding set given by $\mathcal{M}=\{1,2,...,M\}$. At the BS, $M$ antennas are selected to concurrently transmit $M$ individual data streams associated with $M$ UEs²²2In this paper, we assume that $M$ antennas have been selected. The task of simultaneously optimizing antenna selection and wave-based precoding constitutes future research work.. Let $L$ denote the number of equidistance-spaced metasurface layers of the SIM and $N$ denote the number of meta-atoms on each metasurface layer, satisfying $N\geq M$, while the metasurface layer set is represented by $\mathcal{L}=\{1,2,...,L\}$ and the meta-atom sets of each layer are represented by $\mathcal{N}=\{1,2,...,N\}$, respectively. Besides, let $\phi_{n}^{l}=e^{j\varphi_{n}^{l}}$ denote the transmission coefficient of the $n$-th meta-atom on the $l$-th metasurface layer. $\varphi_{n}^{l}$ denotes the corresponding phase shift, which is assumed to be adjustable between $0$ and $2\pi$ continuously [18, 34, 7], i.e., $\varphi_{n}^{l}\in[0,2\pi)$, $n\in\mathcal{N}$, $l\in\mathcal{L}$. Then, the transmission coefficient vector and the corresponding matrix representation of the $l$-th metasurface are denoted by $\boldsymbol{\phi}^{l}=[\phi_{1}^{l},\phi_{2}^{l},...,\phi_{N}^{l}]^{T}\in\mathbb{C}^{N\times 1}$ and $\boldsymbol{\Phi}^{l}=\text{diag}(\boldsymbol{\phi}^{l})\in\mathbb{C}^{N\times N}$, respectively.

Without loss of generality, each metasurface layer of the SIM is modeled as a uniform planar array, with all the metasurface layers being arranged in an identical square configuration. Let $r_{e}$ denote the element spacing between two adjacent meta-atoms within the same layer. Thus, the spacing between the $n$-th meta-atom and the $\tilde{n}$-th meta-atom on the same layer is

where $n_{x}$ and $n_{y}$ are the $x$-axis and $y$-axis indices of the $n$-th meta-atom, respectively, with $\tilde{n}_{x}$ and $\tilde{n}_{y}$ defined similarly for the $\tilde{n}$-th meta-atom. The indices are calculated as $n_{x}=\text{mod}(n-1,\,n_{\mathrm{max}})+1$ and $n_{y}=\lceil n/n_{\mathrm{max}}\rceil$, where $n_{\mathrm{max}}=\sqrt{N}$ represents the maximum number of meta-atoms per row.

Furthermore, the total thickness of the SIM structure is denoted as $D$, where the parallel metasurface layers are equally spaced with a distance of $d_{s}=D/(L-1)$ between adjacent layers. Thus, the spacing between the $n$-th meta-atom on the $l$-th layer and the $\tilde{n}$-th meta-atom on the ($l-1$)-th layer is $r_{n,\tilde{n}}^{l}=\sqrt{r_{n,\tilde{n}}^{2}+d_{s}^{2}}$. Moreover, we assume that the transmit antennas are arranged in a uniform linear array (ULA) with a half-wavelength distance between them. The array’s center aligns with that of the first metasurface. Thus, the distance between the $m$-th transmitting antenna and the $n$-th meta-atom on the first metasurface is defined as

where $\lambda$ denotes the wavelength of the carrier signal.

Furthermore, the propagation coefficient of the EM wave between adjacent metasurface layers can be obtained from the Rayleigh-Sommerfeld diffraction equation [34, 18]. According to these equations, the propagation coefficient from the $\tilde{n}$-th meta-atom on the $(l-1)$-th metasurface layer to the $n$-th meta-atom on the $l$-th metasurface layer is defined by

$\mathrm{for}\,l\in\mathcal{L}/\{1\},\,n,\tilde{n}\in\mathcal{N}$, where $s_{a}$ represents the area of each meta-atom. The propagation coefficient $[\mathbf{W}^{1}]_{n,m}$ from the $m$-th transmitting antenna to the $n$-th meta-atom on the first metasurface is obtained by replacing $r_{n,\tilde{n}}^{l}$ in (2) with $r^{1}_{m,n}$ in (II-A). Furthermore, we assume that non-adjacent layers do not interact with each other [43, 44, 46, 47, 48, 50, 51]. Thus, the overall EM response in the SIM can be expressed as

For the wireless channels shown in Fig. 1, we consider a spatially correlated channel model due to the closely spaced meta-atoms [43, 44, 46, 47, 48, 50, 51]. Thus, the channel spanning from the output metasurface to the $M$ UEs can be expressed as

where $\mathbf{R}\in\mathbb{C}^{N\times N}$ denotes the spatial correlation matrix of the SIM, and $\tilde{\mathbf{G}}\in\mathbb{C}^{M\times N}$ represents the independent and identically distributed (i.i.d.) Rayleigh fading channel satisfying $[\tilde{\mathbf{G}}]_{m,:}\sim\mathcal{CN}(\mathbf{0},\rho_{m}^{2}\mathbf{I})$, where $\rho_{m}^{2}$ denotes the path loss between the output metasurface of the SIM and the $m$-th UE. In particular, the path loss of the $m$-th UE is given by

where $C_{0}\in\mathbb{R}$ denotes the path loss at a reference distance of $1$ meter, $d_{c,m}$ denotes the transmission distance between the SIM and the $m$-th UE, and $\alpha$ denotes the corresponding path loss exponent. Besides, assuming far-field wave propagation in an isotropic scattering wireless environment, the spatial correlation matrix $\mathbf{R}$ is thus defined by [50, 4]

We focus on the downlink data transmission and the signal received at the $m$-th UE is given as

where $z_{m}$ denotes the additive white Gaussian noise (AWGN) with variance $\sigma^{2}_{m}$, i.e., $z_{m}\sim\mathcal{CN}(0,\sigma^{2}_{m})$. $\mathbf{P}=\text{diag}(\sqrt{p}_{1},\sqrt{p}_{2},...,\sqrt{p}_{M})\in\mathbb{C}^{M\times M}$ denotes the transmission power matrix with $p_{m},m\in\mathcal{M}$ representing the transmit power allocated to the $m$-th UE. $\mathbf{x}=[x_{1},x_{2},...,x_{M}]\in\mathbb{C}^{M\times 1}$ denotes a column vector of the data streams transmitted to all the UEs, satisfying $\mathbb{E}[\lVert\mathbf{x}\rVert^{2}]=1$. We also assume that our passive SIM does not introduce any noise to its internally propagating signals [17].

Furthermore, the received signal of (7) can be rewritten as

The second addend of (8) is considered as the co-channel interference caused by other $(M-1)$ users. Thus, the received signal-to-interference-plus-noise ratio (SINR) at the $m$-th UE is defined as

In this paper, we adopt the sum rate of $M$ UEs as a system performance metric, which is defined by

Our objective is to design the SIM phase shifts $\mathbf{\Phi}^{l}$ and the transmit power allocation strategy $\mathbf{P}$ that jointly maximize the system sum-rate $C(\mathbf{P},\boldsymbol{\Phi}^{l})$ under perfect CSI of all involved wireless channels³³3The examination of the effects of imperfect CSI is left for future work.. As such, the joint design can be formulated as the following joint optimization problem:

where $P$ denotes the maximum transmit power budget at the BS. It can be verified that problem (P1) constitutes non-convex optimization due to non-convexities arising from the objective function, constraints (11c), and consequently, in the multilayered architecture of the SIM (11b). It has been shown that traditional techniques to address such problems, such as the AO methods shown, e.g., in [52, 53, 18], often struggle to find satisfactory solutions within moderate complexity. For this reason, in this paper, we propose a novel DRL approach to determine an effective $\mathbf{\Phi}^{l},\,l\in\mathcal{L}$, and power allocation solution $\mathbf{P}$ with low complexity. It is noted that given the availability of CSI, an efficient DRL approach will empower the SIM to persistently interact with the environment, enabling autonomous learning and self-guided exploration through feedback from the environment. This methodology enables DRL to derive optimal configuration without the need for manual data gathering for training, which is in contrast with the conventional deep learning algorithms that necessitate both offline training and online learning [54, 13].

Model-free DRL is a dynamic tool capable of solving a series of decision-making problems, enabling an agent to learn the best strategy in a real-time manner [54, 41, 14]. A general DRL system consists of two components: the agent and the environment. The agent continuously interacts with the environment, applying actions based on its current policy and observing the immediate rewards and new states from the environment. The fundamental concepts of the proposed DRL formulation for handling (P1) are as follows:

• •

  Action space: Let $\mathcal{A}$ denote the action space. According to the state of the environment $\mathbf{s}_{t}$, the agent takes an action $\mathbf{a}_{t}\in\mathcal{A}$ to influence the environment according to a policy. Then, the agent receives the reward $r_{t}$ and observes the new state $\mathbf{s}_{t+1}$ feedback. For our system model, the action at each time step $t$ is constructed by the current SIM phase shifts $\boldsymbol{\phi}_{t}^{l},\,l\in\mathcal{L}$ and transmit power allocation $p_{m,t},\,m\in\mathcal{M}$. Since DNNs require real-valued inputs and outputs, the SIM phase shifts are rearranged into real and imaginary components, i.e.,

  $$\mathbf{a}_{t}=\big{[}\Re(\boldsymbol{\phi}^{1}_{t}),...,\Im(\boldsymbol{\phi}^{1}_{t}),...,\Im(\boldsymbol{\phi}^{L}_{t}),p_{1,t},...,p_{M,t}\big{]}^{T},$$

  (12)

  where $p_{m,t}$ is the transmit power of the $m$-th antenna at time step $t$. Hence, the action is represented as a vector with dimension $D_{a}=2NL+M$.

• •

  State space: Let $\mathcal{S}$ denote the environment state space. The current time-slot environment state $\mathbf{s}_{t}\in\mathcal{S}$ is observed from the environment. The current state incorporates the reward $r_{t-1}$ and action $\mathbf{a}_{t-1}$ from the previous time step along with the CSI for all UEs, i.e.,

  $\displaystyle\mathbf{s}_{t}=\big{[}r_{t-1},\mathbf{a}_{t-1},\Re([\mathbf{G}]_{1,:}),...,\Im([\mathbf{G}]_{1,:}),...,\Im([\mathbf{G}]_{M,:})\big{]}^{T}.$

  (13)

  Hence, the dimension of the state space is $D_{s}=2N(L+M)+M+1$.

• •

  Policy: Let $\pi$ denote the policy of agent and $\pi(\mathbf{s}_{t},\mathbf{a}_{t})$ denote the probability of taking $\mathbf{a}_{t}$ based on the environment state $\mathbf{s}_{t}$, satisfying $\sum_{\mathbf{a}_{t}\in\mathcal{A}}\pi(\mathbf{s}_{t},\mathbf{a}_{t})=1$.

• •

  Reward: The reward $r_{t}$ serves as a measure to quantify the quality of the policy. The agent takes action $\mathbf{a}_{t}$ in the environment and then dynamically adjusts its policy based on the feedback reward, aiming to enhance its action performance and maximize the reward. The objective of this paper is to maximize the sum rate of $M$ UEs. As such, the reward at each time step $t$ is defined as $C(\mathbf{P},\boldsymbol{\Phi}^{l})$, which is determined based on the current power allocation scheme $\mathbf{P}$ and phase shifts of SIM $\{\boldsymbol{\Phi}^{l}\}$, both are defined by the output of the actor training neural network.

• •

  Transitions: Action $\mathbf{a}_{t}$ is sampled from $\pi(\mathbf{s}_{t},\mathbf{a}_{t})$ which determines the phase shifts of the SIM and the power allocation. A transmission then occurs within a coherence block where $\mathbf{G}$ remains constant and $C(\mathbf{P},\boldsymbol{\Phi}^{l})$ is computed. Then, the environment/system proceeds to time step $t+1$ where new channel matrices are sampled/observed in i.i.d. fashion.

Generally, there are two categories of algorithms to determine the policy in DRL, i.e., value-based algorithms and policy-based algorithms [15, 40]. Specifically, the well-known deep Q-network (DQN) is a quintessential value-based algorithm that employs experience replay and target network strategies to enhance the stability and efficiency of its learning process. However, DQN is typically applicable to a discrete action space. Yet, in the SIM-aided scenario, a continuous action space spanned by SIM phase shifts and the transmit power allocation strategy is presented in (11). On the other hand, the policy gradient (PG) method is a basic policy-based algorithm that aims to maximize the long-term expected discounted reward of each episode for continuous action spaces [40]. In fact, the PG directly adjusts the likelihood of selecting particular actions based on the gradient of expected reward. However, the pure PG approach exhibits high variance, leading to unstable learning. Moreover, it demands an extensive sample set. Thus, the convergence performance and stability of PG in the wireless communication environment need to be improved [54]. To address this, we adopt a DDPG-based solution, which can handle continuous action space along with excellent stability and convergence, to solve the problem in (P1).

DDPG also exploits a Q-function to choose an action $\mathbf{a}_{t}$ under the policy $\pi$ according to the current state $\mathbf{s}_{t}$. Specifically, the goal of DDPG is to find the best policy $\pi_{\mathrm{op}}$ that maximizes the long-term expected discounted reward. The Q-function evaluates state-action pairs via a Q-value, and the DDPG employs a separate critic as the Q-function. As shown in Fig. 2, the DDPG algorithm utilizes two neural networks with distinct roles: actor networks (depicted in orange) and critic networks (in green). The DDPG exploits the actor neural network $\pi(\mathbf{s}_{t};\theta_{\pi})$, with $\theta_{\pi}$ representing its parameters, as the policy to generate the SIM phase shifts and transmit power allocation strategy from a continuous action space $\mathcal{A}$. Besides, the critic neural network $Q(\mathbf{s}_{t},\mathbf{a}_{t};\theta_{q})$, having the parameters $\theta_{q}$, is utilized to output a Q-value, which measures the output of the actor network.

During the training phase, DRL operates without target data, which actually differentiates it from classic supervised training [40, 13, 11]. To address this, the proposed DDPG adopts two networks with identical structures but with different parameters as shown in Fig. 2. The critic target network $\tilde{Q}(\mathbf{s}_{t+1},\tilde{\mathbf{a}}_{t};\theta_{\tilde{q}})$ is utilized to update the critic training network $Q(\mathbf{s}_{t},\mathbf{a}_{t};\theta_{q})$, while the actor target network $\tilde{\pi}(\mathbf{s}_{t+1};\theta_{\tilde{\pi}})$ updates the actor training network $\pi(\mathbf{s}_{t};\theta_{\pi})$, where $\tilde{\mathbf{a}}_{t}$ denotes the output of actor target network $\tilde{\pi}(\mathbf{s}_{t+1};\theta_{\tilde{\pi}})$. The training neural networks and the target neural networks are built with an identical structure but different parameters. While the training networks undergo training as classic DNNs, the target networks are untrainable and aim to provide a label for the training network. The parameters of the target neural network are updated periodically. At the end of each step, the updates on the critic and actor target neural networks are expressed as

where $\eta_{c}\in(0,1)$ and $\eta_{a}\in(0,1)$ denote the learning rate for updating the critic and actor target neural network, respectively. The updates on the target neural networks can stabilize and accelerate the convergence of the DRL training process [40].

Moreover, the critic training neural network’s update employs the gradient of a loss function. The label in the loss function consists of the reward $r_{t}$, which is returned by the SIM-aided multi-user MISO wireless environment, and the output of the critic target neural network $\tilde{Q}(\mathbf{s}_{t+1},\tilde{\mathbf{a}}_{t};\theta_{\tilde{q}})$, which can be expressed as

where $\mu$ denotes a discount factor. In this paper, we employ the mean square error (MSE) loss function, which is defined as

where the error between $V(\mathbf{s}_{t},\mathbf{a}_{t})$ and estimated value $Q(\mathbf{s}_{t},\mathbf{a}_{t};\theta_{q})$ is termed the temporal difference error [40]. Thus, the parameter update of the critic training neural network can be expressed as

where $\gamma_{c}>0$ denotes the corresponding learning rate.

The update on the actor training neural network is defined through the policy gradient theorem [40], as follows

where $\gamma_{a}>0$ denotes its learning rate. $\Delta_{\pi}Q(\mathbf{s}_{t},\pi(\mathbf{s}_{t};\theta_{\pi});\theta_{q})$ is the gradient of the critic training network regarding the output of the actor training neural network. $\Delta_{\theta_{\pi}}\pi(\mathbf{s}_{t};\theta_{\pi})$ represents the gradient of the actor training neural network with respect to its parameters $\theta_{\pi}$. Notably, the update process of the actor training neural network depends heavily on the gradient of the critic training neural network according to the current policy. This interaction ensures that the action policy is updated towards a direction that maximizes the long-term expected sum rate, as estimated through the Q-values.

It is noted, however, that value-based algorithms can be trapped in local optima due to the correlation between samples and nonstationary targets [40, 54]. In addition, when considering short timescales, this behavior is particularly significant, since data generated through interactions with wireless environments tend to exhibit high temporal correlation. Therefore, the experience replay technique is utilized to reduce the negative impact of sample correlation on the agent by using a buffer window to store a portion of data. For each training phase, networks adopt a mini-batch randomly sampled from experience replay to calculate the gradient and update parameters. Moreover, in our DRL framework, we introduce noise to the action to prevent the learning process from being trapped in a locally optimal solution [40]. This whitening process is applied to the SIM phase shifts and transmit power allocation prior to utilization within the wireless network. Specifically, truncated random white noise $\mathcal{W}\in[-w_{a},w_{a}]$, where $w_{a}$ defines the truncation value, is integrated into the action processing, which is expressed as

The noise $\mathcal{W}$ follows a Gaussian distribution with zero mean and variance $v$, i.e., $\widetilde{\mathcal{W}}\sim\mathcal{CN}(0,v)$. To enable smooth convergence, $v$ decreases exponentially during training as

where $v_{0}$ denotes the initial $v$ value while $\zeta\in(0,1)$ and $t_{\mathrm{gap}}$ are the decay rate and the gap factor for discounting the whitening process, respectively.

Fig. 3 shows the structures of the proposed networks. The actor neural network operates by inputting a state vector, directly producing the SIM phase shifts and transmit power allocation strategy. This network consists of two branches, each addressing the optimization of SIM phase shifts and power allocation, respectively. The branch assigned to process SIM phase shifts consists of a reshape operation, two blocks, an adaptive average pooling layer, a dense layer, a flattening operation, and a normalization layer. The reshape operation serves to transform the input vector $\boldsymbol{\phi}^{l}\in\mathbb{C}^{N\times 1},\,l\in\mathcal{L}$, into the corresponding matrix in the shape of $n_{\mathrm{max}}\times n_{\mathrm{max}}$, while the flatten operation is employed to streamline the subsequent dense layer processing. Each block is formed by the interweaving connection of two convolutional layers of a $3\times 3$ kernel size, two layer-normalization (LN) layers, and a residual connection with a convolutional layer of a $1\times 1$ kernel size. The output channels of each convolutional layer are consistent. The adaptive average pooling layer restricts the output feature map size of the block to $3\times 3$. The power allocation branch is composed of three dense layers. Finally, the outputs from the two branches are concatenated to yield an action. The critic network starts with two input branches, each comprising a dense layer. Their outputs are summed as input to alternating dense and LN layers.

The dense layer neuron quantity exceeds input and output dimensions, based on the SIM meta-atom count and data streams. The normalization layer is exploited to restrict the output of the actor neural network such that the module of meta-atoms of the metasurface layers is 1, i.e., $|\phi^{l}_{n}|^{2}=1,l\in\mathcal{L},n\in\mathcal{N}$. The LN layer is designed to overcome the influence of the data distribution change caused by the variations in the previous layer on the neural network, which is capable of accelerating the learning speed and improving the robustness of the neural network. In contrast to the batch normalization layer that normalizes data across the batch, the LN layer performs normalization individually for each sample, maintaining independence among samples. This property endows it with superior performance in handling sequential data.

On the other hand, the LeakyReLu activation function is utilized after the hidden dense layer and convolutional layer, while the Tanh function is employed before the normalization layer in the branch processing the SIM phase shifts. The Softmax function is inserted at the end of the power allocation branch for power normalization. The Adam optimizer [13] is adopted to update the parameters of both the actor and critic training neural networks.

Regarding the learning rate, if the model exhibits no improvement even after $\iota_{p}$ iteration steps and the reward ceases to increase, we multiply the learning rate by a decay factor of $\iota_{f}$ to reduce the learning rate. This approach enhances the model’s convergence stability toward a near-optimal solution.

In our DRL approach, an agent is assigned to continuously collect channel coefficients $\mathbf{G}$ and combine action $\mathbf{a}_{t-1}$ in experience replay memory to form current state information $\mathbf{s}_{t}$. At the very beginning of the algorithm, we need to create four neural networks, i.e., the actor training neural network $\pi(\mathbf{s}_{t};\theta_{\pi})$, the actor target neural network $\tilde{\pi}(\mathbf{s}_{t};\theta_{\tilde{\pi}})$, critic training neural network $Q(\mathbf{s}_{t},\mathbf{a}_{t};\theta_{q})$, and critic target neural network $\tilde{Q}(\mathbf{s}_{t},\tilde{\mathbf{a}}_{t};\theta_{\tilde{q}})$, and then initialize the parameters of $\pi(\mathbf{s}_{t};\theta_{\pi})$ and $Q(\mathbf{s}_{t},\mathbf{a}_{t};\theta_{q})$ i.e., $\theta_{q}$ and $\theta_{\pi}$, by randomly sampling in a continuous space. In contrast, the parameters of $\tilde{\pi}(\mathbf{s}_{t};\theta_{\tilde{\pi}})$ and $\tilde{Q}(\mathbf{s}_{t},\tilde{\mathbf{a}}_{t};\theta_{\tilde{q}})$ are given by (14) and (15), respectively. Moreover, the experience replay memory $M_{\mathrm{er}}$ with the capacity $C_{\mathrm{er}}$ is also established.

The algorithm executes for $E$ episodes and iterates $T$ steps in each episode. At episode start, the agent collects CSI, resets the experience buffer and whitening noise, randomly initializes SIM phase shifts from $0$ to $2\pi$ for $\mathbf{B}$, and equally allocates transmit antenna power. During a specific episode, the agent first obtains the initial state $\mathbf{s}_{0}$. Then, taking the state $\mathbf{s}_{t}$ into the actor training network to output corresponding action $\mathbf{a}_{t}$. The agent obtains the reward $r_{t}$ by (10) and the next state $\mathbf{s}_{t+1}$ from the wireless environment. Subsequently, storing a tuple ($\mathbf{s}_{t}$, $\mathbf{a}_{t}$, $r_{t}$, $\mathbf{s}_{t+1}$) into experience replay memory $M_{\mathrm{er}}$. When the number of stored tuples attain the capacity $C_{\mathrm{er}}$ of the $M_{\mathrm{er}}$, signifying a full replay buffer, the critic and actor training networks start to randomly sample mini-batches of size $N_{B}$ from $M_{\mathrm{er}}$ and update their parameters utilizing (18) and (19), respectively. Ultimately, the critic target network and actor target network are updated through a soft update method, as outlined in (14) and (15), respectively.

Finally, the optimal SIM phase shifts $\{\boldsymbol{\phi}^{1}_{\mathrm{op}},\boldsymbol{\phi}^{2}_{\mathrm{op}},...,\boldsymbol{\phi}^{L}_{\mathrm{op}}\}$ and the optimal transmit power allocation strategy $\mathbf{P}_{\mathrm{op}}$ are directly derived from the action corresponding to the maximized sum rate $C_{\mathrm{op}}(\mathbf{P},\boldsymbol{\Phi}^{l})$, which is associated the largest instantaneous reward in the current episode. The proposed design is detailed in Algorithm 1.

The computational complexity analysis of our proposed DRL algorithm can be divided into the training and predicting phases. The former involves experience replay, whitening, and four neural networks, while the latter only relates to the actor network. Denote the complexity of the adaptive pooling layer, the LN layer, and the activation function as $v_{\mathrm{pool}}$, $v_{\mathrm{ln}}$, and $v_{\mathrm{activation}}$, respectively.

1) Training phase: The convolutional layer complexity is encompassed within the actor network calculations. As detailed in [55], the complexity considering the output feature map size can be expressed as $c_{l-1}\times s_{c}^{2}\times c_{l}\times s_{f}$, where $c_{l-1}$ denotes the number of channels in the previous layer, $s_{c}$ is the convolutional kernel size, $c_{l}$ refers to the output channel number, and $s_{f}$ indicates the output feature map size. Thus, the complexity of the two blocks can be expressed as $v_{\mathrm{blocks}}=2\times(L\times 3^{2}\times c\times N)+2\times(c\times 3^{2}\times c\times N)+2\times(c\times 1^{2}\times 2L\times N)+6v_{\mathrm{activation}}+4v_{\mathrm{ln}}\approx\mathcal{O}(NLc+Nc^{2}+v_{\mathrm{ln}})$, where $c$ denotes the adopted input/output channel of the convolutional layer. The adaptive pooling layer maintains the number of channels in the output feature map of the block, while decreasing the dimensions to $3\times 3$ to reduce the computational complexity. The complexity of the normalization layer is proportional to $2NL$, as each meta-atom needs to calculate its phase once and then apply the modulus normalization. Moreover, according to [56], the overall complexity of the actor network is characterized by

where $u_{a,0}$ denotes the node count in the dense layer of the branch that handles phase shifts for the SIM, and $u_{a,i},\,i=\{1,2,3\}$ denotes the node count in the dense layers of the branch that allocates power in the actor network. Similarly, the complexity of the critic network can be expressed as

where $u_{c,i}$ represents the node count in the $i$-th dense layer of critic network. $u_{c,0}$ equals the input size $|\mathbf{a}|+|\mathbf{s}|$, where $|\mathbf{a}|$ and $|\mathbf{s}|$ are the size of the action and state, respectively.

On the other hand, the complexity of the whitening process is indeed related to the size of the action. Each SIM meta-atom’s coefficient is whitened and then renormalized to meet the conditions in (11c). Thus, the computational complexity of the whitening process is $NL+2NL$, where $2NL$ corresponds to the cost of the renormalization operation. Therefore, the overall complexity of the training phase is

2) Predicting phase: Since the critic network and the experience replay are adopted to ensure the actor network a faster and more stable training process, the complexity of the predicting phase only needs to consider the actor network. Therefore, the complexity of the predicting phase is

Note that the computational complexity of the DRL scheme exhibits linear scaling with the number of meta-atoms and layers in the SIM during both the training and predicting phases. This linear complexity ensures high computational efficiency and scalability, which is suitable for solving large-scale joint parameter optimization problems in SIM-aided wireless communications. Moreover, DNN implementations in practice are highly parallelizable and run in suitable hardware (e.g., graphics processing units) therefore the computational complexity is practically reduced by a near-linear factor.

As shown in Fig. 4, we consider a SIM-assisted multi-user MISO system that operates in the downlink at a frequency of 28 GHz, corresponding to a wavelength of $\lambda=10.7$ mm. The thickness of the SIM is set to $D=5\lambda$ [7, 18], consisting of $L=4$ layers of isomorphic square metasurfaces, with each layer consists of $N=49$ meta-atoms. The area of each meta-atom is set to $s_{a}=\lambda^{2}/4$. In addition, we set the distance between meta-atoms to be $r_{e}=\lambda/2$. The considered simulation setup is depicted in Fig. 4, where the BS is at the height $H_{b}=10$ meter and a SIM is integrated with the BS to facilitate wave-based precoding. The $M$ UEs are randomly distributed in an annular region as shown in Fig. 4 at the beginning of every episode, with the center of the SIM projected onto the ground serving as the annular center. The annulus has an inner radius $R_{l,1}=100$ meter and an outer radius $R_{l,2}=250$ meter. Due to the random distribution of the UEs, we calculate the path losses independently using (5), where we set $C_{0}=-35$ dB and $\alpha=3.5$. The transmission distance of the $m$-th UE can be obtained as $d_{c,m}=(H_{b}^{2}+R_{m}^{2})^{\frac{1}{2}}$, where $R_{m}$ represents the horizontal distance from the $m$-th UE to the annular center. At the transmitter, we set the maximum transmit power to $P=10$ dBm and the noise power is set to $\sigma_{m}^{2}=-104\;\text{dBm},\;\forall m\in\mathcal{M}$. Additionally, we adaptively adjust the number of neurons of the dense layer based on the number of SIM meta-atoms. This strategy aims to mitigate the potential degradation in system performance due to insufficient network capacity. The other relevant parameters are presented in Table II.

To verify the performance of the proposed scheme, we consider the following four benchmark schemes for comparison:

1) DRL-UPA: Optimize the SIM phase shifts based on the proposed DRL algorithm considering the uniform transmit power allocation.

2) Random: The phase shifts of the SIM are configured randomly, while the transmit power of different antennas is allocated via the iterative water-filling algorithm, which follows the water-filling principle [57] in iteration $t$ as

where the power allocated to data stream $m$ in iteration $t$ is denoted as $p_{m,w}^{(t)}$ and the water level $p_{w}^{(t)}$ is determined such that $\sum_{m=1}^{M}p_{m,w}^{(t)}=P$.

3) Codebook: The codebook size corresponds to the number of iterations in the proposed DRL scheme, and for each codeword, we utilize the Random scheme. Its performance is identified as the maximum sum rate attained across all the codeword configurations explored during training.

4) AO: Perform iterative water-filling and gradient ascent approaches to separately optimize transmit power allocation and SIM phase shifts as in [18].

All results are obtained through an average of 100 independent random channel realizations.

In Fig. 5, we depict the average sum rate versus the transmit power $P$ at the BS, by considering $M=4$, $L=4$, and $N=49$. For comparison, we also illustrate the ZF and MMSE fully digital precoding schemes [8] adopting water-filling transmit power allocation without a SIM. As observed, the joint optimization of SIM phase shifts and transmit power allocation via DRL effectively harnesses inter-user interference, yielding superior sum rate performance, of about 2 bps/Hz, compared to the considered AO algorithm for transmit power larger than $P=0$ dBm. At low transmit power levels, ZF and MMSE underperform all SIM-assisted precoding methods with the optimized phase shifts; this happens due to the large SIM aperture. This observation validates the superiority of SIM’s multi-layer structure over traditional precoding schemes in multi-user MISO systems. However, the sum rate of the digital precoding scheme grows faster than SIM with the increasing transmit power, surpassing all other schemes at $P\geq 24$ dBm. In fact, the high transmit power induces an extensive dynamic range in the sum rate, resulting in more significant oscillations and poorer convergence for DRL and AO approaches, at the expense of higher implementation cost for the former.

Fig. 6 illustrates the sum rate versus the number $L$ of SIM layers for two transmit power levels. As shown, the proposed DRL scheme consistently outperforms the considered SIM-aided benchmark schemes in both considered setups. In particular, the average sum rate initially increases and then becomes saturated when the number of SIM layers increases. The initial improvement trend stems from the enhanced interference mitigation precoding capability offered by the multi-layer structure of SIM. Then, the gains diminish and eventually become saturated, since optimizing the numerous parameters becomes intractable.

We show in Fig. 7 the average sum rate versus the transmit power $P$ considering three system settings, i.e., $N=81$, $N=49$, and $N=25$. The proposed DRL algorithm demonstrates superior performance compared to the other three benchmark schemes. At low transmit power levels, the proposed DRL scheme achieves substantially higher sum rates compared to Codebook and AO. Specifically, the proposed scheme with $N=25$ achieves a sum rate performance comparable to the AO algorithm with $N=81$, since the AO algorithm tends to trap into some local suboptimal solutions [18]. In particular, DRL, benefiting from DNN’s powerful computational power, achieves a four-fold increase in average sum rate over Codebook at $P=0$ dBm. However, the elevated transmit power induces an extensive dynamic range at the neural network output of DRL, which directly degrades the performance of the scheme, diminishing the gap versus Codebook scheme. Moreover, setting a large codebook size enables the Codebook scheme to outperform the AO at low $N$.

In Fig. 8, we show the sum rate versus the number $N$ of meta-atoms per layer considering three cases for the number $L$ of SIM layers, while the other hyperparameters are identical to those described in Section IV-A. It can be seen that the sum rate increases with the number of meta-atoms per layer, showing diminishing returns beyond a threshold, e.g., $N=80$. In this regime, the spatial gain attained through the SIM approaches saturation. Specifically, the AO algorithm gradually approaches the performance of the DRL scheme in terms of sum rate, at the cost of an increased computational complexity, which is proportional to the number of meta-atoms. Additionally, both Random and Codebook with random SIM phase shifts schemes attain only mild performance gains with the increase of $L$, as they are inefficient in utilizing the degrees of freedom offered by increasing $N$.

To better demonstrate the optimization capability of the proposed DRL algorithm, Fig. 9a shows its convergence process for three transmit power levels. It can be seen that the transmit power levels significantly affect the final convergence result. Specifically, in our proposed scheme, the improvements in average reward are considerable as the transmit power increases. This is because the proposed DRL can jointly adapt both the SIM phase shifts and transmit power allocation strategy, thus improving the system sum rate performance. However, DRL-UPA, which only optimizes the SIM phase shifts, has a relatively small gain from the increase of transmit power. Furthermore, we compare the average reward versus training steps in Fig. 9b under different numbers of data streams, i.e., $M=\{2,3,4\}$, while other parameters remain the same as those described in Section IV-A. Observed from Fig. 9b that Codebook and Random are unable to efficiently optimize such a massive number of parameters as the number of data streams increases. Thus, adjusting the number of data streams has less impact on Codebook and Random in terms of the convergence speed and resultant performance.

Fig. 10a demonstrates the effect of different learning rates, considering $\gamma_{a}=\gamma_{c}=\{0.001,0.0001,0.00001\}$, and the data stream is set to $M=4$. It is demonstrated that the learning rate dictates the system’s optimal performance and ability to converge. Specifically, when adopting a learning rate of 0.0001, DRL can converge to a relatively satisfactory reward level and further decrease the learning rate over the training process exhibiting a trend of continued reward optimization. This is because, after sufficient training, a relatively lower learning rate could effectively constrain the parameter space and avoid the model learning noise with a higher probability. These benefits allow the model to fine-tune the parameters to find a better solution near the local minimum. However, with a lower learning rate, such as 0.00001, the system is more likely to become trapped in a local optimum because of the lack of exploration competence. Therefore, it is crucial to determine an appropriate learning rate to ensure that the system achieves optimal performance. Moreover, the codebook scheme has a fast convergence faster but performs poorly, whereas the DRL approach achieves a superior sum rate by relying on a well-configured SIM.

Next, we show in Fig. 10b the effect of different batch sizes considering $N_{B}=\{8,16,32,64,128\}$, while the transmit power is $P=10$ dBm and the number of the data streams is $M=4$. During the DRL training process, the batch size does not significantly influence the resultant sum rate performance. Instead, it slightly affects the convergence speed and efficacy of the DRL algorithm. For example, when employing a batch size of 32, the DRL algorithm rapidly converges to the maximum. A larger batch size reduces the gradient’s variance while improving the gradient’s stability. In contrast, when the batch sizes are set to 8 or 16, the algorithm exhibits slower convergence behavior. The convergence behaviors of DRL-UPA are less sensitive to the batch size, which, however, suffers from severe performance loss due to the CSI-unaware power allocation.

Finally, in Fig. 11, we investigate the impact of the whitening process in (20) on the convergence behavior of the DRL algorithm. We note that different parameters involved in the whitening process significantly influence the convergence effect of rewards. Setting $v_{0}=2$ and $t_{\mathrm{gap}}=100$ provides sufficient exploration ability in both early and late training stages, enabling a broader range of attempted actions circumventing over-dependence on proven high-return tactics. This parameterization led to a $20\%$ performance increase compared to the suboptimal values of $v_{0}=0.5,t_{\mathrm{gap}}=100$.

Additionally, the whitening process leads to a smoother DRL training process as it simulates the uncertainty in a real wireless communications environment, thus granting a certain level of robustness to the trained model against uncharted interference. Indeed, the whitening process has the effect of regularization to reduce overfitting and improve the generalization ability of the model. Nevertheless, introducing excessive noise also impedes convergence, which requires careful controlling of the noise variance. Accordingly, attenuation of the whitening process over time is implemented, with Fig. 11 illustrating the efficacy of this attenuation strategy. According to the above analysis, the whitening process plays a pivotal role in the proposed DRL algorithm.