1 Introduction↩︎

The 3rd Generation Partnership Project (3GPP) has been investigating the integration between terrestrial networks (TNs) and non-terrestrial networks (NTNs) since Release 15 [@TR38.811; @TR38.821]. By incorporating the wide coverage of NTN elements, such as geostationary orbit (GEO) and low Earth orbit (LEO) satellites, communication services can be extended far beyond the limitations of terrestrial infrastructure, which could enhance global connectivity. The NTNs can also provide coverage to aerial users such as drones, planes, and urban air mobility vehicles. In the forthcoming 6G standard, 3GPP is expected to make a unified standard for TNs and NTNs. As these integrated networks extend connectivity to diverse users and environments, ensuring the confidentiality of transmissions over satellite links against unauthorized interception becomes an increasingly critical design consideration.

Physical layer security (PLS) exploits the inherent randomness of wireless channels to provide information-theoretic secrecy guarantees, as originally established in [@wyner75]. While PLS has been extensively studied in terrestrial networks [@mukherjee14], its application to NTNs has attracted growing attention in recent years [@JSAC_Zhu; @TIFS_Lei; @TWC_Zheng; @lin19_robust_sat; @guo20_pls_sat; @lin18_cog_sat; @WCL_Li]. For instance, the ergodic secrecy capacity in unmanned aerial vehicle (UAV) networks was studied in [@JSAC_Zhu], where a jamming strategy was proposed to confuse eavesdroppers randomly located on the ground. In [@TIFS_Lei] and [@TWC_Zheng], zero-forcing (ZF)-based beamforming schemes were developed for multi-beam satellite systems, aiming to minimize the satellite’s transmit power while maintaining a secrecy rate constraint. More recently, robust secure beamforming under imperfect channel state information (CSI) was investigated in [@lin19_robust_sat] for multibeam satellite systems, and a threshold-based scheduling scheme for multiuser satellite PLS was proposed in [@guo20_pls_sat]. Joint beamforming designs for cognitive satellite-terrestrial networks were also studied in [@lin18_cog_sat] and [@WCL_Li]. However, these studies focus solely on ground-based eavesdroppers, which intercept the signals from aerial nodes.

Beyond conventional ground-based eavesdroppers, aerial eavesdroppers have recently attracted growing attention, where the secrecy performance against UAV-based eavesdroppers was investigated in [@yuan19_uav_eve; @tang19_uav_eve; @bao20_uav_eve]. Satellite-based eavesdroppers have also emerged as an increasingly relevant threat due to the rapid proliferation of LEO mega-constellations, which increases the likelihood of unauthorized interception from space by reconnaissance satellites or compromised LEO satellites operating in adjacent orbital planes. In [@my22TIFS], the secrecy performance of ground-to-satellite uplink transmissions was analyzed when satellites serve as eavesdroppers using a stochastic geometry framework, showing that orbital geometry significantly affects the secrecy capacity. Unlike UAV eavesdroppers, which are characterized by limited operational range and quasi-static hovering positions, satellite eavesdroppers follow deterministic orbital trajectories determined by Keplerian mechanics with rapidly time-varying channel geometry due to high orbital velocities. These two distinctive properties, namely predictable orbits and rapidly changing channels, pose new security design challenges that are fundamentally different from both ground-based and UAV eavesdropper scenarios.

Reinforcement learning (RL) has been widely adopted for optimizing various aspects of wireless communication systems, such as joint power control and beamforming in terrestrial 5G networks [@23_DRL_5G], dynamic resource allocation for multibeam satellites [@23_sat_rl_beam], and energy-efficient beamforming for integrated satellite-aerial-terrestrial networks [@24_sat_rl_beam]. RL can learn adaptive policies without explicit channel models or per-slot iterative optimization. Designing secure beamforming policies against satellite eavesdroppers involves non-convex optimization problems that must be re-solved at every time slot as the orbital geometry evolves. Conventional iterative algorithms, such as semidefinite relaxation and successive convex approximation (SCA), incur high per-slot computational complexity and require careful per-slot initialization, making them impractical for real-time decision-making in dynamic LEO satellite environments. An RL-based framework is therefore needed to learn beamforming policies that exploit orbital predictability under probabilistic outage constraints with multiple satellite eavesdroppers.

Motivated by this, in this paper, we address this gap by adopting an RL-based approach, where the secure beamforming problem is formulated as a constrained Markov decision process (CMDP) [@book99_cmdp; @06_RL_constrained]. We propose a primal-dual soft actor-critic (PD-SAC) algorithm that learns secure beamforming policies under connection and secrecy outage probability constraints in the presence of multiple satellite eavesdroppers. This approach relies on closed-form outage expressions under Nakagami-\(m\) fading. Combined with a closed-form inequality for the incomplete gamma function, it yields conservative closed-form cost functions that upper bound the outage probabilities using only elementary operations. These cost functions enable gradient-based RL to directly handle probabilistic outage constraints. The proposed PD-SAC algorithm employs Lagrangian relaxation techniques to enforce the outage requirements within the CMDP formulation. The main contributions of this paper are summarized as follows:

• To the best of our knowledge, this is the first study to investigate secure uplink beamforming for LEO satellite communication systems in the presence of multiple satellite eavesdroppers. While prior works have mainly considered ground-based eavesdroppers [@JSAC_Zhu; @TIFS_Lei; @TWC_Zheng; @lin19_robust_sat; @guo20_pls_sat; @lin18_cog_sat; @WCL_Li], they do not capture the distinctive characteristics of satellite eavesdroppers, such as orbital trajectories and time-varying channel geometry. Although satellite eavesdroppers were considered in [@my22TIFS], secure beamforming design was not addressed. In contrast, this work explicitly incorporates multiple satellite eavesdroppers into the uplink beamforming problem and accounts for their geometry-driven, time-varying channels.

• We derive closed-form expressions for the connection and secrecy outage probabilities under Nakagami-\(m\) fading. We develop tractable cost functions based on a closed-form gamma function inequality that serve as conservative upper bounds on the outage probabilities. These cost functions enable the use of gradient-based policy optimization for the probabilistic outage constraints.

• We reformulate the non-convex secrecy rate maximization problem as a CMDP with average outage constraints. We propose a PD-SAC algorithm with a multi-head cost critic that jointly optimizes the beamforming policy through Lagrangian relaxation. We also analyze the computational complexity of the proposed algorithm against conventional beamformers and the offline SCA benchmark, confirming its suitability for real-time deployment.

• We provide simulation results demonstrating that the proposed RL policy outperforms maximum ratio transmission (MRT) and surpasses zero-forcing when the number of eavesdroppers is large, approaching an offline per-slot SCA benchmark while requiring only a single forward pass at inference. The off-policy PD-SAC further outperforms an on-policy primal-dual proximal policy optimization (PD-PPO) counterpart under tight outage constraints, highlighting the role of sample efficiency in constrained RL.

Notations: The superscript \(\mathrm{T}\) indicates the transpose operation. The absolute value of a complex number \(x\) is \(|x|\), and the \(\ell_2\)-norm of a vector \(\mathbf{x}\) is \(\|\mathbf{x}\|\). \(\mathbf{0}\) and \(\mathbf{I}\) denote the all-zero vector and the identity matrix of appropriate dimensions, respectively. The first-kind Bessel function of order \(j\) is \(J_j(\cdot)\). The Gamma function is \(\Gamma(\cdot)\), and the Pochhammer symbol is defined as \((x)_n=\Gamma(x+n)/\Gamma(x)\). The lower incomplete gamma function is defined as \(\gamma(a, x)=\int_0^x t^{a-1}\exp(-t)\,\mathrm{d}t\). The ramp function is \([x]^+=\max(0,x)\). The Hermitian (conjugate) transpose is denoted by \((\cdot)^{\mathrm{H}}\). The inner product of two vectors \(\mathbf{x}\) and \(\mathbf{y}\) is \(\mathbf{x} \cdot \mathbf{y}\). The ceiling function of a real value \(x\) is \(\ceil{x}\). The Kronecker product is \(\otimes\). The Hadamard product is denoted by \(\odot\). The basic 3D rotation matrices about the \(x\)-, \(y\)-, and \(z\)-axes are denoted by \(R_x(\cdot)\), \(R_y(\cdot)\), and \(R_z(\cdot)\), respectively. \(\operatorname{Re}(\cdot)\) and \(\operatorname{Im}(\cdot)\) denote the real and imaginary parts of a complex vector, respectively.

2 System Model↩︎

We consider an uplink satellite communication system in which a ground terminal \(\mathrm{t}\) communicates with a serving satellite \(\mathrm{s}\) in the presence of \(E\) satellite eavesdroppers \(\mathrm{e}_j\), \(j \in \{1,2,\ldots,E\}\). The eavesdroppers are assumed to be non-colluding, i.e., each eavesdropper operates independently without cooperation with others [@my17VTC]. The orbital parameters of potential eavesdropper satellites are assumed to be known, since satellites are physical objects whose orbits are determined by deterministic Keplerian mechanics and can be readily observed and predicted. As illustrated in Fig. 1, the satellite orbits are defined in a global Earth-centered coordinate system \((x,y,z)\), where the origin \(\mathrm{O}\) is located at the Earth’s center and \(r_{\mathrm{E}}\) denotes the Earth’s radius. Since the satellites are assumed to follow circular orbits, the orbit of satellite \(k\) is characterized by four orbital elements \(\mathcal{O}_k= (a_k, \Omega_k, i_k, u_k)\), \(k \in \{\mathrm{s}, \mathrm{e}_j\}\), representing the orbital altitude (or equivalently the semi-major axis), right ascension of the ascending node (RAAN), inclination, and argument of latitude, respectively. The ground terminal is assumed to be fixed on the Earth’s surface, and its location is specified by the latitude \(\phi\) and longitude \(\psi\). As shown in Fig. 2, the terminal is equipped with a uniform planar array (UPA), which is defined in a local coordinate system \((\bar{x},\bar{y},\bar{z})\) centered at the terminal. The UPA consists of \(M \triangleq M_{\bar{x}} \times M_{\bar{y}}\) antenna elements, where \(M_{\bar{x}}\) and \(M_{\bar{y}}\) denote the numbers of antennas along the \(\bar{x}\)- and \(\bar{y}\)-axes, respectively.

