# Modeling and Analysis of Staged Constellation Deployment from a Single-Unit System
## 3 Methodology  
As discussed in Section 2, the parameters of interest, \(\theta_{k,\ i_{f}}\) (from \(k=1\)_to_\(k_{max}-1\)) and \(\left(\Delta t_{a},\ \Delta t_{b}\right)_{k_{max}}\), increase in complexity with each deployment stage. However, this complexity can be mitigated through assumptions on satellite uniformity. Specifically, each disengagement stage reduces the satellite mass and surface area by half, while the drag coefficient remains constant. As a result, the maximum and minimum ballistic coefficients \(\beta_{Max}\) and \(\beta_{Min}\) remain unchanged across the early stages, ensuring that the angular acceleration \(\bar{\theta}\) remains constant during those phases. Given this consistency, satellites with identical orbital parameters will phase at similar rates. However, due to the cumulative effects of disengagement \(\Delta v\), satellites in later stages begin with greater relative separation rates. Therefore, initiating disengagements early is favorable, provided sufficient phasing time is allowed to prevent collisions. This assumption significantly simplifies the deployment analysis, reducing it to identifying the values of \(\theta_{\left(k_{max}-1\right),\ i_{f}}\) and timing parameters \(\left(\Delta t_{a},\ \Delta t_{b}\right)_{k_{max}}\) that yield the shortest total phasing time. The values of \(\left(\Delta t_{a},\Delta t_{b}\right)_{k_{max}}\) are determined analytically for each satellite-reference pair using Equations (3) and (4). Since the target constellation configuration predefines the final angular separations and the control authority \(\bar{\theta}\) is updated throughout the simulation, the unknowns reside in the initial conditions at stage \(k_{max}\), which are derived from the preceding stage. This is achieved by generating a set of candidate insertion points for each satellite at the penultimate stage. To resolve this, the phasing problem is reformulated as identifying the set of insertion angles \(\theta_{\left(k_{max}-1\right),\ i_{f}}\) that minimize the total phasing time. Each satellite is phased towards discrete angular targets, \(\theta\in\left\{0,\ \Delta 0,\ 2\Delta 0,\ \ldots,\ 360\right\}\), followed by the application of a disengagement \(\Delta v\). This process is summarized in Algorithm 1, which outputs the resulting initial relative motions and the associated timing values \(t_{\left(k_{max}-1\right),\ i}\). Subsequently, Algorithm 2 computes the phasing durations \(\left(\Delta t_{a},\Delta t_{b}\right)\) for each candidate insertion point, which formed the output of Algorithm 1, using the analytical solution. The combination that yields the shortest sum of \(\Delta t_{a},\ \Delta t_{b},\) and \(t_{\left(k_{max}-1\right),\ i}\) is identified as the most favorable. This simulation process is repeated for all \(n-1\) satellites in the final deployment stage. The longest of these values defines the final phasing time for the constellation, as it represents the slowest-converging satellite relative to the reference. This value is then used to determine the final phasing segment for the last satellite-reference pair. For instance, for a system of eight satellites and \(k_{max}=3\), the output will consist of optimal values for \(i\), \(1\), \(t_{2}\), \(2\), and \(t_{2}\), \(3\), calculated for satellites \(S_{2,1}\), \(S_{2,2}\), and \(S_{2,3}\), covering the branches from Section 2.4. The largest of these values defines the reference timing for the final stage. Throughout this process, \(\Delta v\) is assumed constant across all disengagements and is scaled proportionally to satellite mass. To determine the most effective \(\Delta v\) magnitude, Algorithms 1 and 2 are executed over a predefined range of candidate values. The value yielding the shortest total phasing time among the tested configurations is selected. Top-level flow charts of Algorithms 1 and 2 are illustrated in Figures 4 and 5, respectively.