## 4 Joint User Scheduling and Hybrid Beamforming Design for Maximizing EE
### Penalty Iteration Algorithm
Therefore, the problem \(Q_{8}\) can be converted into  
\[Q_{9}:\max_{\mathbf{\eta},\text{W},\text{SINR}_{l}^{\text{min}}, \zeta}\text{EE} =F\Big{(}\mathbf{\eta},\mathbf{W},\text{SINR}_{l}^{\text{min}},\zeta\Big{)}- \lambda_{3}\sum\limits_{l}^{L}(Tr(\mathbf{W}_{l})-\lambda_{\text{max}}(\mathbf{W}_{l}) ),\] (74) \[s.t.\ C_{2},\ C_{3},\ C_{5},C_{6},C_{8}\ in\ Q_{8}\] (75) \[\text{where}\quad F\Big{(}\mathbf{\eta},\mathbf{W},\text{SINR}_{l}^{\text {min}},\zeta\Big{)} = 2q\bigg{(}B\sum\limits_{l=1}^{L}Y_{l}(\text{SINR}_{l}^{\text{ min}})\bigg{)}^{\frac{1}{2}}\ -\ q^{2}\bigg{(}\sum\limits_{l=1}^{L}Tr(\mathbf{W}_{l})+P_{0} \bigg{)}\] \[+\lambda_{1}\sum\limits_{l=1}^{L}\text{P}(\eta_{k,l})^{te}-\sum \limits_{l}^{L}\lambda_{2}\bigg{(}\sum\limits_{k}^{|U_{l}|}\eta_{k,l}-U_{s} \bigg{)}^{2}\ \text{and}\ \lambda_{3}>0\ \text{is the penalty factor, which is generally larger enough to ensure that a smaller value of \(Tr(\mathbf{W}_{l})-\lambda_{\text{max}}(\mathbf{W}_{l})\) can be obtained.  
According to (74), the iterative calculation of the problem \(Q_{9}\) can maximize function \(F\Big{(}\mathbf{\eta},\mathbf{W},\text{SINR}_{l}^{\text{min}},\zeta\Big{)}\) and minimize function \(Tr(\mathbf{W}_{l})-\lambda_{\text{max}}(\mathbf{W}_{l})\). It should be noted that \(Tr(\mathbf{W}_{l})\) is an affine function and \(\lambda_{\text{max}}(\mathbf{W}_{l})\) is nonsmooth, which can result in the nonsmoothness of the objective function in the problem \(Q_{9}\). To handle the challenge, we replace \(\lambda_{\text{max}}(\mathbf{W}_{l})\) with its first-order Taylor expansion. The subgradient of \(\lambda_{\text{max}}(\mathbf{W}_{l})\) is \(\frac{\partial\lambda_{\text{max}}(\mathbf{W}_{l})}{\partial\mathbf{W}}=\mathbf{w}_{l, \text{max}}\mathbf{w}_{l,\text{max}}^{H}\) and its first-order Taylor expansion can be expressed as follows, which is reflected in line eight of Algorithm 1.  
\[\lambda_{\text{max}}\Big{(}\mathbf{W}_{l}^{(t)}\Big{)}\geq\lambda_{\text{max}} \Big{(}\mathbf{W}_{l}^{(t-1)}\Big{)}+\Big{\langle}\mathbf{w}_{l,\text{max}}\mathbf{w}_{l, \text{max}}^{H}\mathbf{W}_{l}^{(t)}-\mathbf{W}_{l}^{(t-1)}\Big{\rangle}, \tag{76}\]  
\[\mathbf{W}_{l}=\lambda_{\text{max}}(\mathbf{W}_{l})\mathbf{w}_{l,\text{max}}\mathbf{w}_{l, \text{max}}^{H}, \tag{70}\]0We substitute (76) into the objective function in the problem \(Q_{9}\) to replace \(\lambda_{\max}(W_{l})\), and the problem \(Q_{9}\) can be expressed as:  
\[\mathbf{W}_{l}=\lambda_{\text{max}}(\mathbf{W}_{l})\mathbf{w}_{l,\text{max}}\mathbf{w}_{l, \text{max}}^{H}, \tag{70}\]1 \[\mathbf{W}_{l}=\lambda_{\text{max}}(\mathbf{W}_{l})\mathbf{w}_{l,\text{max}}\mathbf{w}_{l, \text{max}}^{H}, \tag{70}\]2  
In conclusion, the robust joint user scheduling and hybrid beamforming design algorithm for the massive MIMO LEO satellite multigroup multicast communication system is shown in Algorithm 1. ```
Input: CCCP algorithm iteration index \(k\), thresholds \(\epsilon_{1}\), penalty iteration algorithm iteration index \(m\), thresholds \(\epsilon_{2}\), penalty factor \(\lambda_{1},\lambda_{2},\lambda_{3}\). 1. Initial: \(\left(\boldsymbol{\eta},\boldsymbol{W},\text{SINR}_{l}^{\min},\zeta\right)^{( k-0)},q^{(k-0)}\). 2. while \(\left|\text{EE}^{(k)}-\text{EE}^{(k-1)}\right|\geq\epsilon_{1}\)
3. Convexification step by (57), (59), (61). 4. Calculation \(q^{(k)}\), substitute \(q^{(k)}\) into (77). 5. Optimization step. 6. Let \(\left(\boldsymbol{\eta},\boldsymbol{W},\text{SINR}_{l}^{\min},\zeta\right)^{( m=0)}=\left(\boldsymbol{\eta},\boldsymbol{W},\text{SINR}_{l}^{\min},\zeta \right)^{(k=0)}\). 7. while \(\left|\left\{Tr\Big{(}\boldsymbol{W}_{l}^{(m)}\Big{)}-\lambda_{\max}\Big{(} \boldsymbol{W}_{l}^{(m)}\Big{)}\right\}_{l=1}^{L}\right|\geq\epsilon_{2}\)
8. Calculate the maximum eigenvalue \(\lambda_{\max}(\boldsymbol{W}_{l})\) of \(\boldsymbol{W}_{l}^{(m)}\) and the corresponding eigenvector \(\boldsymbol{w}_{l,\max}^{(m)}\). 9.