The terminal transmits information signals only when the serving satellite is within the 3-dB beamwidth of the satellite’s receive antenna. This restriction ensures that transmissions occur only when the link quality is sufficiently high, while avoiding unnecessary information leakage when the serving satellite is outside the main beam coverage. To characterize the time-varying satellite geometry, the visible period is discretized into \(N_{\mathrm{vis}}\) time slots with a sufficiently small slot interval \(\delta\). The angular evolution of satellite \(k \in \{\mathrm{s},\mathrm{e}_j\}\) between adjacent time slots is characterized by an angular offset \(\Delta_{k}\) from the initial argument of latitude \(u_{k}\), as shown in Fig. 1. Accordingly, the argument of latitude of satellite \(k\) at time slot \(n \in \{1,2,\ldots,N_{\mathrm{vis}}\}\) is given by \({u_{k,n}}= u_{k}+ (n-1)\Delta_{k}\). Since the angular velocity of satellite \(k\) is \(\omega_k = \sqrt{\frac{G M_{\mathrm{E}}}{(r_{\mathrm{E}}+ a_{k})^3}}\), where \(G\) denotes the gravitational constant and \(M_{\mathrm{E}}\) is the mass of the Earth [@book13OrbitalVelocity], the angular offset between adjacent time slots is given by \(\Delta_{k}= \omega_k \delta\).

The Earth-centered inertial (ECI) coordinate system is a fundamental reference frame commonly used in satellite and space object tracking. This system is inertial with respect to distant celestial objects, and thus its axes remain fixed without rotating with Earth. In contrast, the Earth-centered, Earth-fixed (ECEF) coordinate system rotates together with Earth. The position of the ground terminal on the Earth’s surface, specified by the latitude \(\phi\) and longitude \(\psi\), can be expressed in the ECEF coordinate system as \(\Phi_{\mathrm{t}}^{\mathrm{ECEF}} = \left[ r_{\mathrm{E}}\cos \phi \cos \psi,\; r_{\mathrm{E}}\cos \phi \sin \psi,\; r_{\mathrm{E}}\sin \phi \right]^\mathrm{T}\). To represent the terminal position in the ECI frame, the Earth’s rotation must be taken into account. Let \(\omega_{\mathrm{E}}\) denote the Earth’s rotation rate. In time slot \(n\), the Earth rotates by an angle \(\omega_{\mathrm{E}} n \delta\), which can be modeled as a rotation about the \(z\)-axis using the rotation matrix \[R_z(\omega_{\mathrm{E}} n \delta) = \begin{bmatrix} \cos (\omega_{\mathrm{E}} n \delta) & -\sin (\omega_{\mathrm{E}} n \delta) & 0 \\ \sin (\omega_{\mathrm{E}} n \delta) & \cos (\omega_{\mathrm{E}} n \delta) & 0 \\ 0 & 0 & 1 \end{bmatrix}.\] Accordingly, the terminal position in the ECI frame at time slot \(n\) is given by \[\label{eq:pos95ter} \Phi_{\mathrm{t},n} = R_z(\omega_{\mathrm{E}} n \delta)\Phi_{\mathrm{t}}^{\mathrm{ECEF}} = \begin{bmatrix} r_{\mathrm{E}}\cos\phi \cos \psi_n \\ r_{\mathrm{E}}\cos\phi \sin \psi_n \\ r_{\mathrm{E}}\sin\phi \end{bmatrix},\tag{1}\] where \(\psi_n = \psi + \omega_{\mathrm{E}} n \delta\).

We assume that the satellites maintain the boresight of their receive beams fixed toward the subsatellite point, i.e., the nearest point on Earth to the satellite [@TR38.821]. Let \(\nu_{k,n}\) denote the off-boresight angle of satellite \(k\) toward the ground terminal at time slot \(n\), defined as the angle between the terminal direction and the subsatellite point with respect to satellite \(k\), as illustrated in Fig. 1. Furthermore, let \(\nu_k^{\mathrm{3dB}}\) denote the 3-dB beamwidth angle beyond which the receive beam power drops by more than 3 dB. Based on these definitions, the receive antenna gain at satellite \(k\) in time slot \(n\) is modeled as [@antBessel] \(G_{k,n} = G_{k}^{\mathrm{max}}\left(\frac{J_1(g_{k,n})}{2g_{k,n}} + 36\frac{J_3(g_{k,n})}{g_{k,n}^3}\right)^2\), where \(G_{k}^{\mathrm{max}}\) denotes the maximum antenna gain, \(J_1(\cdot)\) and \(J_3(\cdot)\) are the first- and third-order Bessel functions of the first kind, respectively, and \(g_{k,n} \triangleq 2.07123\frac{\sin \nu_{k,n}}{\sin \nu_k^{\mathrm{3dB}}}\). Let \(d_{k,n}\) denote the distance between the terminal and satellite \(k\) at time slot \(n\). Then, the large-scale path loss between the terminal and satellite \(k\) at time slot \(n\) is given by \(\ell_{k,n} = \left(\frac{c}{4\pi f_c}\right)^2 d_{k,n}^{-\kappa}\) , where \(c\) is the speed of light, \(f_c\) is the carrier frequency, and \(\kappa\) is the path-loss exponent.

The characteristics of satellite channels can be accurately modeled by a shadowed-Rician channel model, which explicitly accounts for random shadowing effects. However, its complicated distribution often leads to analytical intractability. Instead, the Nakagami-\(m\) fading model is widely adopted as an effective alternative for tractable analysis in satellite communications [@my23WCL]. Thus, we employ the Nakagami-\(m\) fading model, which can flexibly capture the dominant line-of-sight (LoS) characteristics between the terminal and satellite \(k\) through the fading parameter \(m_k\). Let \(\tilde{h}_{k,n}\) denote the small-scale fading coefficient between the terminal and satellite \(k\) at time slot \(n\). Under the Nakagami-\(m\) model, the channel power gain follows the Gamma distribution, i.e., \(|\tilde{h}_{k,n}|^2 \sim \Gamma(m_k,\frac{1}{m_k})\). Then, the cumulative distribution function (CDF) of the channel power gain \(|\tilde{h}_{k,n}|^2\) is given by \(F_{|\tilde{h}_{k,n}|^2}(x)=\frac{\gamma(m_k,m_k x)}{\Gamma(m_k)}\).

As illustrated in Fig. 2, we let \(\zeta_{k,n}\) and \(\tau_{k,n}\) denote the zenith and azimuth angles of satellite \(k\) at time slot \(n\), respectively, i.e., \((\zeta_{k,n},\tau_{k,n})\) represents the angle-of-departure (AoD) pair from the terminal toward satellite \(k\). Based on the UPA structure, the corresponding array response vector toward satellite \(k\) at time slot \(n\) is given by \(\mathbf{a}_{k,n}=\mathbf{a}_{k,n}^{\bar{x}} \otimes \mathbf{a}_{k,n}^{\bar{y}} \in \mathbb{C}^{M}\), where \(\mathbf{a}_{k,n}^{\bar{x}} = [1, \:\: e^{-j\pi\sin\zeta_{k,n}\cos\tau_{k,n}}, \:\: \cdots, \:\: e^{-j\pi(M_{\bar{x}}-1)\sin\zeta_{k,n}\cos\tau_{k,n}}]^{\mathrm{T}}\) and \(\mathbf{a}_{k,n}^{\bar{y}} = [1, \:\: e^{-j\pi\sin\zeta_{k,n}\sin\tau_{k,n}}, \:\: \cdots, \:\: e^{-j\pi(M_{\bar{y}}-1)\sin\zeta_{k,n}\sin\tau_{k,n}}]^{\mathrm{T}}\) represent the array response vectors along the \(\bar{x}\)- and \(\bar{y}\)-axes, respectively. Stacking the eavesdropper array responses column-wise yields \(\mathbf{A}_{\mathrm{e},n} \triangleq [\mathbf{a}_{\mathrm{e}_1,n}, \ldots, \mathbf{a}_{\mathrm{e}_E,n}] \in \mathbb{C}^{M \times E}\). Consequently, the overall channel vector between the terminal and satellite \(k\) at time slot \(n\) is modeled as [@my24JSAC_EA] \[\begin{align} \label{eq:channel95vector} \mathbf{h}_{k,n}=\tilde{h}_{k,n}\sqrt{\ell_{k,n}}\,\mathbf{a}_{k,n}, \end{align}\tag{2}\] which incorporates the effects of small-scale fading, large-scale path loss, and array geometry. Let \(\mathbf{w}_{n}\in \mathbb{C}^{M}\) denote the beamforming vector employed by the terminal at time slot \(n\), subject to the transmit power constraint \(\|\mathbf{w}_{n}\|^2 \leq P_{\mathrm{max}}\), where \(P_{\mathrm{max}}\) is the maximum transmit power. Then, the received signal-to-noise ratio (SNR) at satellite \(k\) in time slot \(n\) is expressed as \[\begin{align} \label{eq:SNR} \Gamma_{k,n}= \begin{cases} \dfrac{G_{k,n}\left|\mathbf{h}_{k,n}^{\mathrm{H}}\mathbf{w}_{n}\right|^2}{N_0 W}, & \text{if } d_{k,n}<d_k^{\mathrm{max}}, \\ 0, & \text{otherwise}, \end{cases} \end{align}\tag{3}\] where \(d_k^{\mathrm{max}}\) denotes the maximum distance for which satellite \(k\) remains visible due to Earth blockage, given by \(d_k^{\mathrm{max}}=\sqrt{a_{k}(2r_{\mathrm{E}}+a_{k})}\) [@my22TIFS], \(N_0\) denotes the noise power spectral density, and \(W\) is the system bandwidth. For non-colluding eavesdroppers, the secrecy performance is dominated by the most detrimental eavesdropper, i.e., the one achieving the highest received SNR among all eavesdroppers [@my17VTC]. Accordingly, the instantaneous secrecy rate at time slot \(n\) is given by \[\begin{align} \label{eq:sec95rate} R_n=\left[\log_2\left(1+\Gamma_{\mathrm{s},n}\right)-\log_2\left(1+ \underset{j\in\{1,2,\ldots,E\}}{\max} \Gamma_{\mathrm{e}_j,n}\right)\right]^+. \end{align}\tag{4}\] The secrecy rate in 4 is used to formulate the optimization problem in Section 4.

