Paper deep dive
Multi-Agent Temporal Logic Planning via Penalty Functions and Block-Coordinate Optimization
Eleftherios E. Vlahakis, Arash Bahari Kordabad, Lars Lindemann, Pantelis Sopasakis, Sadegh Soudjani, Dimos V. Dimarogonas
Intelligence
Status: succeeded | Model: Gemma-4-26B-A4B | Prompt: intel-v1 | Confidence: 93%
Last extracted: 7/20/2026, 11:52:40 PM
Summary
The paper proposes a scalable framework for multi-agent planning under Signal Temporal Logic (STL) constraints. It formulates the problem as an optimization program and introduces a penalty-based unconstrained relaxation using smooth STL semantics. The problem is solved efficiently using a Block-Coordinate Gradient Descent (BCGD) method, where each block corresponds to a single agent's decision variables. A two-layer optimization scheme enforces feasibility by increasing the penalty parameter in an outer loop while performing inner BCGD updates. The approach mitigates the high computational complexity associated with collaborative multi-agent tasks.
Entities (6)
Relation Signals (5)
STL Planning Problem → solvedby → Block-Coordinate Gradient Descent
confidence 95% · ...introduce a penalty-based unconstrained relaxation that can be efficiently solved via a Block-Coordinate Gradient Descent (BCGD) method...
Multi-Agent Planning → uses → Signal Temporal Logic
confidence 95% · Multi-agent planning under Signal Temporal Logic (STL) is often hindered by collaborative tasks...
Penalty Function → enables → Unconstrained Relaxation
confidence 92% · ...introduce a penalty-based unconstrained relaxation that can be efficiently solved...
BCGD → handles → Multi-Agent System
confidence 90% · ...where each block corresponds to a single agent's decision variables, thereby mitigating complexity.
Smooth STL Semantics → defines → Robustness Function
confidence 88% · By utilizing a quadratic penalty function defined via smooth STL semantics...
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:Multi-agent planning under Signal Temporal Logic (STL) is often hindered by collaborative tasks that lead to computational challenges due to the inherent high dimensionality of the problem, preventing scalable synthesis with satisfaction guarantees. To address this, we formulate STL planning as an optimization program under multi-agent STL constraints and introduce a penalty-based unconstrained relaxation that can be efficiently solved via a Block-Coordinate Gradient Descent (BCGD) method, where each block corresponds to a single agent's decision variables, thereby mitigating complexity. By utilizing a quadratic penalty function defined via smooth STL semantics, we show that BCGD iterations converge to a stationary point of the penalized problem under standard regularity assumptions. To enforce feasibility, the BCGD solver is embedded within a two-layer optimization scheme: inner BCGD updates are performed for a fixed penalty parameter, which is then increased in an outer loop to progressively improve multi-agent STL robustness. The proposed framework enables scalable computations and is validated through various complex multi-robot planning scenarios.
Tags
Links
- Source: https://arxiv.org/abs/2602.17434v2
- Canonical: https://arxiv.org/abs/2602.17434v2
Trouble viewing inline? Open PDF directly →
Full Text
47,278 characters extracted from source content.
Expand or collapse full text
Multi-Agent Temporal Logic Planning via Penalty Functions and Block-Coordinate Optimization Eleftherios E. Vlahakis1 Arash Bahari Kordabad2 Lars Lindemann3 Pantelis Sopasakis4 Sadegh Soudjani2 and Dimos V. Dimarogonas1 This work was supported by the Swedish Research Council (VR), the Knut & Alice Wallenberg Foundation (KAW), the Horizon Europe Grant SymAware and the ERC Consolidator Grant LEAFHOUND.1 Division of Decision and Control Systems, School of Electrical Engineering and Computer Science, KTH Royal Institute of Technology, 10044, Stockholm, Sweden; vlahakis,dimos@kth.se2 Max Planck Institute for Software Systems, Kaiserslautern, Germany; arashbk, sadegh@mpi-sws.org3 Automatic Control Laboratory, ETH Zürich, Zürich, 8092, Switzerland; llindemann@ethz.ch4 School of Electronics, Electrical Engineering and Computer Science, Queen’s University Belfast, Northern Ireland, BT9 5BN, UK; p.sopasakis@qub.ac.uk Abstract Multi-agent planning under Signal Temporal Logic (STL) is often hindered by collaborative tasks that lead to computational challenges due to the inherent high-dimensionality of the problem, preventing scalable synthesis with satisfaction guarantees. To address this, we formulate STL planning as an optimization program under arbitrary multi-agent constraints and introduce a penalty-based unconstrained relaxation that can be efficiently solved via a Block-Coordinate Gradient Descent (BCGD) method, where each block corresponds to a single agent’s decision variables, thereby mitigating complexity. By utilizing a quadratic penalty function defined via smooth STL semantics, we show that BCGD iterations converge to a stationary point of the penalized problem under standard regularity assumptions. To enforce feasibility, the BCGD solver is embedded within a two-layer optimization scheme: inner BCGD updates are performed for a fixed penalty parameter, which is then increased in an outer loop to progressively improve multi-agent STL robustness. The proposed framework enables scalable computations and is validated through various complex multi-robot planning scenarios. IEEEkeywords Block-coordinate optimization, multi-agent systems, penalty methods, signal temporal logic 1 Introduction Multi-agent systems (MAS) research deals with the task of coordinating collections of autonomous systems, e.g., in logistics, exploration, and smart infrastructure. These applications require agents to satisfy complex spatio-temporal and logical constraints governing both individual and interactive behaviors. Signal Temporal Logic (STL) [10] has emerged as a powerful framework suitable for these requirements, offering an expressive language to encode time-bounded properties over continuous-time signals. Unlike automata-based LTL synthesis [1], STL’s quantitative semantics [4, 3] enable the direct optimization of satisfaction margins over system trajectories. While exact solutions can be obtained via mixed-integer MPC formulations [15, 8], their poor scalability has motivated the development of smooth robustness relaxations [14, 11, 5, 6] for efficient gradient-based optimization. However, in MAS settings, collaborative tasks further amplify computational complexity, rendering scalable planning formidable. Important existing multi-agent approaches, including distributed MPC [2, 22], decentralized feedback control [9], and sequential planning [21], address these challenges but are often limited to restricted STL fragments, rely on the assumption of feasible solutions, or employ heuristic coordination schemes that lack formal guarantees. Developing scalable planning methods that can handle complex collaborative specifications while ensuring rigorous satisfaction thus remains an important open problem. To address this challenge, this paper introduces a scalable, optimization-based framework for multi-agent STL planning under arbitrary collaborative tasks, leveraging smooth STL semantics [5] and the computational efficiency of the Block-Coordinate Gradient Descent (BCGD) method [18]. Specifically, we demonstrate that the original planning problem, featuring an objective function that is separable across agents, yet subject to generally coupled multi-agent STL constraints, can be relaxed into an unconstrained problem via a quadratic penalty defined over smooth robustness metrics. This relaxation is solved efficiently via BCGD, where computations are performed at the block level, with each block corresponding to the decision variables of a single agent. Under standard regularity assumptions, we show that the BCGD iterations converge to a stationary point of the penalized problem, providing a computational architecture that remains invariant to the complexity of the multi-agent specification. To enforce feasible solutions for the original planning problem, the BCGD solver is embedded in a two-layer optimization scheme, forming a penalty method (PM) [12, Chap. 17], where the inner loop optimizes for a fixed penalty parameter, which is then updated in the outer loop to progressively improve multi-agent STL robustness. This framework is the first to systematically integrate smooth STL semantics, block-coordinate optimization, and penalty functions for efficient multi-agent STL synthesis. We validate BCGD-PM across complex multi-robot scenarios, benchmarking it against an LBFGS-based implementation [12, Chap. 9] within the same modular penalty framework. This comparison highlights BCGD’s advantages in handling high-dimensional STL planning. For readability, the main technical proofs are provided in the Appendix. 2 Problem setup Notation: The sets of real numbers and nonnegative integers are IR I\!R and IN I\!N, respectively. Let N∈IN∈ I\!N so that IN[0,N]=0,1,…,N I\!N_[0,N]=\0,1,…,N\. Let x1,…,xnx_1,…,x_n be vectors so that x=(x1,…,xn)=[x1⊺⋯xn⊺]⊺x=(x_1,…,x_n)=[x_1 \;·s\;x_n ] . Let a,b∈INa,b∈ I\!N with a≤ba≤ b. We denote by (a:b)=(x(a),…,x(b)) x(a:b)=(x(a),…,x(b)) an aggregate vector consisting of x(t)x(t), t∈IN[a,b]t∈ I\!N_[a,b], representing a trajectory, or an aggregate trajectory when x(t)=(x1(t),…,xM(t))x(t)=(x_1(t),…,x_M(t)) with M∈INM∈ I\!N. The cardinality of a set V is |||V|. We call x⋆x a stationary point of f:IRn→IRf\,:\, I\!R^n→ I\!R if its gradient at x⋆x ∇f(x⋆)=0∇ f(x )=0. Given a sequence xk\x^k\, x¯ x is called an accumulation point if there exists a subsequence xkj\x^k_j\ such that xkj→x¯x^k_j→ x. 2.1 Signal temporal logic We consider STL formulas in positive normal form (PNF) with syntax φ≔⊤∣π∣¬π∣ϕ1∧ϕ2∣ϕ1∨ϕ2∣□ℐϕ1∣ϕ1ℐϕ2, π π _1 _2 _1 _2 _I _1 _1\,U_I _2, (1) where π≔(μ(x)≥0)π (μ(x)≥ 0) is a predicate, with predicate function μ:IRnx→IRμ: I\!R^n_x→ I\!R, ϕ1 _1 and ϕ2 _2 are STL formulas built recursively using the grammar in (1), ¬ , ∧ , and ∨ are the logical operators denoting negation, conjunction, and disjunction, respectively, and □ℐ _I and ℐU_I are the always and until temporal operators, respectively, defined over the discrete interval ℐ⊂INI⊂ I\!N. We omit the eventually operator (◆ℐ _I) from (1) since ◆ℐϕ1=⊤ℐϕ1 _I _1= _I _1. The above definition has negation appearing only beside atomic predicates. This form of STL specifications in PNF is equivalent to the full class of STL specifications [10], and any STL formula can be transformed into a PNF using usual logical identities [16, Prop. 2]. We denote by (t)⊧ϕ x(t) φ, t∈INt∈ I\!N, the satisfaction of ϕφ, verified over (t)=(x(t),x(t+1),…) x(t)=(x(t),x(t+1),…). The validity of ϕφ can be determined recursively using the Boolean semantics of STL; for details, we refer to [10] due to space limitations. STL is endowed with quantitative semantics [3]: A scalar-valued function ρϕ:×t=0ℋϕIRn→IRρ^φ: _t=0^H^φ I\!R^n→ I\!R of a signal (t) x(t), termed robustness function, where ℋϕH^φ is the horizon of ϕφ [10], indicates how robustly a signal (t) x(t) satisfies a formula ϕφ, and is defined recursively as ρπ((t)) ρ^π( x(t)) =μ(x(t)), =μ(x(t)), ρ¬π((t)) ρ π( x(t)) =−ρπ((t)), =-ρ^π( x(t)), ρϕ1∧ϕ2((t)) ρ _1 _2( x(t)) =min(ρϕ1((t)),ρϕ2((t))), = \! (ρ _1( x(t)),\,ρ _2( x(t)) ), ρϕ1∨ϕ2((t)) ρ _1 _2( x(t)) =max(ρϕ1((t)),ρϕ2((t))), = \! (ρ _1( x(t)),\,ρ _2( x(t)) ), ρ□ℐϕ1((t)) ρ _I _1( x(t)) =minτ∈t⊕ℐρϕ1((τ)), = _τ∈ t ρ _1( x(τ)), ρϕ1ℐϕ2((t)) ρ _1U_I _2 -2.0mu( x(t)) =maxτ∈t⊕ℐ(min(ρϕ1((τ)),minτ′∈IN[t,τ]ρϕ2((τ′)))), = -8.0mu _ subarrayc\\[0.1pt] τ∈ t subarray -7.0mu ( -7.0mu \! (ρ _1 -2.0mu( x(τ)), -17.0mu _ subarrayc\\[0.1pt] τ ∈ I\!N_[t,τ] subarray -17.0muρ _2 -2.0mu( x(τ )) ) -6.0mu ) -3.0mu, (2) where ⊕ denotes the Minkowski sum. The satisfaction of a formula ϕφ by a signal (t) x(t) is indicated by the sign of its robustness function in the sense that ρϕ((t))>0⇒(t)⊧ϕρ^φ( x(t))>0 x(t) φ, ρϕ((t))<0⇒(t)⊧̸ϕρ^φ( x(t))<0 x(t) φ, while ρϕ((t))=0ρ^φ( x(t))=0 does not in general determine satisfiability. 2.2 Multi-agent system 2.2.1 Dynamics We consider a MAS with M agents, with the ithi^th agent following the dynamics xi(t+1)=fi(xi(t),ui(t)),x_i(t+1)=f_i(x_i(t),u_i(t)), (3) where xi(t)∈IRnix_i(t)∈ I\!R^n_i and ui∈IRmiu_i∈ I\!R^m_i are the state and input vectors, respectively, fi:IRni×IRmi→IRnif_i: I\!R^n_i× I\!R^m_i→ I\!R^n_i is continuously differentiable and locally Lipschitz in (xi,ui)(x_i,u_i), t∈INt∈ I\!N, and the initial condition xi(0)x_i(0) is known, with i∈=1,…,Mi =\1,…,M\, where V is the set of all agents. To formally group agents participating in collaborative tasks, we adopt the notion of cliques from graph theory, defined next [13]. Definition 1 Consider an undirected graph =(,ℰ)G=(V,E) potentially containing self-loops and multiple edges with node set V and edge set ℰE. Let ν⊆ν and ℰν⊆ℰE_ν be the set of edges connecting the nodes in ν. Then, (ν,ℰν)(ν,E_ν) is called a clique if ∀νi,νj∈ν∀\, _i, _j∈ν, (νi,νj)∈ℰν( _i, _j) _ν, i.e., it is a complete subgraph of G. Consider a graph =(,ℰ)G=(V,E) with clique set K. With a slight abuse of notation, we denote a clique (ν,ℰν)∈(ν,E_ν) simply by ν, and write |ν||ν| for the number of nodes in ν. Let ν∈ν contain the agents i1,…,i|ν|i_1,…,i_|ν|, i.e, ν=i1,…,i|ν|ν=\i_1,…,i_|ν|\. By collecting individual state and input vectors, as xν(t)=(xi1(t),…,xi|ν|(t))∈IRnνx_ν(t)=(x_i_1(t),…,x_i_|ν|(t))∈ I\!R^n_ν and uν(t)=(ui1(t),…,ui|ν|(t))∈IRmνu_ν(t)=(u_i_1(t),…,u_i_|ν|(t))∈ I\!R^m_ν, respectively, we write the aggregate dynamics of |ν||ν| agents as xν(t+1)=fν(xν(t),uν(t)),x_ν(t+1)=f_ν(x_ν(t),u_ν(t)), (4) with fν(xν(t),uν(t))=(fj(xj(t),uj(t)))j∈νf_ν(x_ν(t),u_ν(t))= (f_j(x_j(t),u_j(t)) )_j∈ν. When ν=ν=V, the aggregate dynamics of the entire MAS are written as x(t+1)=f(x(t),u(t)),x(t+1)=f(x(t),u(t)), (5) with x(t)=(x1(t),…,xM(t))x(t)=(x_1(t),…,x_M(t)), u(t)=(u1(t),…,uM(t))u(t)=(u_1(t),…,u_M(t)), and f=(f1,…,fM)f=(f_1,…,f_M). 2.2.2 STL specification The MAS is subject to ϕ=⋀ν∈ϕν,φ= _ν _φ _ν, (6) which is a conjunctive STL formula, where each conjunct ϕν _ν is defined over ν(t) x_ν(t) and follows the syntax in (1), with the aggregate trajectory ν(t) x_ν(t) collecting the individual trajectories of the agents in the clique ν. The set ϕK_φ collects all these cliques induced by ϕφ, and may include individual agents (|ν|=1)(|ν|=1) or group of agents (1<|ν|≤||)(1<|ν|≤|V|). Note that different cliques may overlap in their agent sets, indicating that some agents participate in multiple collaborative tasks. Let π≔(μ(y)≥0)π (μ(y)≥0) be a predicate in ϕφ, where μ:IRny→IRμ: I\!R^n_y→ I\!R. We assume that all predicate functions μ appearing in ϕφ are continuously differentiable. The vector y∈IRnyy∈ I\!R^n_y may represent an individual state xi∈IRnix_i∈ I\!R^n_i, for i∈i∈V, an aggregate state xν∈IRnνx_ν∈ I\!R^n_ν for a clique ν∈ϕν∈K_φ, or an aggregate state xκ∈IRnκx_κ∈ I\!R^n_κ collecting the states of a subset of agents in a clique ν. For example, in formula ϕν=πνκ1ℐπνκ2 _ν= _ν _1U_I _ν _2, where ν=1,2,3ν=\1,2,3\, the predicate πνκ1≔(x1≥0) _ν _1 (x_1≥0) involves one agent κ1=1 _1=\1\, while πνκ2≔(x2−x3≥0) _ν _2 (x_2-x_3≥0) involves two agents κ2=2,3 _2=\2,3\. 2.3 Problem statement The multi-agent STL planning synthesis problem is formulated as an optimization problem over the multi-agent control sequence =(u(0),…,u(N−1)) u=(u(0),…,u(N-1)), where u(t)=(u1(t),…,uM(t))u(t)=(u_1(t),…,u_M(t)), and N=ℋϕN=H^φ. Given the initial condition x(0)=x0x(0)=x_0 and the dynamics x(t+1)=f(x(t),u(t))x(t+1)=f(x(t),u(t)), the multi-agent trajectory =(x(0),…,x(N)) x=(x(0),…,x(N)), where x(t)=(x1(t),…,xM(t))x(t)=(x_1(t),…,x_M(t)), is explicitly determined by u. We thus denote the multi-agent STL robustness ρϕ()ρ^φ( x) as ρϕ()ρ^φ( u), and the total cost as ℒ()=∑i∈ℒi(i)L( u)= _i L_i( u_i), where ℒi(i)=∑t=0N−1ℓi(xi(t),ui(t),t)+Vf,i(xi(N))L_i( u_i)= _t=0^N-1 _i(x_i(t),u_i(t),t)+V_f,i(x_i(N)), with i=(ui(0),…,ui(N−1)) u_i=(u_i(0),…,u_i(N-1)) denoting the ithi^th agent’s inputs. Function ℓi:IRni×IRmi×IN→IR _i: I\!R^n_i× I\!R^m_i× I\!N→ I\!R is the running cost for the ithi^th agent, penalizing, e.g., state energy and control effort quantities that cannot be directly encoded through the specification ϕφ. The terminal cost Vf,i:IRni→IRV_f,i: I\!R^n_i→ I\!R penalizes deviations from a prescribed terminal condition, which can be chosen, e.g., to enforce horizon-end requirements in ϕφ or to induce cyclic trajectory planning by penalizing the distance xi(N)−xi(0)x_i(N)-x_i(0) for all i∈i . The multi-agent STL planning problem is formulated as Minimizeℒ()=∑i∈ℒi(i) *Minimize_ u\,\,L( u)= _i L_i( u_i) (7a) subject toρϕ()=minν∈ϕρϕν(ν)>0. \;to\,\,ρ^φ( u)= _ν _φρ _ν( u_ν)>0. (7b) Solving (7) is computationally challenging due to the non-smoothness and coupling in the joint STL constraint in (7b), which implies that ν⊧ϕν x_ν _ν, ∀ν∈ϕ∀ν∈K_φ, i.e., ⊧ϕ x φ. We address this challenge by employing a smooth STL approximation and a penalty-based block-coordinate gradient descent method, which enables agent-level computations, mitigating the complexity of the centralized problem. This approach relies on the following assumption. Assumption 1 The function ℒi:IRni→IRL_i: I\!R^n_ u_i→ I\!R, i∈i , is proper, convex, continuous, and level-bounded. Furthermore, the problem in (7) is feasible. In the remainder of the paper, Assumption 1 holds. We note that the first part holds for linear dynamics, or, more generally, when the running cost ℓi _i depends only on the input variables and when the terminal cost Vf,iV_f,i is omitted. To streamline the presentation of the proposed optimization framework, we omit explicit input constraints from (7). 3 Multi-agent STL optimization 3.1 Smooth STL approximations We recall that the STL constraint in (7b) is non-smooth involving min and max operations over predicate functions. We underapproximate min and max operators [5] as min(μ1,…,μq)≈−1Γlog(∑j=1qe−Γμj), ( _1,…, _q ) ^≥≈- 1 ( _j=1^qe^- _j ), (8a) max(μ1,…,μq)≈∑j=1qμjeΓμj∑j=1qeΓμj, ( _1,…, _q ) ^≥≈ _j=1^q _je _j _j=1^qe _j, (8b) where the approximation becomes tighter for larger values of Γ>0 >0. For a fixed Γ , we denote by ϱΓϕ() ^φ_ ( u) the smooth approximation of the robustness function ρϕ()ρ^φ( u), which satisfies ϱΓϕ()≤ρϕ() ^φ_ ( u)≤ρ^φ( u) for all Γ>0 >0, since ϕφ is in PNF [7]. This implies that the STL constraint in (7b), can be replaced by ϱΓϕ()>0 ^φ_ ( u)>0, introducing additional conservatism in satisfying ⊧ϕ x φ controlled by Γ . In fact, due to the multi-agent structure of ϕφ in (6), the strict requirement ϱΓϕ()>0 ^φ_ ( u)>0 can be relaxed to ϱΓϕ()≥0 ^φ_ ( u)≥0, as formally stated next. Proposition 1 Consider the formula ϕφ in (6), where ϕK_φ contains at least two cliques, i.e., |ϕ|≥2|K_φ|≥ 2. Then, for all multi-agent sequences u and Γ>0 >0, the smooth robustness satisfies ϱΓϕ()<ρϕ() _ ^φ( u)<ρ^φ( u). Consequently, ϱΓϕ()≥0⇒ρϕ()>0. _ ^φ( u)≥0 ρ^φ( u)>0. 3.2 Unconstrained STL optimization Due to Proposition 1, a stricter version of the problem in (7) can now be written using smooth STL semantics as Minimize∈ℳℒ(),whereℳ≔∣ϱΓϕ()≥0. *Minimize_ u \,L( u),\;where\;M \ u _ ^φ( u)≥ 0\. (9) By Assumption 1, this problem has a minimizer [12]. We address this constrained program using a suitable penalty function R()R( u) leading to an unconstrained problem: MinimizeFλ()≔ℒ()+λR(), *Minimize_ u\,F_λ( u) ( u)+λ R( u), (10) for some λ>0λ>0. To specify a suitable penalty λR()λ R( u) that relaxes the original constrained formulation into the unconstrained problem in (10) we introduce the quadratic penalty R()≔max0,−ϱΓϕ()2, R( u) \0,- _ ^φ( u)\^2, (11) which satisfies ℳ=∣R()=0M=\ u R( u)=0\, and R()>0R( u)>0 for all ∉ℳ u∉M. This is a differentiable penalty function with gradient ∇R()=−2max0,−ϱΓϕ()∇ϱΓϕ(),∇ R( u)=-2 \0,- _ ^φ( u)\∇ _ ^φ( u), (12) which exists everywhere and is continuous. We shall first discuss how the unconstrained problem in (10) can be solved for a fixed value of λ>0λ>0 using a block coordinate gradient descent method. Subsequently, we will show how a penalty method can be used to solve the constrained problem in (9). 3.3 Block coordinate gradient descent Problem (10) possesses a certain structure that can be exploited to solve it efficiently. Firstly, note that the cost function ℒ()L( u) is separable into agent-level objective functions ℒi(i)L_i( u_i), which by Assumption 1 are proper, convex, and continuous in i u_i for all i∈i . The penalty term λR()λ R( u), which penalizes violations of the smooth robustness condition ϱΓϕ()≥0 _ ^φ( u)≥ 0, is differentiable. These properties allow us to use the block coordinate gradient descent (BCGD) method [18]. In each BCGD iteration, (i) λR()λ R( u) is approximated by a strictly convex quadratic function, enabling the application of block coordinate descent to generate a descent direction, and (i) a sufficient descent direction is computed for an agent-level block of u, as detailed in Alg. 1. Quadratic approximation: We model the variation R(k+)−R(k)R( u^k+ d)-R( u^k) at k u^k by QH(k,)≔∇R(k)⊺+12⊺HkQ^H( u^k, d) ∇ R( u^k) d+ 12 d H^k d, where d is the update direction to k u^k derived from (13), and Hk≻0H^k 0 approximates the Hessian ∇2R(k)∇^2R( u^k). Block selection: In iteration k we choose a nonempty subset of agent-level blocks Jk⊆J^k , so that only the elements (i)i∈Jk( u_i)_i∈ J^k will be updated. Over iterations, all blocks should be updated at least once (generalized Gauss–Seidel), or one may select the most “active” blocks using a Gauss–Southwell type rule for efficiency—see [19, 18] for details. Block update direction: We compute the descent direction k=argminλQH(k,)+ℒ(k+)|j=0,∀j∉Jk. d^k -2.0mu= -2.0mu *argmin_ d -4.0mu \λ Q^H -3.0mu( u^k -3.0mu, d)+L( u^k+ d) | d_j=0,∀ j∉J^k -2.0mu \. (13) Due to the separability of ℒ()L( u), this step decomposes into independent subproblems for each selected agent block when HkH^k is chosen block-diagonal, i.e., it is equivalent to solving Min.(j)j∈Jk∑j∈Jkλ2j⊺Hjkj+λ∇R(k)j⊺j+ℒj(jk+j). *Min._( d_j)_j∈ J^k -3.0mu _j∈ J^k -3.0mu λ2 d_j H_j^k d_j+λ∇ R( u^k)_j d_j+L_j( u_j^k+ d_j). (14) Thus, for j∈Jkj∈ J^k, we have the optimization problems Min.jλ2j⊺Hjkj+λ∇R(k)j⊺j+ℒj(jk+j), *Min._ d_j\, λ2 d_j H_j^k d_j+λ∇ R( u^k)_j d_j+L_j( u_j^k+ d_j), (15) which are convex and can be solved efficiently, where Hjk≻0H_j^k 0 are the diagonal blocks of HkH^k, enabling a parallel computation across all agents in JkJ^k, when |Jk|≥1|J^k|≥1. Matrix HkH^k approximates the Hessian ∇2R(k)∇^2R( u^k), which is costly to compute. In fact, any Hjk≻0H_j^k 0 may be used, with sufficient decrease along jk d_j^k ensured via the inexact line search below. Armijo-type rule for step-size selection: Let σ,γ∈(0,1)σ,γ∈(0,1), ΔFλk=Fλ(k+αkk)−Fλ(k) F_λ^k=F_λ( u^k+α^k d^k)-F_λ( u^k), and Δℒk=ℒ(k+k)−ℒ(k) ^k=L( u^k+ d^k)-L( u^k). We select the largest αk∈1/2jj∈INα^k∈\1/2^j\_j∈ I\!N to ensure sufficient decrease of Fλ()F_λ( u) along k d^k, i.e., ΔFλk≤σαk(λ∇R(k)⊺k+γ(k)⊺Hkk+Δℒk). F_λ^k≤σα^k 0.0mu (λ∇ R( u^k) d^k+γ( d^k) H^k d^k+ ^k ). (16) Initialization rule: Let i=ν∈ϕ∣i∈νT_i=\ν _φ i∈ν\ be the set of cliques containing i. Assuming i∈i\i\∈T_i for all i∈i∈V, implying that for all i∈i∈V there is an individual task ϕi _i, one can obtain initial guesses for i0 u_i^0 by minimizing ℒi(i)+λmax0,−ϱΓϕi(i)2L_i( u_i)+λ\, \0,\,- _ _i( u_i)\^2 over i u_i, for i∈i∈V, using, e.g., the toolbox in [17]. Otherwise, one may start with 0=0 u^0=0. Termination: Algorithm 1 terminates when ‖∇Fλ()‖≤ϵ\|∇ F_λ( u^k)\|≤ε. Algorithm 1 Block Coordinate Gradient Descent (BCGD) 1:Input: 0 u^0 (initial iterate), λ (penalty), K (maximum iterations), ϵ>0ε>0 (tolerance) 2:for k=0,1,…,K−1k=0,1,…,K-1 do 3: If ‖∇Fλ(k)‖≤ϵ\|∇ F_λ( u^k)\|≤ε break 4: Select agent blocks Jk⊆J^k (Gauss–Seidel or Gauss–Southwell) 5: Compute ∇R(k)∇ R( u^k) and select HkH^k 6: Compute block update direction k d^k via (13) 7: Choose step size αkα^k using the Armijo rule 8: Update iterate: k+1=k+αkk u^k+1= u^k+α^k d^k 9:return ⋆,k u ,k Proposition 2 For a fixed penalty parameter λ>0λ>0, let k\ u^k\ be the sequence generated by Alg. 1 applied to (10). Then every accumulation point of k\ u^k\ is a stationary point of the penalized problem (10). The proof of Proposition 2 follows directly from [18, Theorem 1] tailored to the multi-agent STL setting considered here. For fixed λ, the penalized objective FλF_λ is continuously differentiable, consisting of a block-separable cost and a quadratic penalty function. Consequently, the assumptions of [18] are satisfied, and the proposed block-coordinate scheme converges to a stationary point of FλF_λ. We note that the most computationally demanding part of Alg. 1 is the computation of the gradient ∇R(k)∇ R( u^k), followed by the computation of the block update direction k d^k. To enforce the satisfaction of the constraints in (9) in the next section we will embed the BCGD algorithm within a penalty method. 3.4 Penalty method for the constrained problem We can solve the constrained problem (9) with the penalty method [12, Chap. 17]. Specifically, we solve a sequence of unconstrained problems of the form (10), each corresponding to a fixed penalty parameter λ>0λ>0, and iteratively increase λ to progressively enforce feasibility. Following [17, Alg. 1], the quadratic penalty method is given in Alg. 2. Upon termination, Alg. 2 yields an approximate solution ¯⋆ u which is ϵinfeas _ infeas-infeasible satisfying ϱΓϕ(¯⋆)>−ϵinfeas _ ^φ( u )>- _ infeas. In fact, under certain conditions, as stated in Theorem 1, a near-infeasible solution of the penalized problem (10) satisfies the original STL constraints (7b). As a result, one could use the original STL constraints in (7b) as a termination criterion in Alg. 2. Algorithm 2 Penalty method (PM) 1:Input: 0 u^0 (initial iterate), λ0λ^0 (initial penalty), KPMK_ PM (maximum iterations), ηλ>1 _λ>1 (penalty update factor), ϵinfeas>0 _ infeas>0 (max. infeasibility), ϵ0ε^0 (initial tolerance), ηϵ∈(0,1] _ε∈(0,1] (tolerance update parameter) 2:¯0=0 u^0= u^0 3:for k=0,1,…,KPM−1k=0,1,…,K_ PM-1 do 4: if R(¯k)<ϵinfeasR( u^k)< _ infeas break 5: Use Alg. 1 to solve (10) with λ=λkλ=λ^k, initial guess ¯k u^k, and tolerance ϵkε^k; obtain ϵkε^k-approx. solution ⋆,k u ,k 6: λk+1=ηλλkλ^k+1= _λ^k, ¯k+1=⋆,k u^k+1= u ,k, ϵk+1=ηϵϵkε^k+1= _ε^k 7:Output: ¯⋆,k u ,k Theorem 1 Suppose |ϕ|≥2|K_φ|≥2 and the sequence ¯k\ u^k\ of Alg. 2 has a limit point ¯∗ u^*, with R(¯∗)=0R( u^*)=0. For a given problem of the form (9), there is an ϵinfeas∗ _ infeas^* such that Alg. 2 with ϵinfeas≤ϵinfeas∗ _ infeas≤ _ infeas^* returns a feasible solution for (7). As a remark, as it is the case with the penalty method [12], Alg. 2 can get stuck in local minima. However, in Section 4 we demonstrate the effectiveness of this approach, which works well for multi-robot STL planning. Moreover, the “inner” optimization problems are of the form (10) and are solved efficiently with Alg. 1. The penalty method requires only a small number of “outer” iterations to converge. 4 Numerical validation via ten-robot example 05510101515202025253030353540404545055101015152020□ℐ¬1 _I _1□ℐ¬2 _I _2□ℐ¬3 _I _3◆ℐcC1 _I_cC_1◆ℐdD1 _I_dD_1◆ℐcC2 _I_cC_2◆ℐdD2 _I_dD_2◆ℐcC3 _I_cC_3◆ℐdD3 _I_dD_3◆ℐcC4 _I_cC_4◆ℐdD4 _I_dD_4◆ℐcC5 _I_cC_5◆ℐdD5 _I_dD_5◆ℐcC6 _I_cC_6◆ℐdD6 _I_dD_6◆ℐcC7 _I_cC_7◆ℐdD7 _I_dD_7◆ℐcC8 _I_cC_8◆ℐdD8 _I_dD_8◆ℐcC9 _I_cC_9◆ℐdD9 _I_dD_9◆ℐcC10 _I_cC_10◆ℐdD10 _I_dD_10◆ℐmM123 _I_mM_123◆ℐmM15 _I_mM_15◆ℐmM34 _I_mM_34◆ℐmM45 _I_mM_45◆ℐmM47 _I_mM_47◆ℐmM56 _I_mM_56◆ℐmM68 _I_mM_68◆ℐmM69 _I_mM_69◆ℐmM78 _I_mM_78◆ℐmM810 _I_mM_810◆ℐmM910 _I_mM_910i(0:N) z_i(0:N)i(0:N) y_i(0:N)x1(0)x_1(0)x2(0)x_2(0)x3(0)x_3(0)x4(0)x_4(0)x5(0)x_5(0)x6(0)x_6(0)x7(0)x_7(0)x8(0)x_8(0)x9(0)x_9(0)x10(0)x_10(0)x1(100)\;\;x_1(100)x2(100)\;\;x_2(100)x3(100)\;\;x_3(100)x4(100)\;\;x_4(100)x5(100)\;\;x_5(100)x6(100)\;\;x_6(100)x7(100)\;\;x_7(100)x8(100)\;\;x_8(100)x9(100)\;\;x_9(100)x10(100)\;\;x_10(100) Figure 1: RURAMCA motion planning for ten robots with discrete-time unicycle dynamics: zi(t+1)=zi(t)+vi(t)cosθi(t)z_i(t+1)=z_i(t)+v_i(t) _i(t), yi(t+1)=yi(t)+vi(t)sinθi(t)y_i(t+1)=y_i(t)+v_i(t) _i(t), and θi(t+1)=θi(t)+ωi(t) _i(t+1)= _i(t)+ _i(t), where (zi(t),yi(t))∈IR2(z_i(t),y_i(t))∈ I\!R^2 denote the Cartesian position of the ithi^th robot, θi(t)∈IR _i(t)∈ I\!R its heading, and vi(t)∈IRv_i(t)∈ I\!R, ωi(t)∈IR _i(t)∈ I\!R the linear and angular control inputs, respectively, with a sampling period of 11 s. The trajectories are generated from initial conditions (colored squares) using input sequences obtained via Alg. 2 (using Alg. 1 as inner solver). Simulation parameters: Γ=2 =2 for inner smooth operators, Γ=1 =1 for the outer softmin in (7b), σ=0.5σ=0.5 and γ=0.995γ=0.995 for the Armijo rule, 0=0 u^0=0, λ0=1λ^0=1, ηλ=5 _λ=5, ϵinfeas=5.0 10−4 _ infeas=5.0\,10^-4, ϵ=10−3ε=10^-3, and Hk=103IH^k=10^3I. The running cost functions are selected as ℓi=vi(t)2+ωi(t)2 _i=v_i(t)^2+ _i(t)^2, ∀t∈ℐ∀ t , and the terminal cost Vf,i=0V_f,i=0. We evaluate the proposed framework on the ten-robot motion planning problem introduced in [21, 20], considering a multi-agent specification ϕ=⋀ν∈ϕνφ= _ν _φ _ν with horizon ℋϕ=100H^φ=100 over the time window ℐ=IN[0,100]I= I\!N_[0,100]. Formulas ϕν _ν represent either individual or collaborative STL tasks, while ϕK_φ is the set of 21 cliques. A task ϕν _ν, with ν=iν=i, requires the ithi^th robot avoid obstacles at all times (□ℐ¬l _I _l), visit the ithi^th collection region (◆ℐcCi _I_cC_i) within ℐc=IN[10,50]I_c= I\!N_[10,50], and reach the ithi^th delivery region (◆ℐdDi _I_dD_i) within ℐd=IN[70,100]I_d= I\!N_[70,100]. The workspace, including collection and delivery regions for the ten robots, and three obstacles, is in Fig. 1. A task ϕν _ν, with |ν|>1|ν|>1, requires the robots forming the clique ν approach each other (◆ℐmMν _I_mM_ν) within ℐm=IN[0,70]I_m= I\!N_[0,70]. The formulas CiC_i, DiD_i, i∈i∈V, and lO_l, l∈IN[1,3]l∈ I\!N_[1,3], are conjunctions of predicates, each defined by a linear predicate function, while Mν≔(minκ∈νμκ(xκ)≥0)M_ν ( _κ∈ν _κ(x_κ)≥0), where xκ=(zq,yq)q∈κx_κ=(z_q,y_q)_q∈κ, with κ=(i,j)⊆νκ=(i,j) ν being any subset of two robots in ν, and μκ(xκ)=0.25−‖[I−I]xκ‖ _κ(x_κ)=0.25- \|[I\;-I]\,x_κ \|. We call this scenario R2AM (reach-twice-avoid-meet), previously studied in [21, 20] via heuristic coordination of agent-level MILPs. To evaluate the proposed BCGD-PM, we consider three planning scenarios: the R2AM baseline, and two additional test cases: 1) R2AMCA (R2AM with collision-avoidance), which incorporates global inter-agent collision avoidance specified as □ℐ⋀i≠j(‖xi(t)−xj(t)‖≥0.01) _I _i≠ j(\|x_i(t)-x_j(t)\|≥ 0.01); and 2) RURAMCA, which connects the two reach specifications (R2) in R2AMCA (or R2AM) with the until operator specified as (◆ℐcCi)[0,50](◆ℐcDi)( _I_cC_i)\,U_[0,50]( _I_cD_i). In this scenario, the robots must maintain the temporal promise of visiting their collection regions CiC_i until they reach the delivery goals DiD_i. Simultaneously, robots must fulfill their specific collaborative meeting tasks, ◆ℐmMν _I_mM_ν, while ensuring collision-free motion; both of these requirements remain unchanged from R2AMCA plan. To demonstrate the BCGD’s efficiency, we compare its performance to the LBFGS method [12, Chap. 9] implemented in a block-coordinate fashion, as an alternative solver, which uses quasi-Newton curvature information to address the inner subproblems (10) in Alg. 2. Table 1 details the computational runtime and achieved robustness of the BCGD and LBFGS penalty methods across all scenarios for linear and unicycle agent models. This comparison also underscores the modularity of our proposed penalty framework, which can wrap diverse inner solvers. Table 1: Runtime and achieved true robustness (ρ×103ρ× 10^3) for BCGD and LBFGS solvers. Solver (dynamics) R2AM R2AMCA RURAMCA BCGD (linear) 12s (0.5) 13s (0.4) 35s (1.5) LBFGS (linear) 19s (0.9) 31s (0.3) 43s (1.3) BCGD (unicycle) 234s (0.2) 288s (0.2) 480s (0.1) LBFGS (unicycle) 53s (0.7) 112s (0.4) 195s (0.2) Fig. 1 illustrates a successful RURAMCA plan generated by the proposed BCGD-PM, where Alg. 2 terminated after 6 iterations within 8 minutes, returning a solution u such that ρϕ()>10−4ρ^φ( u)>10^-4. Both BCGD and the LBFGS benchmark (implemented via jaxopt.LBFGS) utilize an identical randomized shuffling strategy, updating each agent block exactly once per 10 iterations. In BCGD, we use a simple diagonal Hessian approximation Hik=103IH_i^k=10^3I, yielding the closed-form update: ik=−(103λk+2)−1(2ik+λk∇R(k)i). d_i^k=-(10^3λ^k+2)^-1 (2 u_i^k+λ^k\,∇ R( u^k)_i ). BCGD completes each “outer” iteration in fewer than 62846284 updates on average in RURAMCA scenario. While LBFGS utilizes second-order curvature history, Table 1 shows BCGD-PM’s competitive performance, offering lower computational complexity by utilizing a simplified approximation of the Hessian curvature. The proposed algorithm is implemented in Python using the JAX library for automatic differentiation. Source code and animations for all scenarios in Table 1 are available at github.com/lefterisvl83/MAS-STL-planning. 5 Conclusions We presented a scalable framework for multi-agent planning under collaborative STL tasks. We employed a penalty method to tackle the underlying optimization problem under smooth STL semantics, using unconstrained optimization with quadratic penalty functions. We showed that the structure of the resulting unconstrained subproblems admits an efficient block-coordinate gradient descent solution, enabling agent-level computations with convergence guarantees while mitigating computational complexity. The proposed framework motivates future research into reactive STL planning for time-varying MAS structures and tasks. Appendix: Proofs of Proposition 1 and Theorem 1 Proposition 1: Let aν=ρϕν(ν)a_ν=ρ _ν( u_ν) for all ν∈ϕν _φ, and let amin=minν∈ϕaν=ρϕ()a_ = _ν _φa_ν=ρ^φ( u) denote the true robustness. By definition, the smooth robustness is ϱΓϕ()=−1Γlog∑ν∈ϕe−Γaν. _ ^φ( u)=- 1 _ν _φe^- a_ν. Writing each aνa_ν as aν=amin+(aν−amin)a_ν=a_ +(a_ν-a_ ) and factoring out e−Γamine^- a_ from the summation yields ∑ν∈ϕe−Γaν=e−Γamin∑ν∈ϕe−Γ(aν−amin). _ν _φe^- a_ν=e^- a_ _ν _φe^- (a_ν-a_ ). Substituting this expression back into the definition of ϱΓϕ() _ ^φ( u) gives ϱΓϕ()=−1Γlog(e−Γamin∑ν∈ϕe−Γ(aν−amin))=amin−1Γlog∑ν∈ϕe−Γ(aν−amin). _ ^φ( u)=- 1 \! (e^- a_ _ν _φe^- (a_ν-a_ ) )=a_ - 1 _ν _φe^- (a_ν-a_ ). Defining δΓ():=1Γlog∑ν∈ϕe−Γ(aν−amin), _ ( u):= 1 _ν _φe^- (a_ν-a_ ), we obtain ϱΓϕ()=amin−δΓ(). _ ^φ( u)=a_ - _ ( u). Since aν≥amina_ν≥ a_ for all ν∈ϕν _φ, it follows that e−Γ(aν−amin)≤1e^- (a_ν-a_ )≤ 1. Because at least on term in ∑ν∈ϕe−Γ(aν−amin) _ν _φe^- (a_ν-a_ ) is 1, each aνa_ν is finite by continuity of the system dynamics and the predicate functions, and |ϕ|≥2|K_φ|≥ 2 by assumption, the summation inside the logarithm satisfies 1<∑ν∈ϕe−Γ(aν−amin)≤|ϕ|,1< _ν _φe^- (a_ν-a_ )≤|K_φ|, which implies 0<δΓ()≤log|ϕ|Γ0< _ ( u)≤ |K_φ| for any Γ>0, >0, hence, ϱΓϕ()=amin−δΓ()<ρϕ(). _ ^φ( u)=a_ - _ ( u)<ρ^φ( u). Therefore, ϱΓϕ()≥0⇒ρϕ()>0 _ ^φ( u)≥0 ρ^φ( u)>0. Theorem 1: Let us define δΓ():=ρϕ()−ϱΓϕ()>0 _ ( u):=ρ^φ( u)- ^φ_ ( u)>0 (from Prop. 1), so for every compact set ⊂IRnU⊂ I\!R^n_ u, α≔min∈δΓ()>0 _U _ u _ ( u)>0. Note that for ⊇′U it holds that α≤α′ _U≤ _U . Let (¯⋆)N( u ) be a closed neighborhood of ¯⋆ u . We have δΓ()≥α(¯⋆) _ ( u)≥ _N( u ) for all ∈(¯⋆) u ( u ). By the assumption on the existence of a feasible limit point, there is a subsequence (¯⋆,ιs)s( u , _s)_s which converges to a KKT point, ⋆ u , of (10). So there is s∈INs∈ I\!N so that R(¯⋆,ιs)<α(¯⋆)2R( u , _s)< _N( u )^2. Therefore, ϱΓϕ(¯⋆,ιs)>−α(¯⋆) ^φ_ ( u , _s)>- _N( u ), hence ρϕ(¯⋆,ιs)>0ρ^φ( u , _s)>0. References [1] C. Belta, B. Yordanov, and E. A. Gol (2017) Formal methods for discrete-time dynamical systems. Studies in Systems, Decision and Control, Vol. 89, Springer. Cited by: §1. [2] M. Charitidou and D. V. Dimarogonas (2024) Distributed MPC with continuous-time STL constraint satisfaction guarantees. IEEE Control Systems Letters 8 (), p. 211–216. External Links: Document Cited by: §1. [3] A. Donzé and O. Maler (2010) Robust Satisfaction of Temporal Logic over Real-Valued Signals. In 8th International Conference on Formal Modeling and Analysis of Timed Systems, FORMATS, 2010, Klosterneuburg, Austria, p. 92–106. Cited by: §1, §2.1. [4] G. E. Fainekos and G. J. Pappas (2009) Robustness of temporal logic specifications for continuous-time signals. Theoretical Computer Science 410 (42), p. 4262–4291. Cited by: §1. [5] Y. Gilpin, V. Kurtz, and H. Lin (2020) A smooth robustness measure of signal temporal logic for symbolic control. IEEE Control Systems Letters 5 (1), p. 241–246. Cited by: §1, §1, §3.1. [6] S. Han, J. Verhagen, and J. Tumova (2025) Exact smooth reformulations for trajectory optimization under signal temporal logic specifications. External Links: 2511.07375, Link Cited by: §1. [7] M. Kazemi and S. Soudjani (2020) Formal policy synthesis for continuous-state systems via reinforcement learning. In Integrated Formal Methods: 16th International Conference, IFM 2020, Lugano, Switzerland, November 16–20, 2020, Proceedings 16, p. 3–21. Cited by: §3.1. [8] V. Kurtz and H. Lin (2022) Mixed-Integer Programming for Signal Temporal Logic With Fewer Binary Variables. IEEE Control Systems Letters 6, p. 2635–2640. Cited by: §1. [9] L. Lindemann and D. V. Dimarogonas (2019) Feedback control strategies for multi-agent systems under a fragment of signal temporal logic tasks. Automatica 106, p. 284–293. Cited by: §1. [10] O. Maler and D. Nickovic (2004) Monitoring Temporal Properties of Continuous Signals. In Formal Techniques, Modelling and Analysis of Timed and Fault-Tolerant Systems, Y. Lakhnech and S. Yovine (Eds.), Berlin, Heidelberg, p. 152–166. Cited by: §1, §2.1, §2.1, §2.1. [11] N. Mehdipour, C. Vasile, and C. Belta (2019) Arithmetic-geometric mean robustness for control from signal temporal logic specifications. In 2019 American Control Conference (ACC), p. 1690–1695. Cited by: §1. [12] J. Nocedal and S.J. Wright (2006) Numerical optimization. Springer. Cited by: §1, §1, §3.2, §3.4, §3.4, §4. [13] J. Orlin (1977) Contentment in graph theory: covering graphs with cliques. Indagationes Mathematicae 80 (5), p. 406–424. Cited by: §2.2.1. [14] Y. V. Pant, H. Abbas, and R. Mangharam (2017) Smooth operator: control using the smooth robustness of temporal logic. In 2017 IEEE Conference on Control Technol. and Applic. (CCTA), Vol. , p. 1235–1240. Cited by: §1. [15] V. Raman, M. Maasoumy, A. Donze, R. M. Murray, A. Sangiovanni-Vincentelli, and S. A. Seshia (2014) Model predictive control with signal temporal logic specifications. In Proceedings of the IEEE Conf. on Decis. and Cont., p. 81–87. Cited by: §1. [16] S. Sadraddini and C. Belta (2015) Robust temporal logic model predictive control. In 53rd Annual Allerton Conference on Communication, Control, and Computing, p. 772–779. Cited by: §2.1. [17] P. Sopasakis, E. Fresk, and P. Patrinos (2020) OpEn: code generation for embedded nonconvex optimization. In IFAC WC, Berlin. Cited by: §3.3, §3.4. [18] P. Tseng and S. Yun (2009) A coordinate gradient descent method for nonsmooth separable minimization. Mathematical Programming 117 (1), p. 387–423. Cited by: §1, §3.3, §3.3, §3.3. [19] P. Tseng (2001) Convergence of a block coordinate descent method for nondifferentiable minimization. Journal of Optimization Theory and Applications 109 (3), p. 475–94. Cited by: §3.3. [20] E. E. Vlahakis, L. Lindemann, and D. V. Dimarogonas (2025) Conformal data-driven control of stochastic multi-agent systems under collaborative signal temporal logic specifications. In 2025 IEEE 64th Conference on Decision and Control (CDC), Vol. , p. 624–629. Cited by: §4. [21] E. E. Vlahakis, L. Lindemann, P. Sopasakis, and D. V. Dimarogonas (2024) Probabilistic tube-based control synthesis of stochastic multi-agent systems under signal temporal logic. In 2024 IEEE 63rd Conference on Decision and Control (CDC), Vol. , p. 1586–1592. Cited by: §1, §4. [22] X. Zhou, Y. Zou, S. Li, X. Li, and H. Fang (2022) Distributed model predictive control for multi-robot systems with conflicting signal temporal logic tasks. IET Control Theory & Appl. 16 (5), p. 554–572. Cited by: §1.