Designing Active Tether-Net Systems for Space Debris Capture with Graph-Learning-Aided Mixed-Combinatorial Optimization

# Designing Active Tether-Net Systems for Space Debris Capture with Graph-Learning-Aided Mixed-Combinatorial Optimization
Source: [https://arxiv.org/html/2605.29021](https://arxiv.org/html/2605.29021)
Feng Liu1Achira Boonrath1,∗∗Equal contributionGishnu Madhu2,∗Eleonora M\. Botta32M\.S\. Student, Department of Computer Science and Engineering3Associate Professor, Department of Mechanical and Aerospace Engineering, AIAA memberSouma Chowdhury44Professor, Department of Mechanical and Aerospace Engineering; Professor \(adjunct\), Department of Computer Science and Engineering, AIAA Senior Member, Corr\. authorThis work is accepted to be presented in the 2026 AIAA Aviation Forum\.

###### Abstract

Active tether\-net systems offer a promising solution for capturing large non\-cooperative targets with uncertain dynamics, including space debris\. This system involves deploying a flexible net from a chaser spacecraft and using maneuverable units \(MUs\) to manipulate the trajectory of the net\. However, concurrent systematic explorations of design and control choices of the tether\-net system to understand its full potential remain limited, partly due to the complex, constrained, nonlinear optimization problem that it presents – one that involves a mixture of continuous, integer and categorical variables, with the latter two arising from net connectivity and component choices, respectively\. Classical binary encoding approaches from the combinatorial optimization paradigm do not lend themselves well to solving such Mixed Combinatorial Nonlinear Programmings \(MCNLPs\) problems in the engineering design domain due to the presence of highly nonlinear, often multimodal, functions as is the case here\. On the other hand, integer coding approaches introduce spurious relations among combinations\. Given the graph\-structured characteristics of the combinatorial space, this paper adopts and extends a new graph\-learning\-aided optimization approach to solve this MCNLP problem\. Here, a Graph Neural Network \(GNN\) is trained to score \(as output\) and thereof recommend candidate combinations represented as nodes in a graph, with the continuous variable vector portion of a candidate design given as input\. As a result, the MCNLP optimization reduces to an NLP, which can be solved using standard solvers\. While this reduction approach is agnostic to the choice of the NLP solver, here a state\-of\-the\-art Particle Swarm Optimization \(PSO\) algorithm with gradient\-based fine\-tuning is used as the solver\. Demonstrated on the problem of concurrently designing the morphology of the net, choice of mass and thrusters in the MUs and aiming points used by the controller of the tether\-net system, the GNN\-based recommender is shown to provide significantly faster convergence to similar optimal solutions, compared to direct solution of the MCNLP problem\.

## IIntroduction

Space debris poses an increasing risk to space missions in Earth’s orbit\[[1](https://arxiv.org/html/2605.29021#bib.bib141)\]\. Approaches to solve this problem fall underActive Debris Removal\(ADR\), which can be achieved with various methods, among which tether\-net systems are particularly promising for their ability to capture uncooperative and tumbling large\-scale debris while maintaining the chaser satellite at a safe distance\[[25](https://arxiv.org/html/2605.29021#bib.bib283),[2](https://arxiv.org/html/2605.29021#bib.bib149),[12](https://arxiv.org/html/2605.29021#bib.bib117),[13](https://arxiv.org/html/2605.29021#bib.bib258)\]\. For robustness of space debris capture missions \(e\.g\., improve the capture success rate\), the use of an active and robotic tether\-net system has been proposed in more recent years\[[20](https://arxiv.org/html/2605.29021#bib.bib263),[37](https://arxiv.org/html/2605.29021#bib.bib108),[17](https://arxiv.org/html/2605.29021#bib.bib170),[3](https://arxiv.org/html/2605.29021#bib.bib99),[38](https://arxiv.org/html/2605.29021#bib.bib254)\]\. Such systems have multipleManeuverable Unitsattached to the corners or edges of the net, whose thrusters can transport the net towards the debris\. After the debris is captured, the chaser initiates a deorbiting or relocation process to transport it to a graveyard orbit\. Methods for controlling theMUsthat have been considered in the literature include closed\-loop and open\-loop approaches\[[16](https://arxiv.org/html/2605.29021#bib.bib69),[20](https://arxiv.org/html/2605.29021#bib.bib263),[38](https://arxiv.org/html/2605.29021#bib.bib254),[3](https://arxiv.org/html/2605.29021#bib.bib99)\]\. For example, Meng et al\. proposed a dual\-loop control scheme for the deployment of a robotic tether\-net system\[[20](https://arxiv.org/html/2605.29021#bib.bib263)\]\. Zhu et al\. developed a feedback controller that would allow a robotic tether\-net system to capture multiple small pieces of space debris located in close proximity\[[38](https://arxiv.org/html/2605.29021#bib.bib254)\]\. Boonrath et al\. proposed a robotic net system with a proportional\-integral\-derivative thrust controller that is guided by a reinforcement learning policy for debris\-capture tasks\[[3](https://arxiv.org/html/2605.29021#bib.bib99)\]\. To enable such autonomous control processes, it is important to carefully consider the structure, components \(that support autonomy, e\.g\., sensors/actuators\) and control choices, instead of simply adopting traditional \(baseline\) physical form for the net and corner maneuverable units \(MUs\)\. This paper presents a computational framework that allows efficiently searching over the complex mixed continuous\-combinatorial design spaces presented by such joint consideration of structure, components and control choices; secondly, it provides a proof of concept for how graph neural nets can be used to systematically \(as opposed to ad hoc\) and meaningfully \(as opposed to indexing\) abstract the space of combinatorial options; this could be useful for automated guided design of other autonomous spacecraft systems as well\.

Simultaneously optimizing the physical design and the controlled trajectory of such systems is a challenging task\. First, the hardware selections for the tether\-net systems may inherently involve fixed\-parameter groupings rather than independent variables only\. In the example of selecting thrusters, key parameters such as the thrust saturation limit, specific impulse, and dry mass are predefined by the manufacturers and are fixed for each thruster type\. In recent years, more methods have been proposed to solve this type ofMixed Combinatorial Nonlinear Programming\(MCNLP\) problem in design or optimal control problems\. For example, to efficiently address the combinatorial aspects of the satellite layout problem, Xia et al\.\[[33](https://arxiv.org/html/2605.29021#bib.bib24)\]formulated it as a mixed\-integer program and used a bi\-level optimization approach that separates discrete and continuous variables\. Ringkamp et al\.\[[23](https://arxiv.org/html/2605.29021#bib.bib43)\]implemented the concept of time transformation, which converts discrete, integer\-controlled time intervals into continuous time functions, thereby smoothing out the integer\-control problem\. Though the above methods prove to be capable of handling theMCNLPproblems in some engineering problems, they rely on direct computing of the objective functions for every candidate combination, which is computationally expensive if the simulation involved in the objective functions is computationally expensive; furthermore, the relations of the candidates in the combinatorial space in those methods are not depicted illustratively\.

On the other hand, graphs have been proven to be useful to abstract the combinatorial space in the field of classical and contemporary planning problems \(and some design problems\) that are formed as combinatorial optimization problems, and the approximate solutions can be learned usingGraph Neural Networks\[[32](https://arxiv.org/html/2605.29021#bib.bib70),[28](https://arxiv.org/html/2605.29021#bib.bib3)\]\. Such problems include vehicle routing\[[15](https://arxiv.org/html/2605.29021#bib.bib52)\], multi\-robot task allocation\[[22](https://arxiv.org/html/2605.29021#bib.bib53),[21](https://arxiv.org/html/2605.29021#bib.bib2)\], and molecular design\[[35](https://arxiv.org/html/2605.29021#bib.bib18)\]\. In our previous research,GNNshave successfully been used to optimize both the physical design and the controlled trajectory of the tether\-net system in a challenging environment\[[19](https://arxiv.org/html/2605.29021#bib.bib72)\]: agraph convolutional network\(GCN\)\[[36](https://arxiv.org/html/2605.29021#bib.bib55)\]was utilized to efficiently recommend the best combinatorial components \(e\.g\., selection of the thruster\) for a given vector of non\-combinatorial variables \(e\.g\., control parameters\)\. MoreGNNstructures remain to be explored for the tether\-net co\-design problem, such asgraph attention networks\[[29](https://arxiv.org/html/2605.29021#bib.bib56)\]and capsule networks\[[30](https://arxiv.org/html/2605.29021#bib.bib57)\]\.

Motivated by the above related research, we propose to use an updatedGNN\-aidedParticle Swarm Optimization\(PSO\) framework from our parallel research to minimize the fuel expended in the deployment maneuver by the tether\-net system while ensuring successful capture; this is done by optimizing its physical design and the employed control\. Compared to our previous research\[[19](https://arxiv.org/html/2605.29021#bib.bib72)\], the main differences are as follows:1\)the structure of theGNNhas been switched from aGNNto an encoder\-decoder architecture comprising aGNNencoder and an MLP decoder\. This architecture is referred to as Edge Flow Graph Network\.2\)The graph abstraction has been enhanced to represent a directed graph – with edges encoding pairwise objective differences between nodes –, which is the ground truth vector to learn\. While in earlier work node\-level scores were learned, the current formulation focuses on edge\-level signals\. This abstraction enables the effective use of edge relationships and allows the direct encoding of the difference between fitness scores, which inherently represent the problem’s mixed combinatorial nature, onto the graph\.3\)The constraints of the problem have been formulated into penalty functions\.4\)The design space of the tether\-net system has been increased to include more continuous and combinatorial variables\.5\)The final \(i\.e\., end\-of\-deployment\) desired velocity of eachMUis now determined by the optimization framework, which, while increasing the complexity of the control framework, also widens the range of possibleMUmaneuvers\.

The key contributions of the paper can be summarized as follows:1\) SpecializedGNNLoss Function:We design a specialized loss function to train theGNNand a corresponding metric \(e\.g\., sign accuracy\) to evaluate the performance of the trainedGNN, enabling robust recommendation of combinatorial choices given the continuous portion of a candidate design\.2\) Optimization Framework:We develop a framework for joint structure\-component\-control optimization with this embeddedGNNto improve the maneuvering efficiency of the active tether\-net system, specifically for minimizing the fuel cost of the maneuver, while preserving successful capture of the debris\. In this work, we usePSOas the optimizer, with its results post\-processed by sequential quadratic programming for optimality assurances; however, the proposed framework is also compatible with other population\-based or gradient\-based algorithms\.3\) Tether\-Net Simulation:We updated our active tether\-net system with aModel Predictive Control\(MPC\) system \(instead of the earlier PID system\) serving as the low\-level controller that tracks a minimum\-energy reference path to achieve the aiming points \(high level control decisions, that is optimized jointly with design\)\.

The subsequent sections of this paper are arranged as follows: Section[II](https://arxiv.org/html/2605.29021#S2)introduces the active tether\-net system, its components, and the controller for theMUs\. Section[III](https://arxiv.org/html/2605.29021#S3)introduces the optimization framework and formulations, along with the embeddedGNNmodel\. Section[IV](https://arxiv.org/html/2605.29021#S4)discusses the results of the optimization in the test scenario\. Section[V](https://arxiv.org/html/2605.29021#S5)then concludes the paper and discusses the future work for this research\.

## IITether\-Net Design Optimization

### II\-ASystem Components, Mission Procedure, and Evaluation of Capture Success

The net\-basedADRsystem of interest for this work includes a robotic net with a closing mechanism interlaced with its perimeter, a main tether connecting the net to a chaser vehicle, and an autonomous chaser vehicle\. The chaser is a rigid body that is assumed to be equipped with a winch, all the required sub\-systems \(including propulsion\), and a sensing system that enables estimating the relative location, velocity, and angular states of the target debris\. For the net itself, at the corners of the square\-shaped structure,MUsand the closing mechanism threads are attached; in simulation, the net structure is modeled as a collection of particles linked together via spring\-damper elements that cannot withstand compression\. Here, eachMUrepresents a nanosatellite\-sized spacecraft that can be controlled using small thrusters and inertial measurement units for navigation\. Figure[1\(a\)](https://arxiv.org/html/2605.29021#S2.F1.sf1)displays the components of the system and the reference frames used in this work, including the inertial reference frame𝒪=\{O,𝐢^,𝐣^,𝐤^\}\\mathcal\{O\}=\\\{O,\\mathbf\{\\hat\{i\}\},\\mathbf\{\\hat\{j\}\},\\mathbf\{\\hat\{k\}\}\\\}– fixed to the initial location of the net’s corner knot adjacent toMU1 – and the target body\-fixed frame𝒟=\{D,𝐝^x,𝐝^y,𝐝^z\}\\mathcal\{D\}=\\\{D,\\mathbf\{\\hat\{d\}\}\_\{x\},\\mathbf\{\\hat\{d\}\}\_\{y\},\\mathbf\{\\hat\{d\}\}\_\{z\}\\\}– fixed to the target’sCenter of Mass\(CoM\)\. The target debris to be captured – with its size and mass information listed in Table[I](https://arxiv.org/html/2605.29021#S2.T1)– is a rigid\-body in\-simulation that is modeled after the upper stage of the Zenit\-2 rocket\. Table[I](https://arxiv.org/html/2605.29021#S2.T1)also lists the chaser spacecraft’s properties\. The physical properties of the net \(except for the axial damping ratio of the net threads, set at 0\.106\) will be determined through the optimization process defined later in this work\.

TABLE I:Debris and Chaser Component SettingParameterValueUnitDebrisMass9000kgDiameter3\.9mLength11\.0mVolume125\.3m3Surface Area159\.9m2ChaserMass1600kgSide length1\.5mThe thrusters on theMUsare switched on at the beginning of the mission to launch the net from the chaser satellite and maneuver it towards the target\. After 25 seconds, the thrust is deactivated, and the net is closed by shortening and then locking the closing mechanism threads, interlaced around the net\. The number of locked segments of the closing mechanism,NLN\_\{\\text\{L\}\}, and theCapture Quality Index\(CQI\) are selected as metrics to determine whether a capture is successful or not\. TheCQIis defined as:

J=0\.1​\|Vn−VD\|VD\+0\.1​\|Sn−SD\|SD\+0\.8​\|qn\|LcJ=0\.1\\frac\{\|V\_\{n\}\-V\_\{D\}\|\}\{V\_\{D\}\}\+0\.1\\frac\{\|S\_\{n\}\-S\_\{D\}\|\}\{S\_\{D\}\}\+0\.8\\frac\{\|q\_\{n\}\|\}\{L\_\{c\}\}\(1\)In the equation, variables are defined as follows: theConvex Hull\(CH\) volume of the net isVnV\_\{n\}, the debris’s volume isVDV\_\{D\}, theCHsurface area of the net isSnS\_\{n\}, the surface area of the debris isSDS\_\{D\}, the distance between the debrisCoMand the netCoMisqnq\_\{n\}, and the minimum distance from the debris’sCoMto its surface isLcL\_\{c\}\. The value ofCQIquantifies the similarity between the shape of the debris and the net, as well as the distance between theirCoMs\. For this work, theCQIvalue at the end of the mission \(i\.e\., 10 s after net closing activation\) is indicated withICQI∗I^\{\*\}\_\{\\text\{CQI\}\}, and a capture simulation withICQI∗≤2\.5I^\{\*\}\_\{\\text\{CQI\}\}\\leq 2\.5andNL=12N\_\{\\text\{L\}\}=12\(i\.e\., the maximum number of locked segments\) is defined as successful\. In this work, all the simulations are executed on a Python\-based simulator, where Google JAX’sjax\.numpy\.arraydata structure is employed to enable accelerated numerical computation performance\[[16](https://arxiv.org/html/2605.29021#bib.bib69)\]\. Interested readers may refer to previous works by the authors for more information on the system’s in\-simulation modeling\[[5](https://arxiv.org/html/2605.29021#bib.bib110),[6](https://arxiv.org/html/2605.29021#bib.bib116),[16](https://arxiv.org/html/2605.29021#bib.bib69)\]\.

![Refer to caption](https://arxiv.org/html/2605.29021v1/x1.png)\(a\)Diagram of the tether\-net system\. The crosses are examples aiming points of theMUsfor illustrative purposes\.
![Refer to caption](https://arxiv.org/html/2605.29021v1/x2.png)\(b\)Illustrations of the selected optimization variables, withKclsK\_\{\\text\{cls\}\}at lower bound in eachNkN\_\{k\}case as the example\.

Figure 1:Illustrations of the tether\-net system
### II\-BManeuvering Control Process

In this work, we aim to employ the proposed optimization framework to determine the optimal controlled trajectory of theMUsand the physical design of the debris removal system\. Here, the control formulation for the net deployment is established, based onMUcontrol frameworks previously employed by Boonrath et al\.\[[3](https://arxiv.org/html/2605.29021#bib.bib99),[4](https://arxiv.org/html/2605.29021#bib.bib106)\]\. The control forces for the net deployment \(activating and deactivating at instances mentioned in Section[II](https://arxiv.org/html/2605.29021#S2)\) are determined through a path\-trackingMPCcontroller\.

For eachMU, the tracked reference path is based on aminimum\-energy solution for a referenceMUmodel without the net attachedto travel from a set initial state – that is, a set combination of position and velocity – to a desired final state within a selected maneuver timetft\_\{f\}, as computed using the methodologies detailed in\[[4](https://arxiv.org/html/2605.29021#bib.bib106)\]\. The initial states for the minimum\-energy optimal control problem correspond to the initial states of theMUs, while the final states are defined as follows\. First, an aiming point is established for eachMU, on a plane that is located 5 m behind the debris’sCoMand is parallel to the𝐢^\\hat\{\\mathbf\{i\}\}\-𝐣^\\hat\{\\mathbf\{j\}\}plane \(see Figure[1\(a\)](https://arxiv.org/html/2605.29021#S2.F1.sf1)\)\. Theii\-thMUaiming point defined as:

𝐫final,i=\(xnom,i\+Δ​xi\)​𝐢^\+\(ynom,i\+Δ​yi\)​𝐣^\\displaystyle\\mathbf\{r\}\_\{\\text\{final\},i\}=\(x\_\{\\text\{nom\},i\}\+\\Delta x\_\{i\}\)\\mathbf\{\\hat\{i\}\}\+\(y\_\{\\text\{nom\},i\}\+\\Delta y\_\{i\}\)\\mathbf\{\\hat\{j\}\}\+\(zD−5\)​𝐤^,\\displaystyle\+\(z\_\{\\text\{D\}\}\-5\)\\mathbf\{\\hat\{k\}\},\(2\)i=1,2,3,4\\displaystyle i=1,2,3,4wherexnom,ix\_\{\\text\{nom\},i\}andynom,iy\_\{\\text\{nom\},i\}are the nominal coordinates of theii\-thMUaiming point on the𝐢^\\hat\{\\mathbf\{i\}\}\-𝐣^\\hat\{\\mathbf\{j\}\}plane, which are defined to bexnom,i=xD\+12\.0​\(−1\)ix\_\{\\text\{nom\},i\}=x\_\{\\text\{D\}\}\+12\.0\(\-1\)^\{i\}andxnom,i=xD\+12\.0​\(−1\)⌊\(i\+1\)/2⌋x\_\{\\text\{nom\},i\}=x\_\{\\text\{D\}\}\+12\.0\(\-1\)^\{\\lfloor\(i\+1\)/2\\rfloor\}, whilexDx\_\{\\text\{D\}\},yDy\_\{\\text\{D\}\}, andzDz\_\{\\text\{D\}\}are the debris’sCoMcoordinates\. The variablesΔ​xi\\Delta x\_\{i\}andΔ​yi\\Delta y\_\{i\}denote coordinate offsets of the aiming point of theii\-thMU, to be optimized for by theGNN\-based optimization process proposed in this work\. Next, the final velocity for eachii\-thMU,𝐯final,i\\mathbf\{v\}\_\{\\text\{final\},i\}, is defined as:

𝐯final,i=\(‖𝐫final,i−𝐫0,i‖tf\+vi\)​𝐫final,i−𝐫0,i‖𝐫final,i−𝐫0,i‖\\mathbf\{v\}\_\{\\text\{final\},i\}=\\left\(\\frac\{\|\|\\mathbf\{r\}\_\{\\text\{final\},i\}\-\\mathbf\{r\}\_\{0,i\}\|\|\}\{t\_\{f\}\}\+v\_\{i\}\\right\)\\frac\{\\mathbf\{r\}\_\{\\text\{final\},i\}\-\\mathbf\{r\}\_\{0,i\}\}\{\|\|\\mathbf\{r\}\_\{\\text\{final\},i\}\-\\mathbf\{r\}\_\{0,i\}\|\|\}\(3\)where‖𝐫final,i−𝐫0,i‖tf\\frac\{\|\|\\mathbf\{r\}\_\{\\text\{final\},i\}\-\\mathbf\{r\}\_\{0,i\}\|\|\}\{t\_\{f\}\}is the nominal final speed for theMUsandviv\_\{i\}is a scalar quantity, selected via theGNN\-based optimization process, defining a discrepancy from the nominal final speed\. For eachii\-thMU, solving the optimal minimum\-energy control problem with the initial and final states defined above results in the states𝐱r,i​\(t\)=\[𝐫r,iT​\(t\),𝐯r,iT​\(t\)\]T\\mathbf\{x\}\_\{r,i\}\(t\)=\[\\mathbf\{r\}\_\{\{r,i\}\}^\{T\}\(t\),\\mathbf\{v\}\_\{\{r,i\}\}^\{T\}\(t\)\]^\{T\}defined for the entire reference path, consisting of the reference positions𝐫r,i​\(t\)\\mathbf\{r\}\_\{r,i\}\(t\)and the reference velocities𝐯r,i​\(t\)\\mathbf\{v\}\_\{r,i\}\(t\)\.

To track the reference path, the control force for eachii\-thMUis determined via aMPCcontroller\. To accomplish this, a referenceMUmodel without the net attached \(i\.e\., the same model used to compute the minimum\-energy reference path\) is employed\. This model and reference trajectory are converted from their original continuous\-time form to a discrete\-time form \(with a sampling frequency of 20 Hz\[[3](https://arxiv.org/html/2605.29021#bib.bib99)\]\) using thepython\-controllibrary, yielding the discrete\-timeAAandBBmatrices for the state and control, respectively\. With the discrete\-time model established, the mathematical formulation for theMPCis defined as follows for eachii\-thMU:

min\{𝐱k,i,𝐮k,i\}\\displaystyle\\min\_\{\\\{\\mathbf\{x\}\_\{k,i\},\\mathbf\{u\}\_\{k,i\}\\\}\}\(𝐱Nh,i−𝐱Nh,r,i\)T​QN​\(𝐱Nh,i−𝐱Nh,r,i\)\\displaystyle\(\\mathbf\{x\}\_\{N\_\{h\},i\}\-\\mathbf\{x\}\_\{N\_\{h\},r,i\}\)^\{T\}Q\_\{N\}\(\\mathbf\{x\}\_\{N\_\{h\},i\}\-\\mathbf\{x\}\_\{N\_\{h\},r,i\}\)\(4\)\+∑k=0Nh−1\[\(𝐱k,i−𝐱k,r,i\)TQ\(𝐱k,i−𝐱k,r,i\)\\displaystyle\\quad\+\\sum\_\{k=0\}^\{N\_\{h\}\-1\}\\Big\[\(\\mathbf\{x\}\_\{k,i\}\-\\mathbf\{x\}\_\{k,r,i\}\)^\{T\}Q\(\\mathbf\{x\}\_\{k,i\}\-\\mathbf\{x\}\_\{k,r,i\}\)\+𝐮k,iTR𝐮k,i\]\\displaystyle\\quad\\quad\+\\mathbf\{u\}\_\{k,i\}^\{T\}R\\mathbf\{u\}\_\{k,i\}\\Big\]subject to𝐱k\+1,i=A​𝐱k,i\+B​𝐮k,i,k=0,…,Nh−1\\displaystyle\\mathbf\{x\}\_\{k\+1,i\}=A\\mathbf\{x\}\_\{k,i\}\+B\\mathbf\{u\}\_\{k,i\},\\quad k=0,\\dots,N\_\{h\}\-1𝐮min,i≤𝐮k,i≤𝐮max,i,k=0,…,Nh−1\\displaystyle\\mathbf\{u\}\_\{\\min,i\}\\leq\\mathbf\{u\}\_\{k,i\}\\leq\\mathbf\{u\}\_\{\\max,i\},\\quad k=0,\\dots,N\_\{h\}\-1
where the subscriptkkindicates thekk\-th time step forMPC111It should be noted that the actual net capture simulation is still performed in continuous time, whereas the above discrete\-time formulation is only used for the MPC controller\.\. States of eachii\-thMUat eachkk\-th step,𝐱i,k=\[𝐫i,kT,𝐯i,kT\]T\\mathbf\{x\}\_\{\{i,k\}\}=\[\\mathbf\{r\}\_\{\{i,k\}\}^\{T\},\\mathbf\{v\}\_\{\{i,k\}\}^\{T\}\]^\{T\}, are unconstrained for this work, while the inputs \(i\.e\., the control force\)𝐅T,i,k=𝐮i,k=\[ux,i,k,uy,i,k,uz,i,k\]T\\mathbf\{F\}\_\{T,i,k\}=\\mathbf\{u\}\_\{i,k\}=\[u\_\{\{x,i,k\}\},u\_\{\{y,i,k\}\},u\_\{\{z,i,k\}\}\]^\{T\}are constrained as𝐮min,i=−13​\[FT,max,FT,max,FT,max\]T\\mathbf\{u\}\_\{\\min,i\}=\-\\frac\{1\}\{\\sqrt\{3\}\}\[F\_\{T,\\text\{max\}\},F\_\{T,\\text\{max\}\},F\_\{T,\\text\{max\}\}\]^\{T\}and𝐮max,i=13​\[FT,max,FT,max,FT,max\]T\\mathbf\{u\}\_\{\\max,i\}=\\frac\{1\}\{\\sqrt\{3\}\}\[F\_\{T,\\text\{max\}\},F\_\{T,\\text\{max\}\},F\_\{T,\\text\{max\}\}\]^\{T\}such that the total maximum thrust magnitude available for eachMUisFT,maxF\_\{T,\\text\{max\}\}\. The weight matrices are heuristically defined asQN=Q=Q\_\{N\}=Q=diag\(\[10\.0, 10\.0, 10\.0, 1\.0, 1\.0, 1\.0\]\) andR=R=diag\(\[0\.01, 0\.01, 0\.01\]\) with the finite\-horizon step size,NhN\_\{h\}, chosen forMPCbeing 10\. Using theosqppackage in Python, Eq\. \([4](https://arxiv.org/html/2605.29021#S2.E4)\) is solved iteratively throughout the net’s flight to the target at 20 Hz \(i\.e\., the thrust value applied to the continuous\-time, full\-net system is updated every 0\.05 s\)\. With theMPC\-computed thrust values, the total fuel consumed is determined as:

mprop=∑i=14∫0tf‖𝐅T,i​\(τ\)‖ge​Isp​𝑑τm\_\{\\text\{prop\}\}=\\sum\_\{i=1\}^\{4\}\\int\_\{0\}^\{t\_\{f\}\}\\frac\{\|\|\\mathbf\{F\}\_\{T,i\}\(\\tau\)\|\|\}\{g\_\{e\}I\_\{\\text\{sp\}\}\}d\\tau\(5\)wherege=9\.81g\_\{e\}=9\.81m/s2andIspI\_\{\\text\{sp\}\}is the specific impulse, whose value is dependent on the propulsion system to be chosen\.

### II\-CGNN\-aided Optimization Formulation

TABLE II:Tether\-Net Fuel Consumption Optimization Combinatorial and Integer VariablesContinuous VariablesCategorySub\-categoryVariableBoundControlMUControlΔ​xi,i=1,2,3,4\\Delta x\_\{i,i=1,2,3,4\}\[m\]\[\-5, 5\]Δ​yi,i=1,2,3,4\\Delta y\_\{i,i=1,2,3,4\}\[m\]\[\-5, 5\]vi,i=1,2,3,4v\_\{i,i=1,2,3,4\}\[m/s\]\[\-1, 4\]Physical DesignMUMassmMUm\_\{\\text\{MU\}\}\[kg\]\[2, 3\]\[0\.8pt/2pt\]Net Morphologyrthreadr\_\{\\text\{thread\}\}\[m\]\[0\.0005, 0\.0015\]rcornerr\_\{\\text\{corner\}\}\[m\]\[0\.0001, 0\.0015\]LnetL\_\{\\text\{net\}\}\[m\]\[19, 25\]LctL\_\{\\text\{ct\}\}\[m\]\[0\.5, 2\]

Combinatorial and Integer Variables \(From\[[19](https://arxiv.org/html/2605.29021#bib.bib72)\]\)CategorySub\-categoryVariable\#1\#2\#3\#4\#5CombinatorialThrusterFT,maxF\_\{T,\\text\{max\}\}\[N\]8\.93\.66\.15\.56\.0IspI\_\{\\text\{sp\}\}\[s\]60\.057\.0277\.0253\.0250\.0mTm\_\{T\}\[kg\]0\.370\.0230\.60\.480\.25\[0\.8pt/2pt\]Net MaterialEnE\_\{\\text\{n\}\}\[GPa\]70\.070\.5112\.4\-\-ρn\\rho\_\{\\text\{n\}\}\[kg/m3\]1390\.01440\.01440\.0\-\-\[0\.8pt/2pt\]Net ShapeNkN\_\{k\}91113\-\-KclsK\_\{\\text\{cls\}\}\-2 to 0\-2 to \+1\-3 to \+1\-\-

For tether\-net systems, the physical design of the net and the controlled trajectory of theMUssignificantly affect the probability of capture success\. However, most previous works\[[14](https://arxiv.org/html/2605.29021#bib.bib50),[34](https://arxiv.org/html/2605.29021#bib.bib105)\]that examined design and control variables used continuous or discrete integer values, with limited studies exploring physical designs involving combinatorial variables\. Considering combinatorial variables is equally important, as certain design choices inherently involve fixed parameter groupings rather than independent continuous values\.

The goal of the optimization task is to minimize the fuel cost, in terms of propellant mass consumedmpropm\_\{\\text\{prop\}\}, during the capture mission while accomplishing a successful capture\. The fuel cost is calculated based on the thrust of theMUs, which is controlled by theMPCcontroller, and the maximum thrust is specified by the manufacturer\. All the related optimization variables are listed in Table[II](https://arxiv.org/html/2605.29021#S2.T2)\.

Continuous Variables:Continuous variables include theMUcontrol variables, the mass of eachMU, and the net morphology parameters\. TheMUcontrol variables are defined by the aiming point position illustrated in Fig\.[1\(a\)](https://arxiv.org/html/2605.29021#S2.F1.sf1)and the final speed of theMU\(i\.e\., the speed attft\_\{f\}\)\. The net morphology parameters include the net thread radiusrnetr\_\{\\text\{net\}\}, the corner thread radiusrcornerr\_\{\\text\{corner\}\}, the length of net threadLnetL\_\{\\text\{net\}\}, and the length of the corner threadLctL\_\{\\text\{ct\}\}\. These parameters are depicted in the top\-left part of Fig\.[1\(b\)](https://arxiv.org/html/2605.29021#S2.F1.sf2)\.

Combinatorial and Integer Variables:Five thrusters were picked as the potential choices, and each thruster is defined by its maximum thrust magnitudeFT,maxF\_\{T,\\text\{max\}\}, its specific impulseIspI\_\{\\text\{sp\}\}, and its massmTm\_\{T\}; the considered thrusters include cold\-gas and combustion\-based thrusters in order to explore the design space of the tether\-net system\[[31](https://arxiv.org/html/2605.29021#bib.bib68),[10](https://arxiv.org/html/2605.29021#bib.bib65),[26](https://arxiv.org/html/2605.29021#bib.bib66)\]\. Three possible types of high\-strength fibers are picked for the threads of the net\[[7](https://arxiv.org/html/2605.29021#bib.bib96),[9](https://arxiv.org/html/2605.29021#bib.bib67)\], each defined by the material’s Young’s ModulusEnE\_\{\\text\{n\}\}, and densityρn\\rho\_\{\\text\{n\}\}; the Poisson ratio for all fibers is assumed to be the same\. The variableNkN\_\{k\}is an integer variable that indicates the number of knots on each side of the net, resulting in a net withNk2N\_\{k\}^\{2\}knots in total\. The variableKclsK\_\{\\text\{cls\}\}represents the index of the node on a net edge where the closing mechanism is connected, as shown on the right side of Fig\.[1\(b\)](https://arxiv.org/html/2605.29021#S2.F1.sf2)\. A higher index value corresponds to a connection node closer to the central node of the edge222The set of connection nodes selected by default \(i\.e\., the connection nodes selected whenKcls=0K\_\{\\text\{cls\}\}=0\) is chosen according to the methodology presented in\[[6](https://arxiv.org/html/2605.29021#bib.bib116)\]\.\. Note that 1\) the upper and lower bounds ofKclsK\_\{\\text\{cls\}\}depend onNkN\_\{k\}\(see Table[II](https://arxiv.org/html/2605.29021#S2.T2)\), and 2\) the sameKclsK\_\{\\text\{cls\}\}value results in different offsets from each edge’s central node whenNkN\_\{k\}differs\. For allKclsK\_\{\\text\{cls\}\}, the connection node offsets are defined so the two closing nodes of each net edge are displaced equally away from the center of each edge \(i\.e\., they are symmetric\)\. In total, 180 unique combinations exist for thruster, net material, and net morphology variables\. Based on previous research\[[19](https://arxiv.org/html/2605.29021#bib.bib72)\], for a problem with this many candidate combinations, the use ofGNNimproves optimization efficiency by quickly identifying the best\-suited combinatorial variables\. The optimization is then defined as the following:

min𝐙comb,𝐗cont\\displaystyle\\min\_\{\\mathbf\{Z\}\_\{\\text\{comb\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\}f​\(𝐙comb,𝐗cont\)=fMU\\displaystyle f\(\\mathbf\{Z\}\_\{\\text\{comb\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)=f\_\{\\text\{MU\}\}\(6\)s\.t\.g​\(𝐙comb,𝐗cont\)=ICQI∗≤2\.5\\displaystyle g\(\\mathbf\{Z\}\_\{\\text\{comb\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)=I^\{\*\}\_\{\\text\{CQI\}\}\\leq 5h​\(𝐙comb,𝐗cont\)=NL=12\\displaystyle h\(\\mathbf\{Z\}\_\{\\text\{comb\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)=N\_\{L\}=2wherefMU=mprop\\displaystyle f\_\{\\text\{MU\}\}=m\_\{\\text\{prop\}\}𝐙comb=\[FT,max,Isp,mT,En,ρn,Nk,Kcls\]\\displaystyle\\mathbf\{Z\}\_\{\\text\{comb\}\}=\\begin\{bmatrix\}F\_\{T,\\text\{max\}\},I\_\{\\text\{sp\}\},m\_\{T\},E\_\{\\text\{n\}\},\\rho\_\{\\text\{n\}\},N\_\{k\},K\_\{\\text\{cls\}\}\\end\{bmatrix\}𝐗cont=\[Δ​x1​…​Δ​x4,Δ​y1​…​Δ​y4,v1​…​v4,mMU,rthread,rcorner,Lnet,Lct\]\\displaystyle\\mathbf\{X\}\_\{\\text\{cont\}\}=\\begin\{bmatrix\}\\begin\{aligned\} &\\Delta x\_\{1\}\\dots\\Delta x\_\{4\},\\ \\Delta y\_\{1\}\\dots\\Delta y\_\{4\},\\\\ &v\_\{1\}\\dots v\_\{4\},\\ m\_\{\\text\{MU\}\},\\ r\_\{\\text\{thread\}\},\\ r\_\{\\text\{corner\}\},\\ L\_\{\\text\{net\}\},\\ L\_\{\\text\{ct\}\}\\end\{aligned\}\\end\{bmatrix\}𝐙comb∈ℤ\\displaystyle\\mathbf\{Z\}\_\{\\text\{comb\}\}\\in\\mathbb\{Z\}ℤ=\[𝐙\(1\),𝐙\(2\),…,𝐙\(180\)\]\\mathbb\{Z\}=\[\\mathbf\{Z\}^\{\(1\)\},\\mathbf\{Z\}^\{\(2\)\},\\dots,\\mathbf\{Z\}^\{\(180\)\}\]is the set of all the 180 combinations of the combinatorial and integer variables in Table[II](https://arxiv.org/html/2605.29021#S2.T2)\. However, when dealing with a mixed\-combinatorial problem with this many combinations, it would be too time\-consuming to run the function evaluations on each combination, especially when the objective and constraints can only be calculated after completion of the time\-expensive, relatively high\-fidelity simulation of the system dynamics\. Therefore, it is beneficial to have an approximator that can quickly approximate which combination is the best\-suited one for the given vector of continuous variables, so that the combinatorial and integer variables no longer need to be co\-optimized with the continuous variables\. Assuming the approximator is represented byΘ\\Theta, the objective function can be transformed into the following:

min𝐗cont\\displaystyle\\min\_\{\\mathbf\{X\}\_\{\\text\{cont\}\}\}f​\(𝐙i∗,𝐗cont\)=fMU\\displaystyle f\(\\mathbf\{Z\}^\{\*\}\_\{\\text\{i\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)=f\_\{\\text\{MU\}\}\(7\)s\.t\.g​\(𝐙,𝐗cont\)=ICQI∗≤2\.5\\displaystyle g\(\\mathbf\{Z\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)=I^\{\*\}\_\{\\text\{CQI\}\}\\leq 5h​\(𝐙,𝐗cont\)=NL=12\\displaystyle h\(\\mathbf\{Z\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)=N\_\{L\}=2wherefMU​\(𝐙comb,𝐗cont\)=mprop\\displaystyle f\_\{\\text\{MU\}\}\\left\(\\mathbf\{Z\}\_\{\\text\{comb\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\\right\)=m\_\{\\text\{prop\}\}𝐲=Θ​\(𝒢,𝐗cont\)\\displaystyle\\mathbf\{y\}=\{\\Theta\}\(\\mathcal\{G\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)i∗=argmaxi=1,2,…,ntot​yi,yi∈𝐲\\displaystyle i^\{\*\}=\\underset\{i=1,2,\\dots,n\_\{\\text\{tot\}\}\}\{\\text\{argmax\}\}y\_\{i\},\\qquad y\_\{i\}\\in\\mathbf\{y\}𝐗cont=\[Δ​x1​…​Δ​x4,Δ​y1​…​Δ​y4,v1​…​v4,mMU,rthread,rcorner,Lnet,Lct\]\\displaystyle\\mathbf\{X\}\_\{\\text\{cont\}\}=\\begin\{bmatrix\}\\begin\{aligned\} &\\Delta x\_\{1\}\\dots\\Delta x\_\{4\},\\ \\Delta y\_\{1\}\\dots\\Delta y\_\{4\},\\\\ &v\_\{1\}\\dots v\_\{4\},\\ m\_\{\\text\{MU\}\},\\ r\_\{\\text\{thread\}\},\\ r\_\{\\text\{corner\}\},\\ L\_\{\\text\{net\}\},\\ L\_\{\\text\{ct\}\}\\end\{aligned\}\\end\{bmatrix\}𝐙i∗∈𝒢,and​𝒢∈ℤ\\displaystyle\\mathbf\{Z\}^\{\*\}\_\{i\}\\in\\mathcal\{G\},~\{\\text\{and\}\}~\\mathcal\{G\}\\in\\mathbb\{Z\}Here,𝐲\\mathbf\{y\}is a vector of scores of each node in the sampled subgraph, and𝒢\\mathcal\{G\}is agraphconstructed by all the combinations, which will be discussed more in the following section\. To handle the constraints, we propose using a penalty function to simplify the objective function, which can then be further transformed into the following:

min𝐗cont\\displaystyle\\min\_\{\\mathbf\{X\}\_\{\\text\{cont\}\}\}f​\(𝐙comb∗,𝐗cont\)=fMU\\displaystyle f\(\\mathbf\{Z\}^\{\*\}\_\{\\text\{comb\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)=f\_\{\\text\{MU\}\}\(8\)where𝐙comb∗=𝒵argmax​\(𝐲\)\\displaystyle\\mathbf\{Z\}^\{\*\}\_\{\\text\{comb\}\}=\\mathcal\{Z\}\_\{\\text\{argmax\}\(\\mathbf\{y\}\)\}𝐲=GNN\-NavCo​\(𝒢,𝐗cont\)\\displaystyle\\mathbf\{y\}=\\text\{GNN\-NavCo\}\(\\mathcal\{G\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)fMU​\(𝐙comb,𝐗cont\)=\{mprop,if​ICQI∗≤2\.5​and​NL=12ln⁡\(\(ICQI∗−2\.5\)2\+1\)\+ln⁡\(\(NL−12\)2\+1\)\+β,otherwise\\displaystyle f\_\{\\text\{MU\}\}\\left\(\\mathbf\{Z\}\_\{\\text\{comb\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\\right\)=𝐗cont=\[Δ​x1​…​Δ​x4,Δ​y1​…​Δ​y4,v1​…​v4,mMU,rthread,rcorner,Lnet,Lct\]\\displaystyle\\mathbf\{X\}\_\{\\text\{cont\}\}=\\begin\{bmatrix\}\\begin\{aligned\} &\\Delta x\_\{1\}\\dots\\Delta x\_\{4\},\\ \\Delta y\_\{1\}\\dots\\Delta y\_\{4\},\\\\ &v\_\{1\}\\dots v\_\{4\},\\ m\_\{\\text\{MU\}\},\\ r\_\{\\text\{thread\}\},\\ r\_\{\\text\{corner\}\},\\ L\_\{\\text\{net\}\},\\ L\_\{\\text\{ct\}\}\\end\{aligned\}\\end\{bmatrix\}whereGNN\-NavCo is theGNNrecommender that will be introduced in the following section,𝐲\\mathbf\{y\}is a vector of scores of each node in the sampled subgraph, andβ\\betais a penalty term, equal to 1\.5 times the maximum fuel cost observed in a successful capture333A failed capture was found to generally have a worse value offMUf\_\{\\text\{MU\}\}than a successful capture\.\. The simulator includes a maximum tension tolerance that checks whether the tension in the net thread exceeds the upper limit; readers may refer to\[[2](https://arxiv.org/html/2605.29021#bib.bib149)\]for details on how the tension limit is defined\. If tension in any thread exceeds the maximum, then the simulation is terminated and tagged as unsuccessful capture, and the values ofICQI∗I\_\{\\text\{CQI\}\}^\{\*\}andNLN\_\{L\}are set to be 50\.0 and 0, respectively, which are the worst observed values during preliminary simulations\.

## IIIGraph\-Learning\-Based Recommender for the Optimizer

### III\-APreprocessing and Creation of the Graph Dataset

The goal of theGNNis to act as a recommender, which, given the input of a graph of nodes and a continuous context, can rank them and give the best combinatorial node for that context\. In this research, theGNNtakes an undirected graph and a vector of continuous variables as input and outputs a directed graph whose edges encode the pairwise differences between nodes with respect to the objective function used to train theGNN\. In this concept, graph𝒢\\mathcal\{G\}is defined to be𝒢=\{𝒵,ℰ\}\\mathcal\{G\}=\\\{\\mathcal\{Z\},\\mathcal\{E\}\\\}, where𝒵\\mathcal\{Z\}is the node feature andℰ\\mathcal\{E\}is the edge feature\. In this work, assuming there arekcombk\_\{\\text\{comb\}\}combinatorial variables, then a unique node is defined to be𝐙=\[z1,z2,…,zkcomb∣zj∈ℝ,j=1,2,…,kcomb\]\\mathbf\{Z\}=\[z\_\{1\},z\_\{2\},\\dots,z\_\{k\_\{\\text\{comb\}\}\}\\mid z\_\{j\}\\in\\mathbb\{R\},\\,\\,j=1,2,\\dots,k\_\{\\text\{comb\}\}\], which represents the direct features of a combination in thekcombk\_\{\\text\{comb\}\}dimensional space, andzjz\_\{j\}is thejj\-th combinatorial variable\. Thus, the node feature𝒵∈ℝntot×kcomb\\mathcal\{Z\}\\in\\mathbb\{R\}^\{n\_\{\\text\{tot\}\}\\times k\_\{\\text\{comb\}\}\}is a vector containing all predefined valid combinations of the combinatorial variables\. All the nodes are defined to be connected with each other in this graph to form a fully connected graph\[[24](https://arxiv.org/html/2605.29021#bib.bib71)\], and the edge featureℰ∈ℝntot×ntot\\mathcal\{E\}\\in\\mathbb\{R\}^\{n\_\{\\text\{tot\}\}\\times n\_\{\\text\{tot\}\}\}is then defined as an adjacency matrix between each node encoding the pairwise differences of the objective function between them, parametrized by the continuous context\.

For a complex problem, the combinatorial variable numberkcombk\_\{\\text\{comb\}\}and the number of valid combinationsmmcan be large\. In that case, the dimensions of the node feature𝒵\\mathcal\{Z\}and edge featureℰ\\mathcal\{E\}are large for training theGNNand for further analysis\. Therefore, instead of directly working on the graph𝒢\\mathcal\{G\}with all the nodes and edges, we propose to randomly picknsnn\_\{\\text\{sn\}\}nodes from𝒵\\mathcal\{Z\}and their corresponding edges fromℰ\\mathcal\{E\}to create a set ofpsgp\_\{\\text\{sg\}\}subgraphs\. The set of subgraphsSSis defined as:

S=\{𝒢1,𝒢2,…,𝒢psg∣𝒢q⊂𝒢,∀q=1,2,…,psg\}S=\\\{\\mathcal\{G\}\_\{1\},\\mathcal\{G\}\_\{2\},\\dots,\\mathcal\{G\}\_\{p\_\{\\text\{sg\}\}\}\\mid\\mathcal\{G\}\_\{q\}\\subset\\mathcal\{G\},\\,\\forall~q=1,2,\\dots,p\_\{\\text\{sg\}\}\\\}\(9\)wherepsgp\_\{\\text\{sg\}\}is the number of subgraphs, and each subgraph𝒢q\\mathcal\{G\}\_\{q\}in the set is defined to be:

𝒢q\\displaystyle\\mathcal\{G\}\_\{q\}=\{𝒵q,ℰq∣𝒵q∈𝒵,ℰq∈ℰ,∀q=1,2,…,psg\}\\displaystyle=\\\{\\mathcal\{Z\}\_\{q\},\\mathcal\{E\}\_\{q\}\\mid\\mathcal\{Z\}\_\{q\}\\in\\mathcal\{Z\},\\mathcal\{E\}\_\{q\}\\in\\mathcal\{E\},\\,\\forall~q=1,2,\\dots,p\_\{\\text\{sg\}\}\\\}𝒵q=\(𝐙q,1𝐙q,2⋮𝐙q,nsn\),ℰq=\(𝐝q,1𝐝q,2⋮𝐝q,nsn\),\\displaystyle\\mathcal\{Z\}\_\{q\}=\\begin\{pmatrix\}\\mathbf\{Z\}\_\{q,1\}\\\\ \\mathbf\{Z\}\_\{q,2\}\\\\ \\vdots\\\\ \\mathbf\{Z\}\_\{q,n\_\{\\text\{sn\}\}\}\\end\{pmatrix\},\\,\\,\\mathcal\{E\}\_\{q\}=\\begin\{pmatrix\}\\mathbf\{d\}\_\{q,1\}\\\\ \\mathbf\{d\}\_\{q,2\}\\\\ \\vdots\\\\ \\mathbf\{d\}\_\{q,n\_\{\\text\{sn\}\}\}\\end\{pmatrix\},\(10\)𝐙q,i∈𝒵,𝐝q,i∈ℰ,∀i=1,2,…,nsn\\displaystyle\\quad\\mathbf\{Z\}\_\{q,i\}\\in\\mathcal\{Z\},\\mathbf\{d\}\_\{q,i\}\\in\\mathcal\{E\},\\,\\forall~i=1,2,\\dots,n\_\{\\text\{sn\}\}where the notation\(⋅\)q,i\(\\cdot\)\_\{q,i\}denotes the operation of randomly picking an unrepeated value between\[1,ntot\]\[1,n\_\{\\text\{tot\}\}\]\. The variable𝐝i\\mathbf\{d\}\_\{i\}is the vector of pairwise differences in objective value parametrized by the continuous context of theii\-th node to every node in the graph\.

Then, within a subgraph𝒢q\\mathcal\{G\}\_\{q\}, the objective function will be evaluatednsnn\_\{\\text\{sn\}\}times, with each node providing the combinatorial variables, and𝐗q\\mathbf\{X\}\_\{q\}providing the continuous variables; thus, all nodes use the same vector of continuous values for function evaluations\. The objective values for all the nodes in𝒢q\\mathcal\{G\}\_\{q\}and the continuous variables in𝐗q\\mathbf\{X\}\_\{q\}is then defined as:

𝐅q\\displaystyle\\mathbf\{F\}\_\{q\}=\(fq,1fq,2⋮fq,nsn\),fq,i=f​\(𝐙q,i,𝐗q\)\\displaystyle=\\begin\{pmatrix\}f\_\{q,1\}\\\\ f\_\{q,2\}\\\\ \\vdots\\\\ f\_\{q,n\_\{\\text\{sn\}\}\}\\end\{pmatrix\},\\,f\_\{q,i\}=f\(\\mathbf\{Z\}\_\{q,i\},\\mathbf\{X\}\_\{q\}\)\(11\)∀i=1,2,…,nsn\\displaystyle\\quad\\forall~i=1,2,\\dots,n\_\{\\text\{sn\}\}
To train theGNNas a recommender \(referred to asGNN\-NavCo in the rest of the paper\) to find the best\-suited combination, a pairwise regression loss function such asMean Squared Error\(MSE\) or Huber loss is used during training\. Huber loss is a special loss function that acts asMSEloss function below a threshold andMean Absolute Error\(MAE\) loss above it\. Here, the special cyclical consistency regularizer is also used, which penalizes nonzero circulation around graph cycles, encouraging globally consistent predictions\. The variableλcycle\\lambda\_\{\\text\{cycle\}\}is the regularization coefficient that weights the cyclical consistency term relative to the primary pairwise regression objective\. Letf^​\(i,j;𝐗\)\\hat\{f\}\(i,j;\{\\mathbf\{X\}\}\)denote the predicted edge difference\. Learning proceeds by minimizing the discrepancy between predicted and observed differences over the edges of the graph, whereℓ​\(⋅\)\\ell\(\\cdot\)denotes a regression loss where edge\-wise loss is defined in Eq\. \([12](https://arxiv.org/html/2605.29021#S3.E12)\)\. Mathematically, cycle consistency penalizes nonzero circulation over triangles, which are the smallest cycles in a fully connected graph\. For a triangle\(i,j,k\)\(i,j,k\), the circulation is defined in Eq\. \([13](https://arxiv.org/html/2605.29021#S3.E13)\), and ideallyc​\(i,j,k;𝐗\)=0c\(i,j,k;\{\\mathbf\{X\}\}\)=0\. If𝒯\\mathcal\{T\}denotes the set of all triangles in the graph, then the cycle consistency loss is defined in Eq\. \([14](https://arxiv.org/html/2605.29021#S3.E14)\)\. Thus, the overall objective isL=Ledge\+λcycle​LcycleL=L\_\{\\text\{edge\}\}\+\\lambda\_\{\\text\{cycle\}\}L\_\{\\text\{cycle\}\}\.

Ledge=𝔼𝐗​\[1\|ℰq​\(𝒢q\)\|​∑\(i,j\)∈ℰq​\(𝒢q\)ℓ​\(f^q​\(i,j;𝐗\)−fq​\(i,j\)\)\]\\displaystyle L\_\{\\text\{edge\}\}=\\mathbb\{E\}\_\{\\mathbf\{X\}\}\\left\[\\frac\{1\}\{\|\\mathcal\{E\}\_\{q\}\(\\mathcal\{G\}\_\{q\}\)\|\}\\sum\_\{\(i,j\)\\in\\mathcal\{E\}\_\{q\}\(\\mathcal\{G\}\_\{q\}\)\}\\ell\\big\(\\hat\{f\}\_\{q\}\(i,j;\{\\mathbf\{X\}\}\)\-f\_\{q\}\(i,j\)\\big\)\\right\]\(12\)
c​\(i,j,k;𝐗\)=f^​\(i,j;𝐗\)\+f^​\(j,k;𝐗\)\+f^​\(k,i;𝐗\)\\displaystyle c\(i,j,k;\{\\mathbf\{X\}\}\)=\\hat\{f\}\(i,j;\{\\mathbf\{X\}\}\)\+\\hat\{f\}\(j,k;\{\\mathbf\{X\}\}\)\+\\hat\{f\}\(k,i;\{\\mathbf\{X\}\}\)\(13\)
Lcycle=𝔼𝐗​\[1\|𝒯\|​∑\(i,j,k\)∈𝒯c​\(i,j,k;𝐗\)2\]\\displaystyle L\_\{\\text\{cycle\}\}=\\mathbb\{E\}\_\{\\mathbf\{X\}\}\\left\[\\frac\{1\}\{\|\\mathcal\{T\}\|\}\\sum\_\{\(i,j,k\)\\in\\mathcal\{T\}\}c\(i,j,k;\{\\mathbf\{X\}\}\)^\{2\}\\right\]\(14\)
In terms of the tether\-net optimization problem, the continuous variables in the top part of Table[II](https://arxiv.org/html/2605.29021#S2.T2)are used to define𝐗q\\mathbf\{X\}\_\{q\}; the combinatorial and integer variables in the bottom part of Table[II](https://arxiv.org/html/2605.29021#S2.T2)are used to create the combinations, which then formulated into𝐙q\\mathbf\{Z\}\_\{q\}\.

![Refer to caption](https://arxiv.org/html/2605.29021v1/x3.png)Figure 2:GNN\-aided Optimization Framework![Refer to caption](https://arxiv.org/html/2605.29021v1/x4.png)Figure 3:Sketch of theGNN\-NavCo\-aided optimization using population\-based optimizer as an example\. The graphs are simplified for illustrative purposes\. At the beginning of each iteration, a subgraph is sampled with the optimal node from the previous iteration \(red node\) included\. During the Least Square Fitness Recovery process, the fitness score of each node \(red and orange nodes\) in the subgraph is evaluated based on the edge values\. A greedy filter is applied to find the node with the optimal score \(blue node\)\. In the top figure, between two nodes in the subgraph, the blue arrow indicates the direction that leads to a node with a better score, and the purple arrow indicates the direction that leads to a node with a worse score\. The optimal node, along with the optimized continuous vectorX, will be used in the function evaluation in the optimization\.
### III\-BFramework of GNN\-aided Optimization Algorithm

Figure[2](https://arxiv.org/html/2605.29021#S3.F2)shows the overall structure of the proposed framework\. The framework consists of two major components: theGNN\-NavCo training and the optimization loop\. The framework starts with creating the dataset by randomly sampling valid combinations to create a pool of subgraphs for training theGNN\-NavCo\. To trainGNN\-NavCo, the input consists of the sampled subgraph providingnsnn\_\{\\text\{sn\}\}combinatorial nodes and a common continuous vector\. After this, the fitness scores of the combination are computed and stored on the nodes\. The subgraph is a complete graph, meaning all nodes are connected to each other\. Here, we capture the pairwise difference in fitness scores at nodes and move it into the edges, forming the ground truth of the abstraction\. Then the fitness scores, represented by the objective valuefMUf\_\{\\text\{MU\}\}here, are sent to the loss function to calculate the gradient and updateGNN\-NavCo\. The main objective loss is a pairwise regression loss\. Therefore, theGNN\-NavCo is trained to predict the edge\-wise pairwise difference in fitness scores between nodes in a subgraph given a continuous vector; in other words, theGNN\-NavCo learns a context\-conditioned parametric edge function over the subgraph that models pairwise objective differences between combinatorial choices\. This gives rise to an inference time phenomenon, where we input the set of undirected subgraphs, with nodes containing the combinatorial vector and a continuous context, andGNN\-NavCo outputs a directed graph, with pairwise differences in fitness scores encoded on the edges\. Using this, ideally, the optimizer can make an informed decision to move to the best combinatorial node within the subgraph for the fixed best continuous vector provided by an outer optimizer, such as PSO\. Thus, it achieves a notion of explicit gradient awareness over combinatorial spaces\.

AfterGNN\-NavCo is trained, it helps the optimizer select the candidate combination from a subgraph generated from sampled combinations and the candidate combination from the previous iteration\. The sampled combinations are randomly sampled from the graph with respect to the candidate combination\. Once the directed graph is encoded with pairwise fitness score differences, the node\-level fitness scores are recovered by solving a least\-squares problem that minimizes the discrepancy between the computed edge differences and the induced node fitness scores\. This yields an approximate global ranking over nodes\[[11](https://arxiv.org/html/2605.29021#bib.bib362)\]\. After this, a greedy filter is applied to select the node with the lowest score\. The selected candidate combination,Zq∗\\textbf\{Z\}^\{\*\}\_\{q\}, is used in the function evaluation with the continuous vector provided by the optimizer, and it is also stored in the cache, so it becomes the new initial candidate combination in the next iteration\. The continuous vector will be updated by the optimizer\. The optimization iterates until the maximum number of iterations is reached or the result converges\. Here, theMixed\-Discrete Particle Swarm Optimization\(MDPSO\) algorithm\[[8](https://arxiv.org/html/2605.29021#bib.bib10),[27](https://arxiv.org/html/2605.29021#bib.bib9)\]is used as the optimizer, which is a popular approach for solving complex nonlinear optimization problems\. However, it should be noted that the optimizer can be any gradient\- or sampling\-based algorithm, such as a Genetic Algorithm or sequential quadratic programming\. To validate that the converged designs are local optima, after the optimization converges, the combinatorial variables are frozen, and a sequential quadratic programming solver is applied to the continuous variables for post\-optimization\. The trainedGNNcan efficiently pick the node to run the function evaluations without iterating over all the candidates\. ThisGNN\-NavCo\-aided MDPSO is defined as in Algorithm[1](https://arxiv.org/html/2605.29021#alg1), and the illustrative flowchart is provided in Figure[3](https://arxiv.org/html/2605.29021#S3.F3)\.

In this research, the population of theMDPSOis set to 100, and the max iteration limit is set to 200\. The optimization is terminated early if either of the following conditions is met: i\) the population remains entirely infeasible and the total constraint violation of the global best does not decrease \(within a tolerance of10−510^\{\-5\}\) over 15 consecutive iterations; or ii\) the global best is feasible but its objective value does not improve \(within a tolerance of10−510^\{\-5\}\) over 15 consecutive iterations\. All the other hyperparameters of theMDPSOare set to their default values\.

Algorithm 1GNN\-aided MDPSO1:Initialize the particle population

𝐗1\\mathbf\{X\}\_\{1\}to

𝐗Npop\\mathbf\{X\}\_\{N\_\{\\text\{pop\}\}\}
2:Initialize the optimal node

𝐙1∗\\mathbf\{Z\}^\{\*\}\_\{1\}to

𝐙Npop∗\\mathbf\{Z\}^\{\*\}\_\{N\_\{\\text\{pop\}\}\}
3:Function evaluation

f​\(𝐙i∗,𝐗i\)f\(\\mathbf\{Z\}^\{\*\}\_\{i\},\\mathbf\{X\}\_\{i\}\),

i=1i=1to

NpopN\_\{\\text\{pop\}\}
4:for

t=1t=1to

m​a​x​\_​i​t​e​r​a​t​i​o​nmax\\\_iterationdo

5:Update particles

𝐗1\\mathbf\{X\}\_\{1\}to

𝐗Npop\\mathbf\{X\}\_\{N\_\{\\text\{pop\}\}\}withMDPSO

6:for

q=1q=1to

NpopN\_\{\\text\{pop\}\}do

7:Select

nsn−1n\_\{\\text\{sn\}\}\-1nodes randomly around each

𝐙q∗\\mathbf\{Z\}^\{\*\}\_\{q\}from

𝒢\\mathcal\{G\}
8:Create

𝒢q\\mathcal\{G\}\_\{q\}
9:run

𝒢q∗=\\mathbf\{\\mathcal\{G\}\}^\{\*\}\_\{q\}=GNN\-NavCo\(

𝒢q,𝐗q\\mathcal\{G\}\_\{q\},\\mathbf\{X\}\_\{q\}\)

10:\#recover node score from adjacency matrix

11:

Φi=1nsn​∑j=1nsnℰi​j∗\\Phi\_\{i\}=\\frac\{1\}\{n\_\{\\text\{sn\}\}\}\\sum\_\{j=1\}^\{n\_\{\\text\{sn\}\}\}\\mathcal\{E\}^\{\*\}\_\{ij\}
12:Set

Zq∗=arg⁡mini⁡Φi\\textbf\{Z\}^\{\*\}\_\{q\}=\\arg\\min\_\{i\}\\Phi\_\{i\}
13:endfor

14:Function evaluation

f​\(𝐙i∗,𝐗i\)f\(\\mathbf\{Z\}^\{\*\}\_\{i\},\\mathbf\{X\}\_\{i\}\),

i=1i=1to

NpopN\_\{\\text\{pop\}\}
15:ifConvergethen

16:Return the result

17:endif

18:endfor

19:Return the result

### III\-CMDPSO without the aid of the GNN\-NavCo

To investigate the contribution ofGNN\-NavCo to the optimization, another round ofMDPSOis performed withoutGNN\-NavCo and uses an integer variable \(xint∈\[1,180\]x\_\{\\text\{int\}\}\\in\[1,180\]\) to represent the choice of the combination among the feasible set\. All the settings of thisMDPSOare exactly the same as the GNN\-NavCo\-aidedMDPSO\. The initial populations, including both the continuous and projected combinatorial variables, are also kept the same between the two methods\. The objective function for this round of optimization is defined in Eq\. \([15](https://arxiv.org/html/2605.29021#S3.E15)\)\. The convergence history and optimized results are provided in the following Results section\.

minxint,𝐗cont\\displaystyle\\min\_\{x\_\{\\text\{int\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\}f​\(xint,𝐗cont\)=fMU\\displaystyle f\(x\_\{\\text\{int\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\)=f\_\{\\text\{MU\}\}\(15\)wherefMU​\(𝐙comb,𝐗cont\)=\{mprop,if​ICQI∗≤2\.5​and​NL=12ln⁡\(\(ICQI∗−2\.5\)2\+1\)\+ln⁡\(\(NL−12\)2\+1\)\+β,otherwise\\displaystyle f\_\{\\text\{MU\}\}\\left\(\\mathbf\{Z\}\_\{\\text\{comb\}\},\\mathbf\{X\}\_\{\\text\{cont\}\}\\right\)=𝐙comb=𝐙\(xint\)\\displaystyle\\mathbf\{Z\}\_\{\\text\{comb\}\}=\\mathbf\{Z\}^\{\(x\_\{\\text\{int\}\}\)\}𝐗cont=\[Δ​x1​…​Δ​x4,Δ​y1​…​Δ​y4,v1​…​v4,mMU,rthread,rcorner,Lnet,Lct\]\\displaystyle\\mathbf\{X\}\_\{\\text\{cont\}\}=\\begin\{bmatrix\}\\begin\{aligned\} &\\Delta x\_\{1\}\\dots\\Delta x\_\{4\},\\ \\Delta y\_\{1\}\\dots\\Delta y\_\{4\},\\\\ &v\_\{1\}\\dots v\_\{4\},\\ m\_\{\\text\{MU\}\},\\ r\_\{\\text\{thread\}\},\\ r\_\{\\text\{corner\}\},\\ L\_\{\\text\{net\}\},\\ L\_\{\\text\{ct\}\}\\end\{aligned\}\\end\{bmatrix\}

## IVResults and Discussion

The training and testing of theGNN\-ReCo, and the optimization are performed on a dual\-socket Intel Xeon Gold 6148 system \(40 physical cores, 80 threads\) with 187 GB RAM and a 16 GB VRAM NVIDIA GeForce RTX 4060 Ti GPU\. In this section, we present an analysis of the training performance of GNN recommender, followed by discussion of the optimization history and results, and its comparison to direct optimization applied to the MCNLP problem as is \(without a GNN\)\. The section ends with physical insights into the optimized designs, and their simulated capture performance\.

### IV\-ATraining of theGNNRecommender

The currentGNNarchitecture utilizes a 3\-layerGNNbackbone with a hidden dimension of 128\. The embeddings are scaled to 64 dimensions to provide richer feature representations\. To map graph features to cost predictions, the model employs a deep, 5\-layer residual Multi\-Layer Perceptron \(MLP\)\. Finally, to prevent overfitting given the increased capacity, a moderate dropout rate of 30% is applied throughout the network\. For data collection, 3000 function evaluations are sampled from the high\-fidelity simulator, and 100 subgraphs are constructed, each with 30 nodes representing combinatorial features and one common continuous variable\. For all simulations, the target is initialized at a position ofrt=9​𝐢^\+9​𝐣^−60​𝐤^\\textbf\{r\}\_\{t\}=9\\mathbf\{\\hat\{i\}\}\+9\\mathbf\{\\hat\{j\}\}\-60\\mathbf\{\\hat\{k\}\}m\. Its initial orientation is specified by a 3\-2\-1 Euler angles set \(i\.e\., the set of rotations describing the orientation of the𝒟\\mathcal\{D\}frame with respect to the𝒪\\mathcal\{O\}frame\) of\{60\.0∘,40\.0∘,0\.0∘\}\\\{60\.0^\{\\circ\},40\.0^\{\\circ\},0\.0^\{\\circ\}\\\}, while its initial angular velocity \(expressed in the𝒟\\mathcal\{D\}frame\) is defined as𝝎t=10​𝐝^x\+30​𝐝^y\+10​𝐝^z\\boldsymbol\{\\omega\}\_\{t\}=10\\mathbf\{\\hat\{d\}\}\_\{x\}\+30\\mathbf\{\\hat\{d\}\}\_\{y\}\+10\\mathbf\{\\hat\{d\}\}\_\{z\}deg/s\. This described target condition is considered one of the most challenging target\-state scenarios in our previous research\[[3](https://arxiv.org/html/2605.29021#bib.bib99),[16](https://arxiv.org/html/2605.29021#bib.bib69)\], as it has both the greatest initial distance from the chaser and the highest angular velocity\.

![Refer to caption](https://arxiv.org/html/2605.29021v1/x5.png)Figure 4:Shows the performance on training data used to train the GNN\-NavCo model\. The blue line denotes the min\-max range of fitness values each node had and the red dot shows the median of fitness valueIn Fig\.[4](https://arxiv.org/html/2605.29021#S4.F4), it can be observed that the fitness scorefMUf\_\{\\text\{MU\}\}of the collected data varies from near 0 to below 14, which indicates the data covers various successful and failed capture scenarios\. By definition offMUf\_\{\\text\{MU\}\}in Eq\.[8](https://arxiv.org/html/2605.29021#S2.E8), the maximum value for a successful capture will be 0\.96\. It can be observed that the median fitness values of most of the subgraphs did not exceed the threshold, indicating the objective values were strongly right\-skewed\. But individual subgraphs still exhibit substantial within\-graph variation, as shown by the min\-max ranges and mean values in the figure\. This variation is useful for training because it exposes the network to both relatively low and high fitness value regions within the same graph\. The successful capture rate in the 3,000 samples is 74\.5%\. Therefore, it is still necessary to apply the optimization to ensure the capture is successful\.

With the collected data, theGNN\-NavCo model is trained while tracking edge and cycle\-consistency regularization losses\. Model selection is performed using the minimum validation edge loss, and the corresponding checkpoint is saved at 581 epochs\. In Fig\.[5](https://arxiv.org/html/2605.29021#S4.F5), it is observed that the edge loss remains stagnant early on, which can be attributed to the model’s early\-stage learning of the hidden representation and to the stabilizing effect of the cycle consistency regularization\. When edge loss and cycle loss in Fig\.[5](https://arxiv.org/html/2605.29021#S4.F5)are compared, it is observed that up until 50 epochs, the cycle loss reduces rapidly\. After this point, the edge loss drops sharply, whereas the cycle loss largely plateaus with occasional spikes, although the overall spike is negligible\. A similar edge\-loss pattern is also observed across experiments with different values of the cycle\-loss weightλc​y​c​l​e\\lambda\_\{cycle\}\. The abrupt improvement in edge loss is not solely due to the regularizer, but rather to the model learning a representation for edge prediction\. This can be observed in Fig\.[5](https://arxiv.org/html/2605.29021#S4.F5), where it is seen that before roughly epochs 50\-60, the model has fluctuating sign accuracy, and after the stabilization mark \(green line\), loss drops and the sign accuracy stabilizes at 75\-76%\. Sign accuracy measures whether the model correctly predicts the direction of an edge relation\. It can be inferred that early on, the model prioritizes both cyclical consistency and accurate sign prediction, and once it stabilizes, it shifts its focus to reducing edge loss\. A hyperparameter study, as shown in Table[III](https://arxiv.org/html/2605.29021#S4.T3)overλc​y​c​l​e\\lambda\_\{cycle\}, shows that sign accuracy was broadly similar across the tested values, withλc​y​c​l​e=0\.003\\lambda\_\{cycle\}=0\.003giving a slightly better overall result; therefore,λc​y​c​l​e=0\.003\\lambda\_\{cycle\}=0\.003was selected for subsequent experiments\. Our parallel research suggests that cycle\-consistency regularization can improve generalization, especially in challenging settings such as distribution shift or skewed targets\. In the present study, however, the effect ofλc​y​c​l​e\\lambda\_\{cycle\}was modest, with similar sign accuracy across the tested values\. In summary, the trainedGNNcan serve as aninstantaneousrecommender because of the negligible execution time of the neural network\. Given a continuous vector and a starting combination, theGNNhas a general 75% accuracy in recommending a more optimal combination for the continuous vector\. The accuracy may seem low, but during the optimization framework, the optimal design is cached thus the actual achieved performance of theGNN\-aided optimization is not limited much by this accuracy\. The accuracy can also be further improved by future research on theGNNstructure, loss function, or subgraph sampling method, but for the current paper, this trainedGNNproves to work as Section[IV\-B](https://arxiv.org/html/2605.29021#S4.SS2)shows\.

![Refer to caption](https://arxiv.org/html/2605.29021v1/x6.png)Figure 5:The top figure shows the training and validation loss of the edge loss function, followed by the cycle loss, and the third one shows the sign accuracy curve of the trained model\. The green line shows the stabilization point\.HyperparameterSign AccuracyEdge LossCycle Lossλc​y​c​l​e=1\\lambda\_\{cycle\}=10\.763330\.043211\.16×10−121\.16\\times 10^\{\-12\}λc​y​c​l​e=0\.003\\lambda\_\{cycle\}=0\.0030\.763560\.043232\.32×10−122\.32\\times 10^\{\-12\}λc​y​c​l​e=0\\lambda\_\{cycle\}=00\.76230\.0432241\.85×10−121\.85\\times 10^\{\-12\}TABLE III:Experiment on cyclical consistency regularization\. sign accuracy, edge loss and cycle loss denote the corresponding converged value at least validation loss
### IV\-BOptimization Performance and Comparisons between with/without GNN recommendation

The total computing time to finish 53 iterations ofGNN\-NavCo\-aided MDPSO is 26\.5 hours\. The convergence history of the objective is shown in Fig\.[6](https://arxiv.org/html/2605.29021#S4.F6)\. The sequential quadratic programming optimizer \(implemented by MATLABfmincon\(\)\) is applied on the continuous variables optimized byGNN\-NavCo\-aided MDPSO to further investigate the potential improvement, but failed to find any more optimal results, proving the current result is already a local optimal\. The original MDPSO converged at the 27th iteration, while theGNN\-NavCo\-aided converged at the 8th iteration\. Here, convergence is defined as the point at which the absolute difference between the current objective value falls into the lower 0\.5% of the total range of objective values observed over the entire optimization history\. Though the MDPSO without GNN eventually reached to an objective slightly better than the GNN\-aided one, it requires more iterations to reach the optimized objective of the GNN\-aided MDPSO, which terminated earlier due to the stopping criteria mentioned in Section[III\-B](https://arxiv.org/html/2605.29021#S3.SS2)\. In summary, theGNN\-NavCo used 1,900 fewer function evaluations during the optimization to reach convergence\. It needs to be noted that thetotalfunction evaluations used by theGNN\-NavCo\-aidedMDPSOare 1,100 more \(considering the 3,000 function evaluations used for training the GNN\) than those of theMDPSOwithout the aid ofGNN\. However, theGNN\-NavCo\-aided methods remain promising\.

First, the additional function evaluations are primarily associated with the data\-generation and training phases of theGNNsurrogate, which can be performedoffline\. This process is highly parallelizable and can be efficiently executed on high\-performance computing \(HPC\) platforms or supercomputers, thereby significantly reducing the computing time despite the increased number of evaluations\. In contrast, the baselineMDPSOrelies on sequential updates between iterations as online evaluations during optimization, which are typically more time\-constrained\.

Second, from a scalability perspective, the advantage of theGNN\-NavCo framework is expected to become more pronounced as the size and complexity of the combinatorial space increase\. The ability to reduce data\-generation costs and leverage learned structural information makes the approach more suitable for large\-scale problems, where purely optimization\-based methods may struggle to efficiently explore the search space\. Our ongoing parallel research will further investigate this statement\.

![Refer to caption](https://arxiv.org/html/2605.29021v1/x7.png)Figure 6:Objective convergence history comparison between the optimization with and without the aid ofGNN\-NavCo
### IV\-CPhysical Insights on the Optimized and Baseline Designs and Their Capture Performance

The comparison of the design variables and optimized objectives is shown in Table[IV](https://arxiv.org/html/2605.29021#S4.T4)\. Both optimization methods converged on the same combination: a thruster with the highestIspI\_\{\\text\{sp\}\}and the 2nd\-highest thrust, and a net with the fewest meshes, which together significantly reduce the system’s fuel consumption\. It should be noted that although the combinatorial variables are identical between the two methods, several continuous variables—such as the aiming points and the net side length—differ significantly\. This observation suggests that the underlying optimization problem is highly nonlinear and admits multiple feasible solutions with comparable performance\. In other words, even under the same discrete configuration, the continuous design space likely contains multiple local optima or near\-optimal regions that can yield similar objective values\. The differences observed in the continuous variables indicate that the two methods have converged to distinct solutions within this landscape, rather than a single unique optimum\. This highlights the problem’s inherent multi\-modality, in which different values of the continuous variables can compensate for one another to achieve similar system performance\. As a result, the final solution is not necessarily unique, and the optimizer may converge to different, yet equally valid, designs even though they have exactly the same initial populations\. The rendered simulation screenshots of theGNN\-NavCo\-aided optimized design are shown in Fig\.[7](https://arxiv.org/html/2605.29021#S4.F7)\. The animated capture of the three methods are provided in a GitHub repository\[[18](https://arxiv.org/html/2605.29021#bib.bib222)\]\. Since all three methods achieved successful capture in this scenario, here thefMUf\_\{\\text\{MU\}\}equals the fuel cost\. Compared with the baseline design from our previous study\[[3](https://arxiv.org/html/2605.29021#bib.bib99)\], the optimized design consumes only 14\.35% of the baseline’s fuel, highlighting the critical role of optimization \(in conjunction with the MPC\) in remarkably improving system fuel efficiency\.

TABLE IV:Comparison between optimized variables and objective with and without the aid ofGNN\-NavCoContinuous VariableGNN\-NavCo\-aided MDPSOMDPSOBaselineΔ​x1,…​Δ​x4\\Delta x\_\{1\},\\dots\\Delta x\_\{4\}\[m\]\[2\.843, \-4\.980, 0\.041, \-4\.997\]\[\-1\.132, \-4\.060, 4\.286, \-4\.999\]\[0,0,0,0\]Δ​y1,…​Δ​y4\\Delta y\_\{1\},\\dots\\Delta y\_\{4\}\[m\]\[0\.569, 4\.103, \-1\.799, \-0\.880\]\[1\.794, 1\.963, \-4\.144, \-2\.404\]\[0,0,0,0\]v1,…,v4v\_\{1\},\\dots,v\_\{4\}\[m/s\]\[0\.527, 0\.785, 0\.837, 0\.841\]\[0\.654, 0\.733, 0\.846, 0\.767\]\[0,0,0,0\]mMUm\_\{\\text\{MU\}\}\[kg\]2\.18622\.75692\.0rthreadr\_\{\\text\{thread\}\}\[m\]0\.00110\.00110\.0010rcornerr\_\{\\text\{corner\}\}\[m\]0\.00120\.00090\.0010LnetL\_\{\\text\{net\}\}\[m\]23\.31519\.66722\.0LctL\_\{\\text\{ct\}\}\[m\]1\.5320\.9312\.414Combinatorial VariableFT,maxF\_\{T,\\text\{max\}\}\[N\]6\.16\.18\.9IspI\_\{\\text\{sp\}\}\[s\]277\.0277\.060\.0mTm\_\{T\}\[kg\]0\.60\.60\.37EnE\_\{\\text\{n\}\}\[Gpa\]70\.070\.070\.0ρn\\rho\_\{\\text\{n\}\}\[kg/m3\]1390\.01390\.01390\.0NkN\_\{k\}9911KclsK\_\{\\text\{cls\}\}001ObjectivefMUf\_\{\\text\{MU\}\}0\.02440\.02440\.170![Refer to caption](https://arxiv.org/html/2605.29021v1/x8.png)Figure 7:Rendered simulation screenshots of theGNN\-NavCo\-aided optimized design\. The dashed lines are the minimum\-energy reference trajectories for theMUs\. The following frames are shown in the above figure\.1\) t=0s:The beginning of the entire mission\.2\) t=10s and 18s:The MUs are maneuvered towards the optimized aiming points\.3\) t=25s:The thrusters are switched off, and the closing mechanism is switched on\.4\) t=34\.9s:The capture is stabilized\.

## VConclusion

In this paper, we developed aGNN\-aided nonlinear optimization framework to search the relatively high\-dimensional mixed\-combinatorial space presented by joint consideration of the physical design of the net, component choices and aiming point control decisions for a tether\-net system used in space debris capture\. The presented framework found a converged physical design and control that provided a 7×\\timesreduction in fuel costs of the tether\-net maneuver compared to a baseline design from previously reported work, while maintaining successful capture\. While both optimizations were allowed to run for a bit more than 50 iterations using the PSO algorithm \(post\-processed by SQP\), the new GNN\-aided PSO approach achieved the threshold convergence tolerance of within 0\.0001 kg of the optimum fuel cost in only 800 function evaluations, compared to the 2700 function evaluations taken by the standard PSO \(that solves the entire MCNLP without GNN\) to achieve the same\. In terms of compute time, this was a savings of 9\.5 hours, relative to overall run time in the order of about 25 hours\.

It is also important to note that the integration of theGNN\-NavCo module introduces a learning\-based guidance mechanism that captures structural relationships within the combinatorial design space, enabling more informed exploration compared to purely stochastic search\. Sample generation and training the GNN is therefore a one\-time overhead, which can then be used with multiple different optimization runs with different settings\. This computational overhead can also be largely amortized through offline processing and parallel execution on high\-performance computing platforms\.

In its current form, the design space does not include sensor choices or a variety of net shapes, addition of which would provide more meaningful exploration of the design space for autonomous tether\-net systems\. In addition, application of theGNN\-NavCo module with other optimization algorithms, especially gradient based methods and Bayesian optimization, with the potential to update the GNN insitu during optimization could provide further insights into the generalized applicability of this approach for automated design guidance of mixed\-combinatorial complex systems\.

## Acknowledgments

This work is supported under the CMMI Awards numbered 2128578 and 2048020 from the National Science Foundation \(NSF\)\. The authors’ opinions, findings, and conclusions or recommendations expressed in this material do not necessarily reflect the views of the National Science Foundation\. The University at Buffalo Presidential Fellowship is acknowledged by Achira Boonrath\.

## References

- \[1\]P\. Anz\-Meador and J\. Liou\(July 2010\)Analysis and Consequences of the Iridium 33\-Cosmos 2251 Collision\.In38th COSPAR Scientific Assembly,Bremen, Germany\.Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1)\.
- \[2\]A\. Boonrath and E\. M\. Botta\(2025\)Robustness and safety of net\-based debris capture under deployment and environmental uncertainties\.Journal of Spacecraft and Rockets,pp\. 1–17\.External Links:[Document](https://dx.doi.org/10.2514/1.A36217)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1),[§II\-C](https://arxiv.org/html/2605.29021#S2.SS3.p6.6)\.
- \[3\]A\. Boonrath, F\. Liu, E\. M\. Botta, and S\. Chowdhury\(2024\)Learning\-aided control of robotic tether\-net with maneuverable nodes to capture large space debris\.IEEE 2024 International Conference on Robotics and Automation\.External Links:[Document](https://dx.doi.org/10.1109/ICRA57147.2024.10610721)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1),[§II\-B](https://arxiv.org/html/2605.29021#S2.SS2.p1.1),[§II\-B](https://arxiv.org/html/2605.29021#S2.SS2.p3.4),[§IV\-A](https://arxiv.org/html/2605.29021#S4.SS1.p1.6),[§IV\-C](https://arxiv.org/html/2605.29021#S4.SS3.p1.2)\.
- \[4\]A\. Boonrath, F\. Rossi, I\. A\. Nesnas, and E\. M\. Botta\(2025\)A mission concept for removing multiple debris using formation\-flying slingshot spacecraft\.Acta Astronautica237,pp\. 443–459\.External Links:[Document](https://dx.doi.org/10.1016/j.actaastro.2025.08.035)Cited by:[§II\-B](https://arxiv.org/html/2605.29021#S2.SS2.p1.1),[§II\-B](https://arxiv.org/html/2605.29021#S2.SS2.p2.4)\.
- \[5\]E\. M\. Botta, I\. Sharf, A\. K\. Misra, and M\. Teichmann\(2016\)On the simulation of tether\-nets for space debris capture with Vortex Dynamics\.Acta Astronautica123,pp\. 91–102\.External Links:[Document](https://dx.doi.org/10.1016/j.actaastro.2016.02.012)Cited by:[§II\-A](https://arxiv.org/html/2605.29021#S2.SS1.p2.10)\.
- \[6\]E\. M\. Botta, I\. Sharf, and A\. K\. Misra\(2019\)Simulation of tether\-nets for capture of space debris and small asteroids\.Acta Astronautica155,pp\. 448–461\.External Links:[Document](https://dx.doi.org/10.1016/j.actaastro.2018.07.046)Cited by:[§II\-A](https://arxiv.org/html/2605.29021#S2.SS1.p2.10),[footnote 2](https://arxiv.org/html/2605.29021#footnote2)\.
- \[7\]E\. M\. Botta\(2017\)Deployment and Capture Dynamics of Tether\-Nets for Active Space Debris Removal\.Ph\.D\. Thesis,McGill University,Montreal, Québec\.Cited by:[§II\-C](https://arxiv.org/html/2605.29021#S2.SS3.p4.13)\.
- \[8\]S\. Chowdhury, W\. Tong, A\. Messac, and J\. Zhang\(2013\)A mixed\-discrete particle swarm optimization algorithm with explicit diversity\-preservation\.Structural and Multidisciplinary Optimization47,pp\. 367–388\.External Links:[Document](https://dx.doi.org/10.1007/s00158-012-0851-z)Cited by:[§III\-B](https://arxiv.org/html/2605.29021#S3.SS2.p2.1)\.
- \[9\]DuPont\(2023\)Kevlar® aramid fiber technical guide\.Note:Accessed: 2025\-03\-20External Links:[Link](https://www.dupont.com/products/kevlar)Cited by:[§II\-C](https://arxiv.org/html/2605.29021#S2.SS3.p4.13)\.
- \[10\]ECAPS\(2023\)ECAPS 5n hpgp thruster brochure\.Note:Accessed: 2025\-03\-20External Links:[Link](https://satsearch.co/products/ecaps-5n-hpgp-thruster)Cited by:[§II\-C](https://arxiv.org/html/2605.29021#S2.SS3.p4.13)\.
- \[11\]A\. N\. Hirani, K\. Kalyanaraman, and S\. Watts\(2015\)Graph laplacians and least squares on graphs\.In2015 IEEE International Parallel and Distributed Processing Symposium Workshop,Vol\.,pp\. 812–821\.External Links:[Document](https://dx.doi.org/10.1109/IPDPSW.2015.73)Cited by:[§III\-B](https://arxiv.org/html/2605.29021#S3.SS2.p2.1)\.
- \[12\]Y\. Hou, C\. Liu, H\. Hu, W\. Yang, and J\. Shi\(2021\)Dynamic computation of a tether\-net system capturing a space target via discrete elastic rods and an energy\-conserving integrator\.Acta Astronautica186,pp\. 118–134\.External Links:[Document](https://dx.doi.org/10.1016/j.actaastro.2021.05.029)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1)\.
- \[13\]W\. Huang, H\. Zou, H\. Liu, W\. Yang, J\. Gao, and Z\. Liu\(2023\)Contact dynamic analysis of tether\-net system for space debris capture using incremental potential formulation\.Advances in Space Research\.External Links:[Document](https://dx.doi.org/10.1016/j.asr.2023.05.054)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1)\.
- \[14\]M\. Jang, H\. Shin, C\. Sim, J\. Park, and H\. Cho\(2023\)Active debris removal simulations using spider\-web space\-nets for kitsat\-1 satellite\.International Journal of Aeronautical and Space Sciences24\(5\),pp\. 1311–1322\.External Links:[Document](https://dx.doi.org/10.1007/s42405-023-00623-2)Cited by:[§II\-C](https://arxiv.org/html/2605.29021#S2.SS3.p1.1)\.
- \[15\]K\. Lei, P\. Guo, Y\. Wang, X\. Wu, and W\. Zhao\(2022\)Solve routing problems with a residual edge\-graph attention neural network\.Neurocomputing508,pp\. 79–98\.External Links:ISSN 0925\-2312,[Document](https://dx.doi.org/10.1016/j.neucom.2022.08.005)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1)\.
- \[16\]F\. Liu, A\. Boonrath, E\. M\. Botta, and S\. Chowdhury\(2025\)Surrogate\-aided learning of active tether\-net maneuver to capture rotating space debris\.IEEE Transactions on Aerospace and Electronic Systems62,pp\. 630 – 645\.External Links:[Document](https://dx.doi.org/10.1109/TAES.2025.3624187)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1),[§II\-A](https://arxiv.org/html/2605.29021#S2.SS1.p2.10),[§IV\-A](https://arxiv.org/html/2605.29021#S4.SS1.p1.6)\.
- \[17\]F\. Liu, A\. Boonrath, P\. KrisshnaKumar, E\. M\. Botta, and S\. Chowdhury\(2023\)Learning constrained corner node trajectories of a tether net system for space debris capture\.InAIAA AVIATION 2023 Forum,pp\. 3920\.External Links:[Document](https://dx.doi.org/10.2514/6.2023-3920)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1)\.
- \[18\]F\. Liu, A\. Boonrath, G\. Madhu, E\. M\. Botta, and S\. Chowdhury\(2026\)Data and video repository for the paper – designing active tether\-net system for space debris capture with graph learning aided mixed\-combinatorial optimization\.Note:[https://github\.com/adamslab\-ub/GNN\-ReCo\-Benchmark/tree/GNN\_NavCo](https://github.com/adamslab-ub/GNN-ReCo-Benchmark/tree/GNN_NavCo)Cited by:[§IV\-C](https://arxiv.org/html/2605.29021#S4.SS3.p1.2)\.
- \[19\]F\. Liu, S\. Chowdhury, A\. Boonrath, and E\. M\. Botta\(2025\-08\)Efficient design optimization over mixed\-combinatorial spaces enabled by graph\-learning\.InASME 2025 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference,pp\. V03AT03A021\.External Links:[Document](https://dx.doi.org/10.1115/DETC2025-169719),https://asmedigitalcollection\.asme\.org/IDETC\-CIE/proceedings\-pdf/IDETC\-CIE2025/89220/V03AT03A021/7557501/v03at03a021\-detc2025\-169719\.pdfCited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1),[§I](https://arxiv.org/html/2605.29021#S1.p4.1),[§II\-C](https://arxiv.org/html/2605.29021#S2.SS3.p4.13),[TABLE II](https://arxiv.org/html/2605.29021#S2.T2.16)\.
- \[20\]Z\. Meng, P\. Huang, and J\. Guo\(2017\)Approach modeling and control of an autonomous maneuverable space net\.IEEE Transactions on Aerospace and Electronic Systems53\(6\),pp\. 2651–2661\.External Links:[Document](https://dx.doi.org/10.1109/TAES.2017.2709794)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1)\.
- \[21\]S\. Paul and S\. Chowdhury\(2024\-08\)Learning to allocate time\-bound and dynamic tasks to multiple robots using covariant attention neural networks\.Journal of Computing and Information Science in Engineering24\(9\),pp\. 091005\.External Links:[Document](https://dx.doi.org/10.1115/1.4065883),ISSN 1530\-9827Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1)\.
- \[22\]S\. Paul, P\. Ghassemi, and S\. Chowdhury\(2022\)Learning scalable policies over graphs for multi\-robot task allocation using capsule attention networks\.In2022 International Conference on Robotics and Automation \(ICRA\),Vol\.,pp\. 8815–8822\.External Links:[Document](https://dx.doi.org/10.1109/ICRA46639.2022.9812370)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1)\.
- \[23\]M\. Ringkamp, S\. Ober\-Blöbaum, and S\. Leyendecker\(2017\)On the time transformation of mixed integer optimal control problems using a consistent fixed integer control function\.Mathematical programming161\(1\-2\),pp\. 551–581\(eng\)\.External Links:ISSN 0025\-5610,[Document](https://dx.doi.org/10.1007/s10107-016-1023-5)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p2.1)\.
- \[24\]B\. Sanchez\-Lengeling, E\. Reif, A\. Pearce, and A\. B\. Wiltschko\(2021\)A gentle introduction to graph neural networks\.Distill\.Note:https://distill\.pub/2021/gnn\-introCited by:[§III\-A](https://arxiv.org/html/2605.29021#S3.SS1.p1.11)\.
- \[25\]M\. Shan, J\. Guo, and E\. Gill\(2015\)A review and comparison of active space debris capturing and removal methods\.Progress in Aerospace Sciences80,pp\. 18–32\.External Links:[Document](https://dx.doi.org/10.1016/j.paerosci.2015.11.001)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1)\.
- \[26\]R\. S\. Systems\(2023\)5N low throughput \(lt\) thruster\.Note:Accessed: 2025\-03\-20External Links:[Link](https://www.rubicon.space/product/19/5n-lt-thruster)Cited by:[§II\-C](https://arxiv.org/html/2605.29021#S2.SS3.p4.13)\.
- \[27\]W\. Tong, S\. Chowdhury, and A\. Messac\(2016\)A multi\-objective mixed\-discrete particle swarm optimization with multi\-domain diversity preservation\.Structural and Multidisciplinary Optimization53,pp\. 471–488\.External Links:[Document](https://dx.doi.org/10.1007/s00158-015-1319-8)Cited by:[§III\-B](https://arxiv.org/html/2605.29021#S3.SS2.p2.1)\.
- \[28\]A\. Vaswani, N\. Shazeer, N\. Parmar, J\. Uszkoreit, L\. Jones, A\. N\. Gomez, Ł\. Kaiser, and I\. Polosukhin\(2017\)Attention is all you need\.Advances in neural information processing systems30\.Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1)\.
- \[29\]P\. Veličković, G\. Cucurull, A\. Casanova, A\. Romero, P\. Liò, and Y\. Bengio\(2018\)Graph attention networks\.InInternational Conference on Learning Representations,Vol\.6\.External Links:[Document](https://dx.doi.org/)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1)\.
- \[30\]S\. Verma and Z\. Zhang\(2018\)Graph capsule convolutional neural networks\.InarXiv,Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1)\.
- \[31\]S\. V\. Weston, C\. D\. Burkhard, J\. M\. Stupl, R\. L\. Ticknor, B\. D\. Yost, R\. A\. Austin, P\. Galchenko, L\. K\. Newman, and L\. S\. Soto\(2025\-02\)State\-of\-the\-art small spacecraft technology 2024\.Technical reportTechnical ReportNASA/TP—20250000142,National Aeronautics and Space Administration \(NASA\),Ames Research Center, Moffett Field, California\.External Links:[Link](https://www.nasa.gov/smallsat-institute/sst-soa)Cited by:[§II\-C](https://arxiv.org/html/2605.29021#S2.SS3.p4.13)\.
- \[32\]Z\. Wu, S\. Pan, F\. Chen, G\. Long, C\. Zhang, and P\. S\. Yu\(2019\)A comprehensive survey on graph neural networks\.CoRRabs/1901\.00596\.External Links:[Document](https://dx.doi.org/10.1109/TNNLS.2020.2978386)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1)\.
- \[33\]Y\. Xia, X\. Chen, Z\. Liu, W\. Zhou, W\. Yao, and Z\. Zhang\(2025\)Mixed integer programming modeling for the satellite three\-dimensional component assignment and layout optimization problem\.Chinese Journal of Aeronautics,pp\. 103415\.External Links:ISSN 1000\-9361,[Document](https://dx.doi.org/10.1016/j.cja.2025.103415)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p2.1)\.
- \[34\]C\. Zeng, G\. R\. Hecht, S\. Chowdhury, and E\. M\. Botta\(0\)Concurrent design optimization of tether\-net system and actions for reliable space\-debris capture\.Journal of Spacecraft and Rockets0\(0\),pp\. 1–11\.External Links:[Document](https://dx.doi.org/10.2514/1.A35812),[Link](https://doi.org/10.2514/1.A35812)Cited by:[§II\-C](https://arxiv.org/html/2605.29021#S2.SS3.p1.1)\.
- \[35\]S\. Zhang, J\. S\. Campos, C\. Feldmann, F\. Sandfort, M\. Mathea, and R\. Misener\(2024\)Augmenting optimization\-based molecular design with graph neural networks\.Computers & chemical engineering186,pp\. 108684–\(eng\)\.External Links:ISSN 0098\-1354,[Document](https://dx.doi.org/10.1016/j.compchemeng.2024.108684)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1)\.
- \[36\]S\. Zhang, H\. Tong, J\. Xu, and R\. Maciejewski\(2019\-11\-10\)Graph convolutional networks: a comprehensive review\.Computational Social Networks6\(1\),pp\. 11\.External Links:[Document](https://dx.doi.org/10.1186/s40649-019-0069-y),ISSN 2197\-4314Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p3.1)\.
- \[37\]Y\. Zhao, F\. Zhang, P\. Huang, and X\. Liu\(2020\)Impulsive super\-twisting sliding mode control for space debris capturing via tethered space net robot\.IEEE transactions on industrial electronics \(1982\)67\(8\),pp\. 6874–6882\(eng\)\.External Links:ISSN 0278\-0046,[Document](https://dx.doi.org/10.1109/TIE.2019.2940002)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1)\.
- \[38\]W\. Zhu, Z\. Pang, Z\. Du, G\. Gao, and Z\. H\. Zhu\(2024\)Multi\-debris capture by tethered space net robot via redeployment and assembly\.Journal of Guidance, Control, and Dynamics,pp\. 1–18\.External Links:[Document](https://dx.doi.org/10.2514/1.G007908)Cited by:[§I](https://arxiv.org/html/2605.29021#S1.p1.1)\.