3 Mathematical Preliminaries↩︎

In this section, we first characterize the distance between the terminal and satellite \(k\in\{\mathrm{s},\mathrm{e}_j\}\), \(j\in\{1,2,\cdots,E\}\), in time slot \(n\), i.e., \(d_{k,n}\). Then, we analyze the satellite visibility constraint based on the orbital configuration of the satellites. Additionally, we derive analytical expressions for the off-boresight angle \(\nu_{k,n}\), the zenith angle \(\zeta_{k,n}\), and the azimuth angle \(\tau_{k,n}\).

3.1 Distance Characterization↩︎

The position of satellite \(k\in\{\mathrm{s}, \mathrm{e}_j\}\) in time slot \(n\) is obtained by using three successive intrinsic rotations⁴ with sequence \(z-x^{\prime}-z^{\prime\prime}\) as [@my23WCL] \[\begin{align} \label{eq:pos95sat} & \Phi_{k,n} = R_z(\Omega_{k})R_x(i_{k})R_z({u_{k,n}}) \mu_{x,k}\nonumber \\ & = (r_{\mathrm{E}}+a_{k}) \begin{bmatrix} \cos\Omega_{k}\cos{u_{k,n}}-\sin\Omega_{k}\cos i_{k}\sin{u_{k,n}}\\ \cos\Omega_{k}\cos i_{k}\sin{u_{k,n}}+\sin\Omega_{k}\cos{u_{k,n}}\\ \sin i_{k}\sin{u_{k,n}} \end{bmatrix}, \end{align}\tag{5}\] where \(\mu_{x,k}=[r_{\mathrm{E}}+a_{k}\:\:\: 0 \:\:\: 0]^{\mathrm{T}}\). From 1 and 5 , the distance between the terminal and satellite \(k\in\{\mathrm{s}, \mathrm{e}_j\}\) in time slot \(n\) is given by \[\begin{align} \label{eq:dist} & d_{k,n}= \Vert \Phi_{k,n}-\Phi_{\mathrm{t},n} \Vert \nonumber \\ & =\Big( (r_{\mathrm{E}}+ a_k)^2 \Big(\eta(\cdot;i_k^2,u_{k,n}^2)+ \eta(\Omega_k^2,i_k^2;u_{k,n}^2)\nonumber \\ & + \eta(u_{k,n}^2;\Omega_k^2) + \eta(\Omega_k^2,u_{k,n}^2;\cdot) +\eta(i_k^2;\Omega_k^2,u_{k,n}^2) \Big) \nonumber \\ & - 2(r_{\mathrm{E}}+ a_k) r_{\mathrm{E}}\Big(\eta(\Omega_k,i_k,\phi;u_{k,n},\psi_n) -\eta(i_k,\phi,\psi_n;\Omega_k,u_{k,n}) \nonumber \\ & + \eta(\Omega_k,u_{k,n},\phi,\psi_n;\cdot) + \eta(u_{k,n},\phi;\Omega_k,\psi_n) + \eta(\cdot;i_k,u_{k,n},\phi) \Big) \nonumber \\ & + r_{\mathrm{E}}^2 \Big( \eta(\phi^2,\psi_n^2;\cdot) + \eta(\phi^2;\psi_n^2) + \eta(\cdot;\phi^2) \Big) \Big)^{1/2}, \end{align}\tag{6}\] where \(\eta(p_1^{o_1},\cdots,p_{\mathrm{A}}^{o_\mathrm{A}};q_1^{v_1},\cdots,q_{\mathrm{B}}^{v_\mathrm{B}}) \mathrel{\ensurestackMath{\stackon[1pt]{=}{\scriptstyle\Delta}}}\prod_{t=1}^{\mathrm{A}}\cos^{o_t}{p_t}\times \prod_{t=1}^{\mathrm{B}}\sin^{v_t}{q_t}\) is a multiplication of cosine and sine functions. These distances vary over time due to satellite mobility and directly affect the path loss \(\ell_{k,n}\).

3.2 Satellite Visibility Analysis↩︎

Since the orbital plane of the serving satellite is obtained by two intrinsic rotations with orbital parameters \(\Omega_{\mathrm{s}}\) and \(i_{\mathrm{s}}\), the normal vector of the orbital plane is obtained as \[\begin{align} \mathbf{n}_{\mathrm{s}}^{\perp} = R_z(\Omega_{\mathrm{s}})R_x(i_{\mathrm{s}}) \hat{z} = \begin{bmatrix} \sin \Omega_{\mathrm{s}}\sin i_{\mathrm{s}}\\ -\cos \Omega_{\mathrm{s}}\sin i_{\mathrm{s}}\\ \cos i_{\mathrm{s}} \end{bmatrix}, \end{align}\] where \(\hat{z}=[0 \:\: 0 \:\: 1]^{\mathrm{T}}\). The angle between the terminal and the normal vector is obtained as \[\begin{align} \beta^{\prime} & =\arccos \frac{\mathbf{n}_{\mathrm{s}}^{\perp} \cdot \Phi_{\mathrm{t},n}}{\lVert \mathbf{n}_{\mathrm{s}}^{\perp} \rVert \cdot \lVert \Phi_{\mathrm{t},n} \rVert}\nonumber \\ & = \arccos (\cos \phi \cos \psi_n \sin \Omega_{\mathrm{s}}\sin i_{\mathrm{s}}\nonumber \\ & \quad - \cos \phi \sin \psi_n \cos \Omega_{\mathrm{s}}\sin i_{\mathrm{s}}+ \sin \phi \cos i_{\mathrm{s}}), \end{align}\] and the angle between the terminal and the orbital plane is derived as \[\begin{align} \beta=\Big\vert\beta^{\prime}-\frac{\pi}{2}\Big\vert. \end{align}\] According to [@my24TWC], the orbital plane is visible only when the following criterion holds: \[\begin{align} \label{eq:vis95criterion} \beta < \arccos\left(\frac{r_{\mathrm{E}}}{r_{\mathrm{E}}+a_{\mathrm{s}}}\right), \end{align}\tag{7}\] and the length of the visible arc of the orbit is given by \[\begin{align} \label{eq:len95vis95arc} l =2(r_{\mathrm{E}}+a_{\mathrm{s}})\arcsin\left(\sqrt{1-\frac{\sec^2\beta}{(1+a_{\mathrm{s}}/r_{\mathrm{E}})^2}}\right). \end{align}\tag{8}\]

Using 8 , the angle between the two endpoints of the visible arc is obtained as \[\begin{align} \label{eq:thetas} \Theta = \frac{l}{r_{\mathrm{E}}+a_{\mathrm{s}}}=2\arcsin\left(\sqrt{1-\frac{\sec^2\beta}{(1+a_{\mathrm{s}}/r_{\mathrm{E}})^2}}\right). \end{align}\tag{9}\] From 9 and the angular offset \(\Delta_{\mathrm{s}}\) between consecutive time slots, the maximum number of visible time slots is obtained as \[\begin{align} N_{\mathrm{vis}} & = \left\lceil \frac{\Theta}{\Delta_{\mathrm{s}}} \right\rceil = \left\lceil \frac{l}{(r_{\mathrm{E}}+a_{\mathrm{s}})\omega_{\mathrm{s}}\delta} \right\rceil \nonumber \\ & = \left\lceil \sqrt{\frac{4(r_{\mathrm{E}}+ a_{\mathrm{s}})^3}{\delta^2 G M_{\mathrm{E}}}}\arcsin\left(\sqrt{1-\frac{\sec^2\beta}{(1+a_{\mathrm{s}}/r_{\mathrm{E}})^2}}\right) \right\rceil. \end{align}\]

Now, we derive the visible range of the argument of latitude \(u_{\mathrm{s}}\) based on the fact that the maximum elevation angle is achieved at the position in orbit with the minimum distance to the terminal.

3.3 Geometric Pointing Angle Characterization↩︎

In this subsection, we derive analytical expressions for the off-boresight angle \(\nu_{k,n}\), the zenith angle \(\zeta_{k,n}\), and the azimuth angle \(\tau_{k,n}\) of satellite \(k\). The off-boresight angle \(\nu_{k,n}\) is defined as the angle between the boresight direction and the terminal with respect to satellite \(k\). To derive the off-boresight angle of satellite \(k\), we apply the law of cosines to the triangle \(\triangle \mathrm{O} \Phi_{\mathrm{t},n} \Phi_{k,n}\), i.e., \(\overline{\mathrm{O}\Phi_{\mathrm{t},n}}^2 =\overline{\Phi_{\mathrm{t},n}\Phi_{k,n}}^2+\overline{\mathrm{O}\Phi_{k,n}}^2-2\overline{\Phi_{\mathrm{t},n}\Phi_{k,n}}\cdot\overline{\mathrm{O}\Phi_{k,n}} \cdot\cos \nu_{k,n}\). Since \(\overline{\mathrm{O}\Phi_{\mathrm{t},n}} = r_{\mathrm{E}}\), \(\overline{\Phi_{\mathrm{t},n}\Phi_{k,n}} = d_{k,n}\), and \(\overline{\mathrm{O}\Phi_{k,n}} = r_{\mathrm{E}}+ a_k\), the off-boresight angle is expressed as \[\begin{align} \nu_{k,n} = \arccos\left(\frac{d_{k,n}^2 + (r_{\mathrm{E}}+ a_k)^2 - r_{\mathrm{E}}^2}{2 d_{k,n} (r_{\mathrm{E}}+ a_k)}\right). \end{align}\]

The zenith angle \(\zeta_{k,n}\) is defined as the angle between the local vertical (zenith) direction at the terminal and the direction toward satellite \(k\). It is calculated as \[\begin{align} \zeta_{k,n} = \arccos\left(\frac{\Phi_{\mathrm{t},n} \cdot (\Phi_{k,n} - \Phi_{\mathrm{t},n})}{\lVert \Phi_{\mathrm{t},n} \rVert \cdot \lVert \Phi_{k,n} - \Phi_{\mathrm{t},n} \rVert}\right). \end{align}\]

To derive the azimuth angle \(\tau_{k,n}\) with respect to the terminal, we first define the equation of the horizontal plane at the terminal’s position \(\Phi_{\mathrm{t},n}\) as \(\Phi_{\mathrm{t},n} \cdot ([x, y, z]^{\mathrm{T}} - \Phi_{\mathrm{t},n}) = 0\), which is equivalent to \[\begin{align} \label{eq:horizontal95plane} & (\cos \phi \cos \psi_n) x + (\cos \phi \sin \psi_n) y + (\sin \phi) z = r_{\mathrm{E}}. \end{align}\tag{10}\] This plane corresponds to the \(\bar{x} \bar{y}\)-plane in the local coordinate system (\(\bar{x}\), \(\bar{y}\), \(\bar{z}\)) where the terminal’s antennas are located. Let \(\Phi^{\prime}_{k,n}\) denote the projection of \(\Phi_{k,n}\) onto this horizontal plane, which is given by \[\begin{align} \Phi^{\prime}_{k,n} = \Phi_{k,n} - \frac{\Phi_{\mathrm{t},n} \cdot (\Phi_{k,n} - \Phi_{\mathrm{t},n})}{\| \Phi_{\mathrm{t},n} \|^2} \Phi_{\mathrm{t},n}. \end{align}\] To determine the reference direction for the planar array, we use an arbitrary point \(\mathrm{H}\) on the horizontal plane 10 , as illustrated in Fig. 2. This reference direction indicates the orientation of the planar array on the horizontal plane at the terminal. Without loss of generality, we set \(\mathrm{H} = \left[\frac{r_{\mathrm{E}}}{\cos \phi \cos \psi_n}, 0, 0\right]^{\mathrm{T}}\). Thus, the azimuth angle is obtained as \[\begin{align} \tau_{k,n} = \arccos\left(\frac{(\Phi^{\prime}_{k,n} - \Phi_{\mathrm{t},n}) \cdot (\mathrm{H} - \Phi_{\mathrm{t},n})}{\|\Phi^{\prime}_{k,n} - \Phi_{\mathrm{t},n}\| \cdot \|\mathrm{H} - \Phi_{\mathrm{t},n}\|}\right). \end{align}\] The above \(\arccos\)-based formulation illustrates the geometric relationship for general terminal latitudes where the reference point \(\mathrm{H}\) is well-defined. The angles \(\nu_{k,n}\), \(\zeta_{k,n}\), and \(\tau_{k,n}\) determine \(G_{k,n}\) and \(\mathbf{a}_{k,n}\), both of which vary with the time slot \(n\).

4 Outage-Constrained Secrecy Rate Maximization Problem↩︎

In this section, we first formulate the secrecy rate maximization problem subject to transmit power and probabilistic outage constraints over the transmission slots. We then derive closed-form expressions for the connection and secrecy outage probabilities under the Nakagami-\(m\) fading. Finally, we transform the probabilistic constraints into deterministic beam gain constraints, yielding a tractable problem reformulation that serves as the basis for the RL-based solution developed in Section 5.

As discussed in Section 2, the terminal transmits uplink signals only when the serving satellite is within the 3-dB beamwidth, i.e., \(\nu_{\mathrm{ser},n} < \nu_{\mathrm{ser}}^{\mathrm{3dB}}\). Therefore, the number of transmission slots, denoted by \(N\), must always be less than the number of visible slots \(N_{\mathrm{vis}}\), i.e., \(N \leq N_{\mathrm{vis}}\). The objective is to maximize the expected secrecy rate while satisfying the transmit power constraint and average outage probability constraints over the \(N\) transmission slots. The secrecy rate maximization problem is formulated as follows: \[\begin{align} (\mathrm{P1}) \quad \underset{\{\mathbf{w}_n\}_{n=1}^{N}}{\text{maximize}} \quad & \mathbb{E}\left[\sum_{n=1}^{N} R_n(\mathbf{w}_n)\right] \tag{11} \\ \text{subject to} \quad & \|\mathbf{w}_n\|^2 \leq P_{\mathrm{max}}, \quad \forall n, \tag{12} \\ & \frac{1}{N}\sum_{n=1}^{N} P_n^{\mathrm{co}}(\mathbf{w}_n) \leq \epsilon_{\mathrm{co}}, \tag{13} \\ & \frac{1}{N}\sum_{n=1}^{N} P_n^{\mathrm{so}}(\mathbf{w}_n) \leq \epsilon_{\mathrm{so}}, \tag{14} \end{align}\] where the expectation in 11 is over the small-scale fading realizations, \(P_n^{\mathrm{co}}(\mathbf{w}_n)\) denotes the connection outage probability at slot \(n\), which is the probability that the serving satellite cannot decode the message, \(P_n^{\mathrm{so}}(\mathbf{w}_n)\) denotes the secrecy outage probability at slot \(n\), defined as the probability that at least one eavesdropper can intercept the message, \(\epsilon_{\mathrm{co}} \in (0,1)\) is the maximum tolerable average connection outage probability, and \(\epsilon_{\mathrm{so}} \in (0,1)\) is the maximum tolerable average secrecy outage probability. The average formulation 13 and 14 is well-suited to the transmission slots because the terminal’s objective is to maintain reliable and secure communication on average over the transmission slots. This allows the beamforming policy to allocate resources adaptively across time slots with varying channel conditions. Closed-form expressions for \(P_n^{\mathrm{co}}\) and \(P_n^{\mathrm{so}}\) are derived in the following subsection to enable tractable evaluation of the constraints 13 and 14 .

4.1 Outage Probability Analysis↩︎

To characterize the per-slot outage probabilities appearing in 13 and 14 , we derive closed-form expressions under the Nakagami-\(m\) fading assumption.

4.1.1 Connection Outage Probability↩︎

The connection outage probability is defined as the probability that the instantaneous data rate at the serving satellite \(\mathrm{s}\) falls below a target rate \(R_{\mathrm{s},n}\). Mathematically, it is given by \[\begin{align} \label{eq:Pco1} P_n^{\mathrm{co}} & = \mathbb{P}\left(\log_2(1+\Gamma_{\mathrm{s},n}) < R_{\mathrm{s},n}\right)\nonumber \\ & =\mathbb{P}\left( \big| \mathbf{h}_{\mathrm{s},n}^{\mathrm{H}} \mathbf{w}_n \big|^2 < \frac{(2^{R_{\mathrm{s},n}} - 1) N_0 W}{G_{\mathrm{s},n}} \right)\nonumber \\ & \mathop{=}^{(a)}\mathbb{P}\left( |\tilde{h}_{\mathrm{s},n}|^2 < \frac{(2^{R_{\mathrm{s},n}} - 1) N_0 W}{G_{\mathrm{s},n} \ell_{\mathrm{s},n} \big| \mathbf{a}_{\mathrm{s},n}^{\mathrm{H}} \mathbf{w}_n \big|^2} \right) \nonumber \\ & \mathop{=}^{(b)} \frac{1}{\Gamma(m_{\mathrm{s}})}\gamma\left(m_{\mathrm{s}}, \frac{m_{\mathrm{s}}(2^{R_{\mathrm{s},n}} - 1) N_0 W}{G_{\mathrm{s},n} \ell_{\mathrm{s},n} \big| \mathbf{a}_{\mathrm{s},n}^{\mathrm{H}} \mathbf{w}_n \big|^2}\right), \end{align}\tag{15}\] where (\(a\)) follows from the channel decomposition \(|\mathbf{h}_{\mathrm{s},n}^{\mathrm{H}} \mathbf{w}_n|^2 = \ell_{\mathrm{s},n} |\tilde{h}_{\mathrm{s},n}|^2 |\mathbf{a}_{\mathrm{s},n}^{\mathrm{H}} \mathbf{w}_n|^2\), and (\(b\)) follows from the CDF of the Gamma distribution.

4.1.2 Secrecy Outage Probability↩︎

The individual secrecy outage probability for the \(j\)-th eavesdropper, denoted \(P_{j,n}^{\mathrm{so}}\), is defined as the probability that the instantaneous achievable rate at that eavesdropper exceeds a secrecy threshold \(R_{\mathrm{e}_j,n}\). Here, the threshold \(R_{\mathrm{e}_j,n}\) represents the redundancy rate allocated through wiretap coding to protect the confidential information, rather than the data rate intended for the eavesdropper or the directly tolerable leakage rate. Therefore, a secrecy outage occurs when the eavesdropper’s instantaneous channel capacity exceeds this redundancy margin, indicating that the confidential message may no longer be fully protected from information leakage. Similar to the derivation of 15 , the closed-form expression is obtained as \[\begin{align} \label{eq:Pso95j} P_{j,n}^{\mathrm{so}} & =\mathbb{P}(\log_2(1+\Gamma_{\mathrm{e}_j,n}) > R_{\mathrm{e}_j,n})\nonumber \\ & =1-\frac{1}{\Gamma(m_{\mathrm{e}_j})}\gamma\left(m_{\mathrm{e}_j}, \frac{m_{\mathrm{e}_j}(2^{R_{{\mathrm{e}_j},n}} - 1) N_0 W}{G_{{\mathrm{e}_j},n} \ell_{{\mathrm{e}_j},n} \big| \mathbf{a}_{{\mathrm{e}_j},n}^{\mathrm{H}} \mathbf{w}_n \big|^2}\right). \end{align}\tag{16}\] Under the assumption that the eavesdroppers operate independently without collusion, a secrecy outage event occurs when at least one of the \(E\) eavesdroppers successfully intercepts the transmitted message. Since the fading channels across different eavesdroppers are statistically independent, due to their physical separation across distinct orbital planes specified by \(\mathcal{O}_{\mathrm{e}_j} = (a_{\mathrm{e}_j}, \Omega_{\mathrm{e}_j}, i_{\mathrm{e}_j}, u_{\mathrm{e}_j})\), the probability that no eavesdropper succeeds is given by the product of individual complement probabilities. The overall secrecy outage probability is therefore expressed as \[\begin{align} \label{eq:Pso95all} P^{\mathrm{so}}_n & = 1 - \mathbb{P}\left( \bigcap_{j=1}^{E} \left\{ \log_2\left(1+\Gamma_{{\mathrm{e}_j},n}\right) \leq R_{{\mathrm{e}_j},n} \right\}\right) \nonumber \\ & = 1 - \prod_{j=1}^{E} \mathbb{P}\left(\log_2(1+\Gamma_{{\mathrm{e}_j},n}) \leq R_{{\mathrm{e}_j},n}\right) \nonumber \\ & = 1 - \prod_{j=1}^{E} \left(1 - P^{\mathrm{so}}_{j,n}\right), \end{align}\tag{17}\] where the second equality holds due to the independence of the fading channels across different eavesdroppers. From 16 and 17 , \(P_n^{\mathrm{so}}\) is monotonically increasing in the beamforming gain \(|\mathbf{a}_{\mathrm{e}_j,n}^{\mathrm{H}} \mathbf{w}_n|^2\) toward each eavesdropper \(j\). This observation motivates the design of beamforming vectors that suppress signal leakage toward potential eavesdroppers.

4.2 Probabilistic Constraint Reformulation↩︎

The outage probability constraints in 13 and 14 involve the lower incomplete gamma function, which does not admit a closed-form expression amenable to direct optimization. To address this difficulty, a conservative bound on the regularized lower incomplete gamma function is employed, as stated in the following lemma.

Applying the upper bound from Lemma 2 to 15 , we obtain the upper bound of the connection outage probability as \[\begin{align} \label{eq:Pco95bound} P_n^{\mathrm{co}} = \frac{\gamma(m_{\mathrm{s}}, x_{\mathrm{s},n})}{\Gamma(m_{\mathrm{s}})} \leq \left(1 - e^{-x_{\mathrm{s},n}}\right)^{m_{\mathrm{s}}} \triangleq c_{\mathrm{co},n}, \end{align}\tag{18}\] where \(x_{\mathrm{s},n} \mathrel{\ensurestackMath{\stackon[1pt]{=}{\scriptstyle\Delta}}}\frac{m_{\mathrm{s}} (2^{R_{\mathrm{s},n}} - 1) N_0 W}{G_{\mathrm{s},n} \ell_{\mathrm{s},n} |\mathbf{a}_{\mathrm{s},n}^{\mathrm{H}} \mathbf{w}_n|^2}\). For the secrecy outage, we apply the lower bound in ?? to each eavesdropper’s CDF in 16 , which gives \(1 - P_{j,n}^{\mathrm{so}} \geq (1 - e^{-(m_{\mathrm{e}_j}!)^{-1/m_{\mathrm{e}_j}} x_{\mathrm{e}_j,n}})^{m_{\mathrm{e}_j}}\). Substituting this into 17 yields the upper bound of the secrecy outage probability, i.e., \[\begin{align} \label{eq:Pso95bound} P_n^{\mathrm{so}} \leq 1 - \prod_{j=1}^{E} \left(1 - e^{-(m_{\mathrm{e}_j}!)^{-1/m_{\mathrm{e}_j}} x_{\mathrm{e}_j,n}}\right)^{m_{\mathrm{e}_j}} \triangleq c_{\mathrm{so},n}. \end{align}\tag{19}\] where \(x_{\mathrm{e}_j,n} \mathrel{\ensurestackMath{\stackon[1pt]{=}{\scriptstyle\Delta}}}\frac{m_{\mathrm{e}_j} (2^{R_{\mathrm{e}_j,n}} - 1) N_0 W}{G_{\mathrm{e}_j,n} \ell_{\mathrm{e}_j,n} |\mathbf{a}_{\mathrm{e}_j,n}^{\mathrm{H}} \mathbf{w}_n|^2}\), analogous to \(x_{\mathrm{s},n}\). Since \(P_n^{\mathrm{co}} \leq c_{\mathrm{co},n}\) and \(P_n^{\mathrm{so}} \leq c_{\mathrm{so},n}\) at every time slot, the averages satisfy \(\frac{1}{N}\sum_n P_n^{\mathrm{co}} \leq \frac{1}{N}\sum_n c_{\mathrm{co},n}\) and \(\frac{1}{N}\sum_n P_n^{\mathrm{so}} \leq \frac{1}{N}\sum_n c_{\mathrm{so},n}\). Therefore, imposing constraints on the average of \(c_{\mathrm{co},n}\) and \(c_{\mathrm{so},n}\) constitutes a sufficient condition for the original average outage constraints 13 and 14 . The connection and secrecy outage bounds exhibit asymmetric tightness, yet both remain conservative. The resulting cost functions are closed-form expressions of elementary functions, readily amenable to gradient-based policy optimization in the CMDP framework developed in the following section.

5 Proposed RL-Based Beamformer Design↩︎

The optimization problem (P1) is non-convex due to the non-concave secrecy rate objective and the dependence between beamforming vectors across time slots. Iterative methods such as SCA face two issues in this setting. First, the per-slot cost is high, since each slot requires repeatedly linearizing the non-concave secrecy rate and solving a dimension-\(2M\) quadratic program from multiple initializations. Second, the constraints are handled per-slot rather than as time-averaged budgets. We therefore use SCA only as an offline benchmark in Section 6.

To address these limitations, we adopt an RL-based approach in which the beamforming vector is obtained by a single forward pass of the policy network, eliminating the per-slot iterative optimization. The problem is reformulated as a CMDP and solved via a PD-SAC algorithm.

5.1 CMDP Formulation↩︎

As discussed in Section 2, the terminal transmits only when it is located within the serving satellite’s coverage region, i.e., \(\nu_{\mathrm{s},n} < \nu_{\mathrm{s}}^{\mathrm{3dB}}\). In the RL framework, the agent acts only during these transmission slots, and the remaining slots are skipped without agent interaction. The sequential beamforming optimization is modeled as a CMDP defined by the tuple \((\mathcal{S}, \mathcal{A}, P, r, \mathbf{c}, \boldsymbol{\epsilon})\). Here, \(\mathcal{S}\) is the state space, \(\mathcal{A}\) is the action space, \(P: \mathcal{S} \times \mathcal{A} \times \mathcal{S} \to [0,1]\) is the transition probability function, and \(r: \mathcal{S} \times \mathcal{A} \to \mathbb{R}\) is the per-step reward function. The vector \(\mathbf{c} = [c_{\mathrm{co}}, c_{\mathrm{so}}]^{\mathrm{T}}\) collects the per-step cost functions for connection and secrecy outage, with realization \(\mathbf{c}_n = [c_{\mathrm{co},n}, c_{\mathrm{so},n}]^{\mathrm{T}}\) at slot \(n\), and \(\boldsymbol{\epsilon} = [\epsilon_{\mathrm{co}}, \epsilon_{\mathrm{so}}]^{\mathrm{T}}\) is the vector of corresponding cost thresholds.

5.1.1 State Space↩︎

The state \(s_n \in \mathcal{S}\) at time slot \(n\) captures the geometry-derived channel characteristics, the previous-slot secrecy rate, and temporal position within the transmission slots, i.e., \[\begin{align} \label{eq:state95space} s_n = \bigl\{ & \bar{R}_{n-1},\, \tfrac{n}{N},\, \tilde{\ell}_{\mathrm{s},n},\, \tilde{G}_{\mathrm{s},n},\, \tilde{\boldsymbol{\ell}}_{\mathrm{e},n},\, \tilde{\mathbf{G}}_{\mathrm{e},n},\, \nonumber \\ & \operatorname{Re}(\mathbf{a}_{\mathrm{s},n}),\, \operatorname{Im}(\mathbf{a}_{\mathrm{s},n}),\, \operatorname{Re}(\mathbf{A}_{\mathrm{e},n}),\, \operatorname{Im}(\mathbf{A}_{\mathrm{e},n}) \bigr\}. \end{align}\tag{20}\] Here, \(\bar{R}_{n-1} \approx [\log_2(1 + \bar{\Gamma}_{\mathrm{s},n-1}) - \log_2(1 + \max_j \bar{\Gamma}_{\mathrm{e}_j,n-1})]^+\) is the long-term secrecy rate from the previous time slot with \(\bar{\Gamma}_{k,n-1} = G_{k,n-1} \ell_{k,n-1} |\mathbf{a}_{k,n-1}^{\mathrm{H}} \mathbf{w}_{n-1}|^2 / (N_0 W)\). Since accurate instantaneous CSI is difficult to obtain in practice for satellite links, \(\bar{R}_{n-1}\) is constructed solely from the satellite orbital geometry and the previously applied beamforming vector, both of which are known to the terminal, thereby providing a deterministic estimate of the achievable secrecy performance without requiring real-time channel estimation. The term \(n/N \in [0,1]\) is the normalized time-slot index, indicating its relative position in time. The normalized path losses are defined as \(\tilde{\ell}_{k,n} = (10\log_{10}(\ell_{k,n}) + 200)/20\) to scale the values near zero for improved neural network training stability, with \(\tilde{\boldsymbol{\ell}}_{\mathrm{e},n} = [\tilde{\ell}_{\mathrm{e}_1,n}, \ldots, \tilde{\ell}_{\mathrm{e}_{E},n}]^{\mathrm{T}}\). The normalized antenna gains \(\tilde{G}_{k,n} = 10\log_{10}(G_{k,n})/G_{\max}^{\mathrm{dBi}}\), where \(G_{\max}^{\mathrm{dBi}} \triangleq 10\log_{10}(G_{k}^{\mathrm{max}})\) is the peak antenna gain in dBi, reflect the angular gain attenuation toward each satellite, with \(\tilde{\mathbf{G}}_{\mathrm{e},n} = [\tilde{G}_{\mathrm{e}_1,n}, \ldots, \tilde{G}_{\mathrm{e}_{E},n}]^{\mathrm{T}}\). Concatenating all components and splitting the complex array responses into real and imaginary parts yields a real-valued state vector of dimension \(D_\mathrm{s} = 4 + 2E + 2M + 2EM\).

5.1.2 Action Space↩︎

To enforce the power constraint \(\|\mathbf{w}_n\|^2 \leq P_{\mathrm{max}}\) while keeping the actor network differentiable, the actor outputs a latent action that consists of the beamforming direction \(\tilde{q}_{n}^{\mathrm{dir}} \in \mathbb{R}^{2M}\) and the transmit power \(\tilde{q}_{n}^{\mathrm{pow}} \in \mathbb{R}\), and is deterministically mapped to the beamforming vector. Here, \(\tilde{q}_{n}^{\mathrm{dir}}\) stacks the real and imaginary parts of \(\mathbf{w}_n \in \mathbb{C}^{M}\), a Cartesian representation chosen to avoid the discontinuity of cyclic phase variables. These two components are concatenated into the latent action vector \(\tilde{q}_n = [(\tilde{q}_{n}^{\mathrm{dir}})^{\mathrm{T}},\, \tilde{q}_{n}^{\mathrm{pow}}]^{\mathrm{T}} \in \mathbb{R}^{2M+1}\). The direction components are bounded by \(\tanh(\cdot)\) and normalized to a unit vector, while the power component is mapped to \((0, P_{\mathrm{max}})\) via the sigmoid function \(\sigma(\cdot)\) as follows: \[\begin{align} \hat{\mathbf{d}}_n & = \frac{\tanh(\tilde{q}_{n}^{\mathrm{dir}})}{\|\tanh(\tilde{q}_{n}^{\mathrm{dir}})\|}, \tag{21} \\ P_n & = P_{\mathrm{max}} \cdot \sigma(\tilde{q}_{n}^{\mathrm{pow}}), \tag{22} \\ \mathbf{w}_n & = \sqrt{P_n}\, \hat{\mathbf{d}}_n. \tag{23} \end{align}\] With a slight abuse of notation, \(\mathbf{w}_n\) denotes either the complex beamforming vector \(\mathbf{w}_n \in \mathbb{C}^{M}\) from Section 2 or its real-valued representation \([\mathrm{Re}(\mathbf{w}_n)^{\mathrm{T}}, \mathrm{Im}(\mathbf{w}_n)^{\mathrm{T}}]^{\mathrm{T}} \in \mathbb{R}^{2M}\) used as the neural network input. The critic networks receive \(\mathbf{w}_n \in \mathbb{R}^{2M}\). The policy samples the latent action \(\tilde{q}_n \in \mathbb{R}^{2M+1}\) from a diagonal Gaussian, and its squashed counterpart \((u_n, v_n) = (\tanh(\tilde{q}_n^{\mathrm{dir}}),\, P_{\mathrm{max}}\sigma(\tilde{q}_n^{\mathrm{pow}}))\) defines the policy density and entropy used in the SAC update.

5.1.3 Reward and Cost Functions↩︎

The objective is to maximize the expected cumulative secrecy rate. The per-step reward is the instantaneous secrecy rate defined in 4 , which includes the small-scale fading realization: \(r_n(s_n, \mathbf{w}_n) = R_n(\mathbf{w}_n).\) Note that \(\bar{R}_{n-1}\) uses only the deterministic geometry, whereas the reward \(R_n\) reflects the actual channel realization including small-scale fading. Since the policy is trained offline through simulation, the reward can be computed using the actual channel realization generated in the simulator, including small-scale fading. The per-step cost functions are the upper-bound cost functions \(c_{\mathrm{co},n}\) and \(c_{\mathrm{so},n}\) derived in Section 4.2. The power consumption is not modeled as a cost function, since the sigmoid parameterization 22 ensures \(\|\mathbf{w}_n\|^2 = P_n \leq P_{\mathrm{max}}\) at every transmission slot.

5.1.4 Optimization Objective↩︎

The CMDP objective is to find a policy \(\pi: \mathcal{S} \to \mathcal{P}(\mathcal{A})\) that maximizes the expected cumulative reward \(J_{\mathrm{R}}\) while satisfying the average cost constraints \(J_{\mathrm{co}}\) and \(J_{\mathrm{so}}\) [@book99_cmdp]: \[\begin{align} (\mathrm{P2}) \quad \underset{\pi}{\text{maximize}} \quad & J_{\mathrm{R}}(\pi) = \mathbb{E}_{\pi} \left[ \sum_{n=1}^{N} r_n \right] \tag{24} \\* \text{subject to} \quad & J_{\mathrm{co}}(\pi) = \mathbb{E}_{\pi} \left[ \frac{1}{N} \sum_{n=1}^{N} c_{\mathrm{co},n} \right] \leq \epsilon_{\mathrm{co}}, \tag{25} \\* & J_{\mathrm{so}}(\pi) = \mathbb{E}_{\pi} \left[ \frac{1}{N} \sum_{n=1}^{N} c_{\mathrm{so},n} \right] \leq \epsilon_{\mathrm{so}}, \tag{26} \end{align}\] where \(\mathbb{E}_{\pi}\) denotes the expectation over the trajectories induced by the policy \(\pi\). As discussed in Section 4.2, satisfying 25 and 26 also ensures the original average outage constraints 13 and 14 . The formulation is undiscounted, i.e., \(\gamma_{\mathrm{d}} = 1\). The finite horizon \(N\) bounds the cumulative reward and cost, and the undiscounted objective 24 matches the original problem 11 exactly.

5.2 PD-SAC Algorithm↩︎

We solve the CMDP using Lagrangian relaxation combined with the SAC framework [@achiam17_cpo; @19_RCPO; @paternain23_safe; @yang21_wcsac; @haarnoja18_sac; @haarnoja18b_sac].

5.2.1 Lagrangian Relaxation↩︎

The constrained problem 24 –26 is transformed into a min-max problem via Lagrangian duality \[\begin{align} \label{eq:lagrangian} \min_{\boldsymbol{\lambda} \succeq 0} \max_{\pi} \mathcal{L}(\pi, \boldsymbol{\lambda}), \end{align}\tag{27}\] where the Lagrangian is \[\begin{align} \mathcal{L}(\pi, \boldsymbol{\lambda}) = J_{\mathrm{R}}(\pi) - \sum_{i \in \{\mathrm{co}, \mathrm{so}\}} \lambda_{i}(J_{i}(\pi) - \epsilon_{i}), \end{align}\] and \(\boldsymbol{\lambda} = [\lambda_{\mathrm{co}}, \lambda_{\mathrm{so}}]^{\mathrm{T}} \succeq 0\) is the vector of Lagrange multipliers for the average outage constraints.

5.2.2 Maximum Entropy Framework↩︎

Following the SAC framework [@haarnoja18_sac], the objective is augmented with entropy regularization to enhance exploration. Since the secrecy rate is non-concave in the beamforming vector, this regularization promotes diverse action selection and mitigates convergence to local optima. The stochastic policy \(\pi_{\boldsymbol{\theta}}\) is parameterized by the actor neural network weights \(\boldsymbol{\theta}\); hence, the maximization with respect to \(\pi\) in 27 is performed by optimizing the actor-network parameters \(\boldsymbol{\theta}\). The entropy-regularized objective becomes \[\begin{align} \label{eq:entropy95obj} J_{\mathrm{ent}}(\boldsymbol{\theta}) = \mathbb{E}_{\pi_{\boldsymbol{\theta}}} \left[ \sum_{n=1}^{N} \left( r_n - \frac{1}{N} \sum_{i \in \{\mathrm{co}, \mathrm{so}\}} \lambda_i c_{i,n} + \alpha \mathcal{H}(\pi_{\boldsymbol{\theta}}(\cdot|s_n)) \right) \right], \end{align}\tag{28}\] where \(\alpha > 0\) is the temperature parameter scaling the entropy term, with larger \(\alpha\) encouraging exploration and smaller \(\alpha\) encouraging exploitation, and \(\mathcal{H}(\pi_{\boldsymbol{\theta}}(\cdot|s_n)) = -\mathbb{E}_{\tilde{q}_n \sim \pi_{\boldsymbol{\theta}}(\cdot|s_n)}[\log \pi_{\boldsymbol{\theta}}(u_n, v_n\,|\,s_n)]\) is the entropy of the policy, which quantifies the randomness of the sampled actions. Throughout, \(\mathbb{E}_{x \sim p}[\cdot]\) denotes the expectation with respect to \(x\) sampled from the distribution \(p\).

5.2.3 Network Architecture↩︎

As summarized in Table 1, the proposed SAC agent employs the following neural networks.

• Actor \(\pi_{\boldsymbol{\theta}}(\tilde{q}_n|s_n)\): A neural network that outputs the mean \(\boldsymbol{\mu}_{\boldsymbol{\theta}}(s_n)\) and log standard deviation \(\log \boldsymbol{\sigma}_{\boldsymbol{\theta}}(s_n)\) of a \((2M{+}1)\)-dimensional diagonal Gaussian over the latent \(\tilde{q}_n\), which is deterministically mapped to the squashed action \((u_n, v_n)\) via the tanh and sigmoid transformations in 21 and 22 . Latent actions are sampled via the reparameterization trick [@haarnoja18_sac] as \(\tilde{q}_n = \boldsymbol{\mu}_{\boldsymbol{\theta}}(s_n) + \boldsymbol{\sigma}_{\boldsymbol{\theta}}(s_n) \odot \boldsymbol{\xi}\), \(\boldsymbol{\xi} \sim \mathcal{N}(\mathbf{0}, \mathbf{I})\), then mapped to \(\mathbf{w}_n\) via 21 –23 .

• Reward Critics \(Q^{\mathrm{R}}_b(s_n,\mathbf{w}_n)\), \(b \in \{1,2\}\): Two networks, parameterized by \(\boldsymbol{\theta}^{\mathrm{R}}_b\), whose minimum is taken as a pessimistic estimate that mitigates the overestimation bias of the value function [@fujimoto18_td3] for the policy update.

• Cost Critics \(\mathbf{Q}^{\mathrm{C}}_b(s_n,\mathbf{w}_n) \in \mathbb{R}^2\), \(b \in \{1,2\}\): Two networks, parameterized by \(\boldsymbol{\theta}^{\mathrm{C}}_b\), each with two output heads; the components \([\mathbf{Q}^{\mathrm{C}}_b]_{\mathrm{co}}\) and \([\mathbf{Q}^{\mathrm{C}}_b]_{\mathrm{so}}\) estimate the connection and secrecy outage costs, respectively. With \(\gamma_{\mathrm{d}} = 1\), \([\mathbf{Q}^{\mathrm{C}}_b]_i = \frac{1}{N}\sum_{n'=n}^{N}\mathbb{E}_{\pi_{\boldsymbol{\theta}}}[c_{i,n'} \mid s_n, \mathbf{w}_n]\), \(i\in\{\mathrm{co},\mathrm{so}\}\), which matches the scale of the average constraints \(J_i\) in 25 and 26 . The element-wise maximum is then taken as a pessimistic cost estimate in the policy update 31 . In contrast to the minimum used for the reward critics, this maximum prevents underestimation of the costs [@yang21_wcsac].

5.2.4 Training Updates↩︎

The complete training procedure is summarized in Algorithm 4, where the specific network updates are performed as follows:

• Critic Update: The reward critics are trained by minimizing the mean-squared error between their output and the regression target: \[\begin{align} \label{eq:critic95loss} L(\boldsymbol{\theta}^{\mathrm{R}}_b) = \mathbb{E}_{(s_n,\mathbf{w}_n,r_n,s_{n+1}) \sim \mathcal{R}} \left[ \left( Q^{\mathrm{R}}_b(s_n,\mathbf{w}_n) - y_{\mathrm{R}} \right)^2 \right], \end{align}\tag{29}\] where the expectation is taken over transition samples drawn from the replay buffer \(\mathcal{R}\), and the target is \[\begin{align} &y_{\mathrm{R}} = r_n + (1\!-\!d) \,\nonumber\\ &\times \mathbb{E}_{\tilde{q}_{n+1} \sim \pi_{\boldsymbol{\theta}}(\cdot\,|\,s_{n+1})} \!\bigl[ \min_{b\in\{1,2\}} \bar{Q}^{\mathrm{R}}_b(s_{n+1},f(\tilde{q}_{n+1})) \nonumber\\ &\qquad\quad -\, \alpha \log \pi_{\boldsymbol{\theta}}(u_{n+1}, v_{n+1}\,|\,s_{n+1}) \bigr]. \end{align}\] Here, \(f(\tilde{q}) = \mathbf{w}\) denotes the deterministic mapping from the latent action to the transmit beamforming vector via 21 –23 , \(d \in \{0,1\}\) is the episode-termination indicator, and \(\bar{\boldsymbol{\theta}}^{\mathrm{R}}_b\) denotes the target parameters maintained as an exponential moving average (EMA) of the online parameters \(\boldsymbol{\theta}^{\mathrm{R}}_b\). The cost critic is updated similarly by minimizing the mean-squared error between its output and the regression target as \[\begin{align} \label{eq:cost95target} \mathbf{y}_{\mathrm{C}} = \frac{1}{N}\mathbf{c}_n + (1 - d) \, \max_{b\in\{1,2\}} \bar{\mathbf{Q}}^{\mathrm{C}}_b(s_{n+1}, f(\tilde{q}_{n+1})), \end{align}\tag{30}\] where \(\bar{\boldsymbol{\theta}}^{\mathrm{C}}_b\) denotes the corresponding target parameters for the cost critic, and the element-wise maximum over the two cost critics acts as a pessimistic target to avoid underestimating the costs.

  

• Actor Update: The policy is updated by minimizing the following loss, which corresponds to maximizing the entropy-regularized objective \(J_{\mathrm{ent}}\) in 28 with the cumulative reward and cost terms estimated by the critics \[\begin{align} \label{eq:actor95loss} &L(\boldsymbol{\theta}) = \mathbb{E}_{s_n \sim \mathcal{R},\, \tilde{q}_n \sim \pi_{\boldsymbol{\theta}}(\cdot\,|\,s_n)} \Big[ \alpha \log \pi_{\boldsymbol{\theta}}(u_n, v_n\,|\,s_n) \nonumber \\ & - \min_{b\in\{1,2\}} Q^{\mathrm{R}}_b(s_n, f(\tilde{q}_n)) + \boldsymbol{\lambda}^{\mathrm{T}} \max_{b\in\{1,2\}} \mathbf{Q}^{\mathrm{C}}_b(s_n, f(\tilde{q}_n)) \Big], \end{align}\tag{31}\] where the element-wise maximum over the two cost critics matches the pessimistic cost target in 30 . Following the standard SAC formulation [@haarnoja18_sac], the action log-density is \[\begin{align} \label{eq:log95prob} &\log \pi_{\boldsymbol{\theta}}(u_n, v_n\,|\,s_n) \!=\! \sum_{i=1}^{2M} \Bigl[\log p(\tilde{q}_{n,i}^{\mathrm{dir}})\! -\! \log\bigl(1\!-\!\tanh^2(\tilde{q}_{n,i}^{\mathrm{dir}})\bigr)\Bigr] \nonumber\\ &+ \log p(\tilde{q}_n^{\mathrm{pow}}) - \log \Bigl[P_{\mathrm{max}}\, \sigma(\tilde{q}_n^{\mathrm{pow}})\bigl(1-\sigma(\tilde{q}_n^{\mathrm{pow}})\bigr)\Bigr], \end{align}\tag{32}\] where \(p(\cdot)\) is the latent Gaussian density of the pre-activation samples, and the remaining terms are the log-Jacobian correction terms induced by the tanh and sigmoid transformations, which convert the latent density into the density of the transformed action \((u_n,v_n)\). The \(\ell_2\)-normalization \(\hat{\mathbf{d}}_n = u_n/\|u_n\|\) in 21 that produces the unit-norm direction is a deterministic post-processing step; its gradient with respect to \(\boldsymbol{\theta}\) propagates through the actor loss 31 via the chain rule, as \(f\) enters both the reward and cost critic terms.⁵ The temperature \(\alpha\) is adjusted via 35 .

• Dual Update: The Lagrange multipliers are updated via gradient descent on the dual loss. To ensure non-negativity, i.e., \(\lambda_i \geq 0\), we adopt a log-space parameterization \(\lambda_i = \exp(\tilde{\lambda}_i)\) with \(\tilde{\lambda}_i \in \mathbb{R}\), stacked as \(\tilde{\boldsymbol{\lambda}} = [\tilde{\lambda}_{\mathrm{co}}, \tilde{\lambda}_{\mathrm{so}}]^{\mathrm{T}}\). The log-space parameters are initialized to a small value \(\tilde{\lambda}_0 < 0\), so that the multipliers begin near zero and the policy first learns basic beamforming before constraint pressure is gradually applied.

  To reduce the bias from buffer samples drawn under earlier policies, we maintain an on-policy EMA of the per-step costs as \[\begin{align} \label{eq:cost95ema} \hat{c}_i \leftarrow \chi\, \hat{c}_i + (1 - \chi)\, \bar{c}_{i,n}, \quad i \in \{\mathrm{co}, \mathrm{so}\} \end{align}\tag{33}\] where \(\bar{c}_{i,n}\) is the average cost across parallel environments at step \(n\), and \(\chi \in (0,1)\) is the decay factor. The dual loss is then \[\begin{align} L_{\mathrm{dual}}(\tilde{\boldsymbol{\lambda}}) = -\sum_{i \in \{\mathrm{co}, \mathrm{so}\}} \lambda_i \left(\hat{c}_i - \epsilon_i\right), \label{eq:dual95update95co} \end{align}\tag{34}\] which is minimized with respect to \(\tilde{\boldsymbol{\lambda}}\) via gradient descent with learning rate \(\eta_{\lambda}\).

• Temperature Update: The temperature \(\alpha\) is adjusted to maintain a target entropy \(\bar{\mathcal{H}} = -\dim(\tilde{q}_n) = -(2M+1)\), as is conventional for continuous action spaces [@haarnoja18b_sac], by minimizing \[\begin{align} \label{eq:temp95loss} L(\alpha) = \mathbb{E}_{s_n \sim \mathcal{R},\,\tilde{q}_n \sim \pi_{\boldsymbol{\theta}}(\cdot\,|\,s_n)} \Big[ &-\alpha \big( \log \pi_{\boldsymbol{\theta}}(u_n, v_n\,|\,s_n) + \bar{\mathcal{H}} \big) \Big]. \end{align}\tag{35}\] When the policy entropy falls below \(\bar{\mathcal{H}}\), \(\alpha\) increases to encourage more stochastic actions, and vice versa.

5.3 Computational Complexity Analysis↩︎

Let \(D_{\tilde{q}} = 2M+1\) denote the latent action dimension, and \(H\) the hidden layer size. For the proposed algorithm, during inference, the actor network requires a single forward pass with complexity \(\mathcal{O}(D_\mathrm{s} H + H^2 + H D_{\tilde{q}})\). Since \(D_{\tilde{q}} \ll H\) in practice, this simplifies to \(\mathcal{O}(D_\mathrm{s} H + H^2)\). During training, the dominant cost per update step arises from the forward and backward passes through the actor and critic networks, yielding a complexity of \(\mathcal{O}(|\mathcal{B}|((D_\mathrm{s} + 2M) H + H^2))\), where \(|\mathcal{B}|\) is the minibatch size. The dual variable updates use only the scalar EMA costs \(\hat{c}_i\) and require \(\mathcal{O}(1)\) operations. Thus, the overall per-step training complexity is \(\mathcal{O}(|\mathcal{B}|((D_\mathrm{s} + 2M)H + H^2))\).

For MRT, the beamforming vector is computed as \(\mathbf{w}_n = \sqrt{P_{\mathrm{max}}} \mathbf{a}_{\mathrm{s},n} / \|\mathbf{a}_{\mathrm{s},n}\|\), which requires \(\mathcal{O}(M)\) operations per time slot. For ZF, the beamforming vector is obtained by projecting the serving array response onto the null space of the eavesdropper array response matrix, with complexity \(\mathcal{O}(M^2 E + M E^2 + E^3)\) per slot. For the per-slot SCA benchmark, each transmission slot is solved by a multistart SCA with \(R_{\mathrm{SCA}}\) initializations, each refined over \(K_{\mathrm{SCA}}\) outer iterations; every outer iteration solves a convex second-order cone program of dimension \(2M\) by an interior-point method whose \(K_{\mathrm{IP}}\) iterations are each dominated by an \(\mathcal{O}(M^3)\) KKT factorization, yielding a per-slot cost of \(\mathcal{O}(K_{\mathrm{SCA}} R_{\mathrm{SCA}} K_{\mathrm{IP}} M^3)\).

While the proposed RL approach incurs an offline training stage, its online execution complexity is strictly \(\mathcal{O}(D_\mathrm{s} H + H^2)\), scaling linearly with the number of antennas \(M\). By contrast, ZF requires algebraic projections that scale as \(\mathcal{O}(M^2 E + M E^2 + E^3)\) and becomes infeasible as \(E\) approaches \(M\), while per-slot SCA scales as \(\mathcal{O}(K_{\mathrm{SCA}} R_{\mathrm{SCA}} K_{\mathrm{IP}} M^3)\) and is unsuitable for real-time execution. The trained RL policy outputs the beamformer via a single forward pass, well-suited for real-time deployment under high satellite mobility.

6 Simulation Results↩︎

The simulation parameters listed in Table 2 are used unless otherwise stated. A separate PD-SAC policy is trained for each \(E\) on the deployment scenario using 50 parallel environment instances.⁶ We compare the proposed PD-SAC against MRT, ZF, per-slot SCA, and PD-PPO. For the per-slot SCA benchmark, the non-convex problem at each transmission slot, which maximizes the secrecy rate evaluated with the geometry-derived average SNRs subject to the power constraint and the per-slot counterparts \(c_{\mathrm{co},n} \leq \epsilon_{\mathrm{co}}\) and \(c_{\mathrm{so},n} \leq \epsilon_{\mathrm{so}}\) of the average outage constraints, is solved by a multistart SCA, where the secrecy rate is iteratively linearized and each convex subproblem is solved by an interior-point method [@markswright78]. The resulting beamformers are evaluated in the same simulation environment as the learned policies. The on-policy alternative PD-PPO shares the same CMDP formulation, network architecture, and dual variable structure as PD-SAC [@schulman2017proximal], with PPO-specific hyperparameters tuned for stable convergence.

Several measures are adopted to stabilize the primal-dual training. The Lagrange multipliers are clamped to \([\lambda_{\min}, \lambda_{\max}]\) to keep the constraint penalty active while avoiding oscillatory updates. The dual variables are further frozen during an initial warmup of \(T_{\mathrm{warm}}\) training steps. For evaluation and deployment, the EMA-averaged policy \(\pi_{\bar{\boldsymbol{\theta}}}\) is used in place of the training policy to smooth short-term oscillations induced by the primal-dual dynamics. The corresponding hyperparameters are listed in Table 2.

Fig. 5 shows the training convergence of PD-SAC for \(E=3\). In Fig. 5 (a), PD-SAC surpasses the MRT baseline within the first few episodes and the ZF baseline after around \(50\) episodes, converging to approximately \(2.8\)  bps/Hz, within roughly \(7\%\) of the offline SCA benchmark. As shown in Figs. 5 (b) and 5 (c), PD-SAC rapidly drives the connection outage below the threshold, while the secrecy outage is steered toward \(\epsilon_{\mathrm{so}}=0.3\) as the dual variable actively enforces the constraint boundary. The primal-dual updates remove the need for manual penalty tuning.

Fig. 6 shows the per-slot floating-point operation (FLOP) count for each method across \(E \in \{1,\ldots,7\}\). The dominant FLOP counts are approximately given by \(8M\) for MRT, \(8M^2 E + 8M E^2 + 2 E^3 + 8M\) for ZF, \(2[D_\mathrm{s} H + H^2 + H \cdot 2(2M+1)]\) for PD-SAC, and \(R_{\mathrm{SCA}} K_{\mathrm{SCA}} K_{\mathrm{IP}} [8M^3 + \mathcal{O}(M^2 E)]\) for SCA, consistent with the complexity analysis in Section 5.3. Using \(R_{\mathrm{SCA}}=10\) restarts per slot and the solver-measured averages \(K_{\mathrm{SCA}}\approx 9.7\) and \(K_{\mathrm{IP}}\approx 12\), at \(E=3\) PD-SAC inference requires approximately \(2.4\times 10^5\) FLOPs whereas the per-slot SCA solve requires about \(4.2\times 10^7\) FLOPs, over two orders of magnitude larger.

Table 3 compares the average secrecy rate and outage probabilities for the five schemes with \(E=3\). The “Bound” uses the upper bounds from Lemma 2, and “Exact” uses the true incomplete gamma function. Among the deployable policies, PD-SAC attains an average secrecy rate of \(2.8\)  bps/Hz while satisfying both average outage constraints under the upper bound, namely \(\frac{1}{N}\sum P_n^{\mathrm{co}}=0.011 \leq \epsilon_{\mathrm{co}}\) and \(\frac{1}{N}\sum P_n^{\mathrm{so}}=0.232 \leq \epsilon_{\mathrm{so}}\). The corresponding exact values, \(0.006\) and \(0.222\), also lie below the thresholds. This corresponds to a \(177\%\) gain over MRT and a \(13\%\) gain over ZF. SCA attains \(3.01\) bps/Hz as an offline benchmark; the proposed PD-SAC closes this gap to within \(7\%\) through a single forward pass at inference, with much less complexity as shown in Fig. 6. PD-PPO attains \(2.18\) bps/Hz, \(22\%\) below PD-SAC. MRT achieves the lowest connection outage since it maximizes the serving link gain, but it ignores the eavesdroppers entirely. Conversely, ZF completely eliminates the secrecy outage by nulling the eavesdropper channels, but this comes at the cost of a high connection outage, as the null constraint significantly reduces the beamforming gain toward the serving satellite. The evaluations using the derived bounds demonstrate close agreement with the exact formulations, confirming that the bound-based approximation incurs almost no practical performance loss.

Fig. 7 illustrates the beam gain \(|\mathbf{h}_{k,n}^{\mathrm{H}} \mathbf{w}_{n}|^2\) as a function of the azimuth angle \(\tau_{k,n}\) and the zenith angle \(\zeta_{k,n}\) at time slot \(n=376\). As expected, MRT consistently directs its main lobe toward the serving satellite but provides no suppression toward the eavesdroppers, resulting in high eavesdropper gain. ZF places nulls at the eavesdropper directions, but reduces the gain toward the serving satellite. The proposed PD-SAC steers the main lobe toward the serving satellite while partially nulling the eavesdropper directions to satisfy the average secrecy outage constraint.

Fig. 8 shows the secrecy performance as the number of eavesdroppers \(E\) increases. The proposed PD-SAC performs comparably to ZF for \(E=\{1,2\}\), but surpasses it with an increasing margin as \(E\) grows, and remains close to the offline SCA benchmark across the entire range of \(E\). In contrast, ZF’s connection outage rises sharply as the null-space dimension shrinks with each additional eavesdropper, exceeding the \(0.3\) threshold for \(E \geq 6\), while MRT’s secrecy outage rises steeply from \(E \geq 2\) and saturates near one for \(E \geq 3\), since it does not actively suppress eavesdroppers. PD-PPO achieves a lower secrecy rate than PD-SAC for \(E \geq 3\), with a gap that becomes pronounced as the number of eavesdroppers increases. Consequently, PD-SAC, PD-PPO, and SCA satisfy both outage constraints across the entire range, whereas ZF becomes infeasible at \(E \geq 6\) due to connection outage, and MRT violates the secrecy outage constraint for \(E \geq 2\). In addition, as the number of eavesdroppers increases, the secrecy rate eventually saturates because, in the non-colluding scenario, the secrecy performance is mainly determined by the most dominant eavesdropper rather than by all eavesdroppers collectively. In the considered along-track deployment, the large satellite-to-ground distances make the link geometry the dominant factor. Hence, adding more satellites beyond the strongest geometric eavesdropping positions causes only marginal additional degradation.

7 Conclusions↩︎

This paper investigated secure uplink beamforming for LEO satellite networks against multiple satellite eavesdroppers. We derived exact outage probabilities under Nakagami-\(m\) fading and developed tractable upper-bound cost functions to manage their intractability. By formulating the secrecy rate maximization as a CMDP, we proposed a PD-SAC algorithm to optimize the precoder. Simulations demonstrated that the proposed algorithm outperforms MRT and ZF, particularly as the number of eavesdroppers increases, while satisfying the outage constraints. It closely approaches an offline per-slot SCA benchmark while outperforming an on-policy PD-PPO counterpart. Although the SCA benchmark attains a higher secrecy rate, its iterative per-slot optimization is unsuitable for real-time deployment. In contrast, the proposed policy provides a deterministic mapping from the geometry-derived channel characteristics to the beamforming vector, enabling practical secure communications on dynamic satellite links. Future work includes addressing colluding eavesdroppers, imperfect channel state information, and massive antenna arrays.

---

1. J. Seo is with Samsung Electronics, Suwon 16677, South Korea (e-mail:juhwanseo92@gmail.com).↩︎

2. H. Cho is with the Department of Electrical and Electronic Engineering, Inha University, Incheon 22212, South Korea (e-mail: hyesang@inha.ac.kr)↩︎

3. D.-H. Jung is with the School of Electronic Engineering, Soongsil University, Seoul 06978, South Korea (e-mail: dhjung@ssu.ac.kr).↩︎

4. In intrinsic rotations, successive rotations are conducted about the axes rotated by the last rotation matrix. The superscript prime (\(\prime\)) is added to indicate the new axes after an elemental rotation.↩︎

5. The reference implementation pre-scales \(u\) by \(\sqrt{P_{\mathrm{max}}}\) before \(\ell_2\)-normalization for numerical conditioning; since the normalization erases any positive scalar, \(\mathbf{w}_n\) is unchanged and the scaling contributes only a state-independent additive constant to \(\log\pi_{\boldsymbol{\theta}}\), which has no effect on the actor gradient.↩︎

6. While the instantaneous CSI of eavesdroppers is unavailable, their cardinality \(E\) is observable from the tracked LEO constellation, since adversarial satellites can be enumerated via ephemeris-based space situational awareness even when their channels remain uncertain.↩︎