# Patent

## Title
Method and device for three-weight message-passing optimization scheme

## Abstract

A method and device determines an optimization solution for an optimization problem. The method includes receiving the optimization problem having cost functions and variables in which each of the cost functions has a predetermined relationship with select ones of the variables. The method includes generating a first message for each of the cost functions for each corresponding variable based upon the respective predetermined relationship and a second message for each of the variables for each corresponding cost function based upon the respective predetermined relationship. The method includes generating a disagreement variable for each corresponding pair of variables and cost functions measuring a disagreement value between the first and second beliefs. The method includes repeating steps (b), (c), and (d) until a consensus is formed between the first and second messages until the optimization solution is determined based upon the consensus.

## Background
RELATED APPLICATION
The present application claims priority to U.S. Provisional Appln. Ser. No. 61/816,305 entitled “Three-Weight Message-Passing Algorithm based on Alternating Direction Method of Multipliers (ADMM)” filed on Apr. 26, 2013, the entirety of which is incorporated by reference herein.


BACKGROUND INFORMATION
An optimization algorithm may be used to determine a solution for a problem in which the solution considers all variables and constraints related to the problem and provides a lowest cost configuration of the values for all the variables. For example, in a televised event, a plurality of cameras may be used to capture the event. The cameras may move along planned trajectories to fully capture the event. Accordingly, the problem associated with such a scenario may be to determine the planned trajectories for the cameras that incorporate the variables and constraints associated therewith. Such variables and constraints may include no overlap of space as the cameras cannot physically be co-located at a common time, a trajectory path that provides a sufficient video capture of the event, a tracking algorithm for the occurrences during the event, etc. Through incorporation of all these considerations, the optimization algorithm may determine the trajectories of the cameras to provide the lowest cost solution (e.g., least distance to be covered by the cameras, best coverage of the event, least energy requirements, etc.). However, most optimization algorithms require significantly high processing requirements for all these considerations to be assessed as well as provide the lowest cost solution to the problem.
A conventional algorithm used for convex optimization is the Alternating Direction Method of Multipliers (ADM). Conventionally, the ADMM algorithm is a variant of an augmented Lagrangian scheme that uses partial updates for dual variables. This algorithm has also been utilized and shown to be well-suited for distributed implementations. Further developments in using the ADMM algorithm have included being used in applications such as Linear Programming (LP) decoding of error-correcting codes and compressed sensing.
However, the ADMM algorithm has in practice been used exclusively for convex optimization problems. Convex optimization problems have a special form for the cost objective function, and in particular, only one local optimum exists for a given problem (e.g., a landscape including one valley has a single solution as to what is the lowest point). Many optimization problems of interest, such as the problem of finding optimal trajectories introduced earlier, are not convex. The basic principles behind the ADMM algorithm including the convexity assumptions underlying the derivation of ADMM, appear to provide a gross contradiction to its application to non-convex problems. That is, conventional uses of the ADMM algorithm are not designed for use with general optimization functions. In particular, the ADMM algorithm was not designed as a heuristic optimization algorithm for non-convex Non-Deterministic Polynomial-Time (NP) hard optimization problems. In such problems, multiple local optima exist for a given problem (e.g., a landscape including various hills and valleys has more than one solution as to what is the lowest point as well as the highest point).
Another conventional approach for optimization uses graphical models in which a graph denotes conditional dependencies between random variables. The graphical models often are combined with a Belief Propagation (BP) message-passing algorithm in which inferences are performed on the graphical models. The BP message-passing algorithm calculates the marginal distribution for each unobserved node, conditional on any observed nodes within the graphical model.



BRIEF DESCRIPTION OF THE DRAWINGS
FIG. 1 shows a device for generating an optimization solution according to an exemplary embodiment.
FIG. 2 shows a graphical model according to an exemplary embodiment.
FIG. 3A shows an unsolved puzzle.
FIG. 3B shows a solved puzzle using the device of FIG. 1 according to an exemplary embodiment.
FIG. 4 shows a packing problem and solution using the device of FIG. 1 according to an exemplary embodiment.
FIG. 5A shows a graphical model as applied to a multi-agent trajectory planning problem according to an exemplary embodiment.
FIG. 5B shows input and output variables for each minimizer block of the graphical model of FIG. 5A according to an exemplary embodiment.
FIG. 6 shows a method of generating an optimization solution according to an exemplary embodiment.



DETAILED DESCRIPTION
The present invention relates to a device and method for determining an optimization solution for an optimization problem. The method comprises (a) receiving the optimization problem, the optimization problem including a plurality of cost functions and a plurality of variables, the cost functions representing possible costs for values of the variables in the optimization problem, each of the cost functions having a predetermined relationship with select ones of the variables; (b) generating a first message for each of the cost functions for each corresponding variable based upon the respective predetermined relationship, the first message indicating a first belief that the corresponding variable has a first value when the optimization problem is solved, the first message having a respective first weight indicating a certainty of the first message; (c) generating a second message for each of the variables for each corresponding cost function based upon the respective predetermined relationship, the second message indicating a second belief that the corresponding variable has a second value when the optimization problem is solved the second message having a respective second weight indicating a certainty of the second message; (d) generating a disagreement variable for each corresponding pair of variables and cost functions measuring a disagreement value between the first and second beliefs; (e) repeating steps (b), (c), and (d) until a consensus is formed between the first and second messages, a subsequent first message being modified based upon the second message, its corresponding second weight, and the corresponding disagreement variable, and a subsequent second message being modified based upon the first message, its corresponding first weight, and the corresponding disagreement variable; and (f) determining an optimization solution based upon the consensus.
The exemplary embodiments may be further understood with reference to the following description and the appended drawings, wherein like elements are referred to with the same reference numerals. The exemplary embodiments describe a method and device for generating an optimization solution. Specifically, the method and device utilize a message-passing version of the ADMM algorithm. The message-passing algorithm includes a weighting system in which a message passed between nodes of a graphical model has an infinite weight, a standard weight, or an indifferent weight. The optimization solution, the generation thereof, the message-passing algorithm, the ADMM algorithm, the confidence weighting system, and a related method will be described in further detail below.
The exemplary embodiments provide how a Divide and Concur (DC) algorithm for constraint satisfaction may be derived as a specific case of a message-passing version of the ADMM algorithm for convex optimization. The message-passing algorithm is improved upon through introduction of three distinct weights for the messages between nodes of a graphical model with “certain” and “no opinion” weights as well as the standard weight used in both ADMM and DC. The “certain” messages allow the optimization algorithm to implement constraint propagation as a special case while the “no opinion” messages speed convergence for some problems by making the optimization algorithm focus only on active constraints. The three-weight version of message-passing in the ADMM/DC further provides greatly improved performance in non-convex problems while retaining the exact performance of ADMM for convex problems.
The exemplary embodiments overcome the conventional limits of the ADMM algorithm and are used as a well-defined algorithm for general optimization functions and provide a powerful heuristic algorithm even for NP-hard non-convex problems. As will be described in further detail below, when the objective function being optimized by the ADMM algorithm consists entirely of “hard” constraints, the ADMM algorithm reduces to a message-passing version of the DC constraint satisfaction algorithm which is effective for a very wide variety of non-convex constraint satisfaction problems.
As will be described in further detail below, the use of the message-passing algorithm between variables and cost functions in the graphical model enables messages to be sent between these nodes until a consensus is reached; that is, until an optimization solution has been generated. Specifically, messages are transmitted from right side nodes representing a first set of functions to left side nodes representing a second set of functions. The left side nodes process those messages and compute new messages that are transmitted back from left side nodes to the right side nodes. The right side nodes then use the messages that have come from the left side nodes to compute new messages that are sent to the left side nodes, and the cycle continues. The weights associated with the messages are determined such that the messages have an appropriate impact in the processing of the messages to compute new messages. The weights for the messages have one of three different weights. A “certain” weight indicates that a cost function or variable (i.e., one of the nodes) is certain (e.g., 100% confidence) about a value of a variable and all other opinions should defer thereto; a “no opinion” weight indicates that the node has no opinion about the value (e.g., may be omitted in the considerations for the consensus); and a “standard” weight indicates that the message is weighted equally with all other standard weight messages.
FIG. 1 shows a device 100 for generating an optimization solution according to an exemplary embodiment. The device 100 may be any electronic component that is configured to receive and process data such that the generating is performed. The device 100 may include a processor 105, a memory arrangement 110, a transceiver 115, and an input/output (I/O) device 120. The processor 105 may be configured to execute an optimization application. The memory 110 may store data that is received as well as data that is determined. The transceiver 115 and the I/O device 120 may provide a communicative connection to other electronic devices or users such that data may be received for processing. For example, a user may manually enter constraints involved in a particular optimization scenario or problem.
The processor 105 may further execute a plurality of engines that are configured to perform a functionality associated with generating the optimization solution. As illustrated, the processor 105 may execute a variable engine 125, a cost engine 130, a factor graph engine 135, a constraint engine 140, an optimization engine 145, and a signal generating engine 150. The engines will be described in further detail below.
The optimization functionality executed by the processor 105 utilizes a graphical model in order to generate the optimization solution. When relating to optimization, the graphical model includes a first plurality of cost functions and a second plurality of variables, as well as edges connecting cost functions to variables. The cost functions may represent either “soft” costs or “hard” constraints. An edge exists between a cost function and a variable if the variable is involved in the cost function. In a “soft” cost function, the cost function encodes a preference for the connected variables to have particular values, but that preference may be violated, while in a “hard” constraint version of a cost function, the preferences must not be violated. Each cost function may be connected to one or more of the variables and vice versa.
In view of the manner in which the optimization functionality is utilized, the variable engine 125 may receive the set of variables while the cost engine 130 may receive the set of cost functions. The factor graph engine 135 may generate the graphical model for the variables and the cost functions. Through the connections between the cost functions and the variables, the constraint engine 140 may provide the lines of the graphical model to the factor graph engine 135. All this data may be received by the optimization engine 145 that generates the optimization solution. The signal generating engine 150 may generate a signal based upon the generated optimization solution such that a particular application that requires this solution may receive it. For example, a particular application may be trajectories for a plurality of cameras. The signal generating engine 150 may generate signals for each camera that indicate the trajectory at which the camera is to move.
Initially, a message-passing version for the ADMM algorithm is derived. The message-passing algorithm used for the ADMM algorithm as described below may be applied to a completely general optimization problem and not only to convex optimization problems. The general optimization problem may be finding a configuration the values of some variables that minimize some objective function subject to some constraints. The variables may be continuous or discrete. For purposes of this description, the variables are represented as continuous. However, the variables may also be represented as discrete through addition of constraints that enforce the discrete nature of the variables. For example, with a variable that takes on Q discrete labels, Q binary indicator variables may be used to represent whether the original variable is in one of the Q states and each binary variable is represented as a continuous variable but subject to the constraint that exactly one of the Q indicator variables is equal to 1 while the rest are equal to 0.
The processor 105 may receive N continuous variables via the variable engine 125 and represent these variables as a vector r∈RN. The general optimization problem becomes one of minimizing an objective function E(r) subject to some constraints on r. The objective function may be, for example, a cost function. All the constraints are considered to be part of the objective function by introducing a cost function Ea(r) for each constraint such that Ea(r)=0 if the constraint is satisfied and Ea(r)=∞ if the constraint is not satisfied. The original objective function E(r) may also be decomposable into a collection of local cost functions. For example, for the problem of minimizing an Ising Model objective function, there is a collection of local cost functions representing the “soft” pairwise interactions and local fields. There are also other cost functions representing the “hard” constraints that each Ising variable may only take on two possible values.
The general problem of minimizing an objective function may be written as Σa=1MEa(r) where there are M cost functions that are either “soft” or “hard” and r∈Rn. Such a problem may have a standard “factor graph” representation. FIG. 2 shows a graphical model 200 according to an exemplary embodiment. Specifically, the graphical model 200 may be the standard “factor graph” for the above described optimization problem.
The graphical model 200 includes a set of hard or soft cost functions E1, E2, E3, and E4 as well as a set of variables r1, r2, r3, r4, r5, and r6. When a line connects a cost function node with a variable node, this indicates that the cost function depends upon the corresponding variable. In order to derive the message-passing algorithm for the ADMM algorithm, the problem may be manipulated into a series of equivalent forms before actually minimizing the objective. The first manipulation is a conversion of the problem over the variables rj into an equivalent problem that depends on variables xij that sit on the edges of a “normalized” Forney-style factor graph. The variables in the standard factor graph are replaced with equality constraints and each of the edge variables is a copy of the corresponding variable that was on its right. Thus, the edge variables attached to the same equality constraint must ultimately equal each other but they may temporarily be unequal while they separately attempt to satisfy different cost functions on the left. The below description uses the notation x to represent a vector consisting of the entire collection of xij edge variables while the notation z represents a vector consisting of the entire collection of zij edge variables. It is noted that x normally has a higher dimensionality than r.
As illustrated in FIG. 2, the cost function E1 may be connected to the variables r1 and r2. Thus, the cost function E1 has variables x11 and x12. The cost function E2 may be connected to variables r1 and r3. Thus, the cost function E2 has variables x21 and x23. The cost function E3 may be connected to variables r2, r3, r4, r5, and r6. Thus, the cost function E3 has variables x32, x33, x34, x35, and x36. The cost function E4 may be connected to variables r5 and r6. Thus, the cost function E4 has variables x45 and x46. In this manner, the variable r1 may have variables z11 and z21; the variable r2 may have variables z12 and z32; the variable r3 may have variables z23 and z33; the variable r4 may have variable z34; the variable r5 may have variables z35 and z45; and the variable r6 may have variables z36 and z46. It is noted that the graphical model 200 is only one exemplary embodiment of a graphical model that may represent the interconnections between the cost functions and the variables as well as the ensuing copied variables. It is also noted that the number of nodes on each side of the graphical model 200 is only exemplary as more or fewer nodes may exist on each side.
Due to the bipartite structure of a Formey factor graph such as the graphical model 200, the cost functions may be split into two groups: the nodes on the left may represent the original soft cost functions and hard constraints and the nodes on the right may represent equality constraints. Each edge variable xij sits on a left side of the edge and a copy may be made denoted as zij that sits on the right side of the edge. The xij and zij variables are referred to as “belief” variables because the xij variables may intuitively be thought of as “beliefs” that the left-hand side of the factor graph has about the values of the original variables rj, while the zij variables may be thought of as “beliefs” that the right-hand side of the factor graph has about the values of the original variables rj. The objective function E(x) may then be split into a sum of the left cost functions denoted as f(x) and the right cost functions denoted as g(z), where z is a vector made from all the zij belief variables.
The constraint that each edge belief variable xij equals its copy belief variable zij is enforced by a Lagrange multiplier yij in a Langrangian represented as L(x, y, z)=f(x)+g(z)+y·(x−z). As will be shown below, it is useful to add another term to augment the Lagrangian. Specifically, the term (ρ/2)(x−z)2 may be added. Since x=z at the optimum, this added term is zero at the optimum and non-negative elsewhere. Therefore, the optimum is not changed from addition of this term. The parameter ρ may be considered a scalar value. However, this parameter may also be generalized to be a vector with a different ρij for each edge. In summary, the original problem of minimizing E(r) has become equivalent to finding the minimum of the augmented Lagrangian with the following:




  
    
      
        
          L
          ⁡
          
            (
            
              x
              ,
              y
              ,
              z
            
            )
          
        
        =
        
          
            f
            ⁡
            
              (
              x
              )
            
          
          +
          
            g
            ⁡
            
              (
              z
              )
            
          
          +
          
            y
            ·
            
              (
              
                x
                -
                z
              
              )
            
          
          +
          
            
              ρ
              2
            
            ⁢
            
              
                (
                
                  x
                  -
                  z
                
                )
              
              2
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          1
        
        )
      
    
  




Accordingly, solving constrained optimization problems is equivalent to a minimization problem that naturally splits into three pieces: (1) minimizing the original soft and hard cost functions f(x) on the left side of the graphical model 200, (2) minimizing the equality cost functions g(z) on the right side of the graphical model 200, and (3) ensuring that the x=z with the Lagranage multipliers y.

To further define the message-passing algorithm for the ADMM algorithm, it may be assumed that each of the local cost functions ƒa(x) on the left side of the graphical model 200 is convex. However, it is again noted that the final algorithm is well-defined even when this assumption is violated (i.e., non-convex). All the equality cost functions on the right side of the graphical model 200 are clearly convex and the augmenting of the quadratic terms are convex. Accordingly, the overall function is convex as well. Therefore, the minimum of the Lagrangian may be determined through maximizing the dual function of the following:

h(y)=L(x*y,z*)  (Equation 2)

where (x*,z*) are the values of x and z that minimize L for a particular choice of y:





  
    
      
        
          (
          
            
              x
              *
            
            ,
            
              z
              *
            
          
          )
        
        =
        
          
            
              arg
              ⁢
              
                
              
              ⁢
              min
            
            
              x
              ,
              z
            
          
          ⁢
          
            
          
          ⁢
          
            L
            ⁡
            
              (
              
                x
                ,
                y
                ,
                z
              
              )
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          3
        
        )
      
    
  




In order to maximize h(y), a gradient ascent algorithm may be used. Thus, given values of yt at some iteration t, an iterative computation may be performed to obtain




  
    (
    
      
        x
        *
      
      ,
      
        z
        *
      
    
    )
  
  =
  
    
      
        arg
        ⁢
        
          
        
        ⁢
        min
      
      
        x
        ,
        z
      
    
    ⁢
    
      
    
    ⁢
    
      L
      ⁡
      
        (
        
          x
          ,
          
            y
            t
          
          ,
          z
        
        )
      
    
  




and then move in the direction of the gradient h(y) according to the following:





  
    
      
        
          y
          
            t
            +
            1
          
        
        =
        
          
            
              y
              t
            
            +
            
              α
              ⁢
              
                
                  ∂
                  h
                
                
                  ∂
                  y
                
              
            
          
          =
          
            α
            ⁡
            
              (
              
                
                  x
                  *
                
                ⁢
                
                  z
                  
                    *
                    
                      
                    
                  
                
              
              )
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          4
        
        )
      
    
  




where α is a step-size parameter.

The bipartite structure of the graphical model may be used to decouple the minimizations of the belief variables x and z. Introducing a scaled Lagrange multiplier and the “messages” mt=xt+ut and nt=zt−ut, the following iterative equations which define the message-passing version of ADMM may be obtained.




  
    
      
        
          x
          t
        
        =
        
          
            
              arg
              ⁢
              
                
              
              ⁢
              min
            
            x
          
          ⁡
          
            [
            
              
                f
                ⁡
                
                  (
                  x
                  )
                
              
              +
              
                
                  (
                  
                    ρ
                    /
                    2
                  
                  )
                
                ⁢
                
                  
                    (
                    
                      x
                      -
                      
                        n
                        t
                      
                    
                    )
                  
                  2
                
              
            
            ]
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          5
        
        )
      
    
  
  
    
      
        
          z
          
            t
            +
            1
          
        
        =
        
          
            
              arg
              ⁢
              
                
              
              ⁢
              min
            
            z
          
          ⁡
          
            [
            
              
                g
                ⁡
                
                  (
                  z
                  )
                
              
              +
              
                
                  (
                  
                    ρ
                    /
                    2
                  
                  )
                
                ⁢
                
                  
                    (
                    
                      z
                      -
                      
                        m
                        t
                      
                    
                    )
                  
                  2
                
              
            
            ]
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          6
        
        )
      
    
  
  
    
      
        
          u
          
            t
            +
            1
          
        
        =
        
          
            u
            t
          
          +
          
            
              α
              ρ
            
            ⁢
            
              (
              
                
                  x
                  t
                
                -
                
                  z
                  
                    t
                    +
                    1
                  
                
              
              )
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          7
        
        )
      
    
  




The algorithm is initialized by selecting u0=0 and starting with some initial z0. A repeating process may be performed in which (1) Equation 5 determines x0, Equation 6 determines z1, and Equation 7 determines u1, (2) Equation 5 determines x1, Equation 6 determines z2, and Equation 7 determines u2, (3) etc. Intuitively, the z and x variables are analogous to the single-node and multi-node “beliefs” in BP while m and n are messages from the nodes on the left side of the graphical model 200 to those on the right side of the graphical model 200 and vice-versa, respectively, much as in the “two-way” version of the standard BP algorithm. In each iteration, the algorithm computes beliefs on the left side of the graphical model 200 based upon the messages coming from the right side of the graphical model 200, then beliefs on the right side of the graphical model 200 based upon the messages from the left side of the graphical model 200. It may intuitively be said that the messages m indicate the beliefs that the left hand side of the graph has for the values of the r variables that minimize the cost functions, while the messages n indicate the beliefs that the right hand side of the graph has for the values of the r variables that minimize the cost functions. Then the algorithm attempts to equalize the beliefs on the left and right sides of the graphical model 200 using the u variables. It is noted that the u variables keep a running total of the differences between the beliefs on the left and right sides of the graphical model 200 and are similar to control variables tracking a time-integrated difference from a target value. Because they track and try to correct the disagreement between the beliefs of the left and right hand side of the graph, the u variables may be referred to as “disagreement variables” below.
It is also noted that all the updates in the above equations are local computations and may be performed in parallel. Thus, if a function cost a is connected to a small set of edges with belief variables {x}a, then it only needs to look at the messages {n}at on the same edges to perform its local computation using the following:




  
    
      
        
          
            {
            x
            }
          
          a
          
            t
            +
            1
          
        
        =
        
          
            
              arg
              ⁢
              
                
              
              ⁢
              min
            
            
              
                {
                x
                }
              
              a
            
          
          ⁡
          
            [
            
              
                
                  f
                  a
                
                ⁡
                
                  (
                  
                    
                      {
                      x
                      }
                    
                    a
                  
                  )
                
              
              +
              
                
                  (
                  
                    ρ
                    /
                    2
                  
                  )
                
                ⁢
                
                  
                    (
                    
                      
                        
                          {
                          x
                          }
                        
                        a
                      
                      -
                      
                        
                          {
                          n
                          }
                        
                        a
                        t
                      
                    
                    )
                  
                  2
                
              
            
            ]
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          8
        
        )
      
    
  




These local computations are usually easy to implement with small “minimizing” subroutines (also called “minimizers” below) specialized to the particular ƒa({x}a). Such a subroutine balances the desire to minimize the local ƒa({x}a) with the desire to agree with the {n}at messages coming from other nodes. The ρ parameter enables the relative strength to be varied for these competing influences.
The minimizations in Equation 6 are similarly all local computations that may be performed in parallel. Because the ga({z}a) functions on the right side of the graphical model 200 all represent hard equality constraints, these minimizations reduce to a particularly simple form (the output belief variables z are given as the mean of the incoming m messages), as will be described in further detail below.
As discussed above, when all the ƒa({x}a) functions are convex, the overall problem is convex and the ADMM algorithm provably converges to the correct global minimum. However, it is noted that no guarantees are made about the speed of convergence. Starting with an arbitrary optimization problem and using Formey factor graphs guarantees that ga({z}a) functions on the right side of the graphical model 200 are all equality constraints that are convex. More generally, the ADMM algorithm may be considered to be an algorithm that operates on any functions f(x) and g(z) to provide exact answers so long as f(x) and g(z) are both convex. However, the algorithm is perfectly well defined even for problems where the ƒa({x}a) functions are not convex, as will be described in further detail below.
The message-passing algorithm used with the ADMM algorithm may also be a special case that is a message-passing version of the DC algorithm. Specifically, the special case is for problems where the ƒa({x}a) functions all represent hard constraints. In this case, the legal values of the belief variables {x}a may be written as a constraint set D and require that ƒa({x}a)=0 for {x}a∈D and ƒa({x}a)=∞ for {x}a∉D. Then each local minimization to compute {x}at+1 as given by Equation 8 reduces to a projection of the incoming messages {n}a onto the constraint set with the following form:

{x}t+1=PD({n}at)  (Equation 9)

The reason that Equation 8 reduces to a projection in this case is that the ƒa({x}a) term in Equation 8 enforces that {x}a must be in the constraint set while minimizing the (ρ/2)({x}a−{n}at)2 term enforces that the computed {x}a values are as close as possible (using a Euclidean metric) to the messages {n}at which is the definition of a projection.

Similarly, the cost functions ga({z}) represent hard equality constraints such that the {z}at+1 updates may be written as projections of the {m}at messages onto the equality constraint sets C with the following:

{z}=at+1=PC({m}at)  (Equation 10)

Assuming that all the ρ weights are equal, this may be further simplified. Specifically, the {z}at+1 values on the edges connected to an equality constraint node are all equal to the mean of the messages incoming to that equality node. However, as will be described in further detail below, the algorithms may go beyond this assumption.

If the step-size parameter α is selected to equal the weight ρ, a further simplification may be performed to show that instead of requiring updates of all the variables xt, zt, ut, mt, and nt, the algorithm reduces to an iteration of a single state variable: the messages mt. That is, straightforward algebraic manipulations of Equations 5, 6, and 7 along with the definitions of mt and nt permit the elimination of the xt, zt, and ut variables and reduce to the following message-update equations:

nt+1=2PC(mt)−mt  ((Equation 11)

mt+1=mt+PD(nt+1)−PC(mt)  (Equation 12)

This equivalently leaves a single update equation for the mt messages which utilizes a substantially similar process of the “difference-map” iteration of the DC algorithm as the following:

mt+1=PD(2PC(mt)−mt)−(PC(mt)−mt)  (Equation 13)

It is noted that certain cases exist where the DC algorithm uses the dual version of the difference-map where nt is updated according to a rule obtained from Equation 13 by swapping the PC and PD projections.

It is again noted that the DC algorithm described above is a special case when used with the ADMM algorithm. However, during the derivation of the above equations, the weights ρ in the ADMM algorithm have disappeared in the DC update equations. The weights for the algorithms described herein have an intuitive meaning in the message-passing algorithm in the ADMM algorithm. Specifically, this intuitive meaning reflects how strongly the messages are to be adhered to in comparison to local function costs. That is, the messages have a “reliability” or “certainty” related thereto. Accordingly, a generalized version of the message-passing algorithm for the ADMM algorithm is determined where each edge (ij) connecting a function cost i to an equality node j is given its own weight of ρij reflecting the certainty of messages on that edge. Again, as described above, the message-passing algorithm utilizes three weights ρ: (1) a “certain” weight, (2) a standard weight, and (3) a “no opinion” or indifferent weight.
The message-passing algorithm may be further generalized with different weights for messages transmitted from the nodes on the left side of the graphical model 200 to the nodes on the right side of the graphical model 200 in which the weights may be changed with each iteration (i.e., each transmission of the messages). The vector of weights going to the left nodes of the graphical model 200 at a time t may be denoted as t while the vector of weights going to the right nodes of the graphical model 200 at a time t may be denoted as {right arrow over (ρ)}t. When a particular weight on an edge (ij) is to be denoted, the weights may be modified as ijt or {right arrow over (ρ)}ijt. It is noted that for a convex problem, the leftward and rightward weights may eventually be equal to each other and may be constant to converge to a global optimum.
Therefore, the modification of the message-passing algorithm is performed with a relatively simple modification that allows for only three possible choices for the weights on each edge. A first weight is the standard weight messages with some weight ρ0 that is greater than zero and less than infinity. The exact choice of ρ0 is important for the rate of convergence for problems with soft cost functions but is irrelevant for problems consisting entirely of hard constraints as is the case for standard DC algorithms. For simplicity, the standard weight ρ0 may be chosen to equal ρ0=1 for such problems. A second weight is for infinite-weight messages which represent the message's value is certainly correct. A third weight is for zero-weight messages which represent that a function cost node or equality node is completely uncertain about the value that a variable is to have and its opinion may be ignored. Therefore, the message-passing algorithm introduces the infinite weight and the zero weight while the standard weight is modified in view of the introduction of the infinite and zero weights.
To modify the message-passing algorithm used with the ADMM algorithm to allow for zero weights and infinite weights, the updating of the uij variables must be properly handled. The uij disagreement variables track the “disagreement” between the beliefs of the left and right nodes of the graphical model 200 on an edge (ij). The uij disagreement variable on an edge grows in magnitude over time if the left belief xij is persistently less than or persistently greater than the right belief zij. Because the uij disagreement variables are added or subtracted to the beliefs to form the messages, the message values become different from the actual belief values as the uij disagreement variables attempt to resolve the disagreement. With infinite and zero weight messages, the algorithms utilize a “reset” for the uij disagreement variables to zero if there is no disagreement on that edge. For example, when an infinite weight message is sent on an edge, it may indicate that the message is certainly correct such that any previous disagreement recorded in the uij variable is to be ignored.
In order to utilize these weights in the message-passing algorithm in the ADMM algorithm, all the uij disagreement variables may be set to equal 0 initially. Thus, uij=0. Accordingly, zero weights (e.g., ij9=0 on each edge) may normally be used for the initial n0 messages from the right node to the left node for the variables in which no information is received. For any variables for which there is a certainty, these messages may be accompanied with infinite weights. The x0 belief may be computed using the standard update equations with the following:




  
    
      
        
          x
          t
        
        =
        
          
            
              arg
              ⁢
              
                
              
              ⁢
              min
            
            x
          
          ⁡
          
            [
            
              
                f
                ⁡
                
                  (
                  x
                  )
                
              
              +
              
                
                  (
                  
                    
                      
                        ρ
                        ←
                      
                      t
                    
                    /
                    2
                  
                  )
                
                ⁢
                
                  
                    (
                    
                      x
                      -
                      
                        n
                        t
                      
                    
                    )
                  
                  2
                
              
            
            ]
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          14
        
        )
      
    
  




Any ties in the xt updates may be broken randomly. The left-to-right messages mt are computed using mt=xt+ut but since u0=0, the initial messages to the right equal the initial x0 beliefs.

The outgoing {right arrow over (ρ)}t weights are computed using an appropriate logic for the function cost on the left that depend on an analysis of the function. For example, as will be described in further detail below, with a mathematical puzzle such as a Sudoku puzzle, only standard or infinite weights may be transmitted depending on a logic that sends out infinite weights when there is a certainty about the corresponding x value. Whenever an infinite weight is used, whether for a right-going message or a left-going message, the uijt variable for that edge is immediately reset to zero.
The zt+1 right beliefs may be computed by taking a weighted average of the mt messages weighted by the {right arrow over (ρ)}t weights. Thus, if any message has an infinite weight, it controls the average while any zero-weight message contributes nothing to the average. If the logic used to send infinite weights is correct, there are no contradictions between infinite weight messages.
To compute the weights returning to the left from an equality node (right nodes), the following logic may be utilized. First, if any edge is sending in an infinite weight {right arrow over (ρ)}ijt, all edges out of the equality node return an infinite weight ijt+1. Otherwise, all edges return a standard weight ijt+1 as long as at least one of the incoming weights is non-zero. Finally, if all incoming weights {right arrow over (ρ)}ijt are zero, the outgoing weights ijt+1 are also set to zero.
The logic continues with all u disagreement variables being updated. Any uij disagreement variable on an edge that has an infinite weight in either direction is reset to zero. Also, any edge that has a zero weight {right arrow over (ρ)}ijt has its u variable reset to zero (the reasoning is that there is no contribution to the average and should agree with the consensus of the rest of the system). Any edge that has a standard weight {right arrow over (ρ)}ijt while all other edges into its equality node have zero weight has its u disagreement variable reset to zero (the reasoning again is that there was no disagreement so there is no longer any need to modify the right belief). Any other edge that has a standard weight {right arrow over (ρ)}ijt and a standard weight {right arrow over (ρ)}ijt+1 will have its uijt+1 disagreement variable updated according to the formula uijt+1=uijt+(α/ρ0)(xjt−zijt+1). Once the u variables are updated, all the right-to-left messages nt+1 may be updated according to the formula nijt+1=zijt+1−uijt+1.
The above description describes one iteration of the optimization algorithm. This process may repeat until a stopping criterion is satisfied in which all the m and n messages are identical from iteration to iteration to some specified numerical tolerance.
As will be described in further detail below, the exemplary three-weight messaging algorithm exhibits several properties that make it ideal for cognitive systems. In a first example, the exemplary three-weight messaging algorithm operates on arbitrary objective functions (e.g., non-linear, non-convex), constraints (e.g., hard, soft, mixed), and variables (e.g., discrete, continuous). The three-weight message-passing algorithm converges to a global minimum for convex problems and, if it converges, the three-weight message-passing algorithm arrives at a feasible solution (e.g., all hard constraints are met). In a second example, the exemplary three-weight message-passing algorithm is an iterative algorithm and, for many problems, intermediate results serve as heuristic input for warm-starting complementary approaches. In a third example, the exemplary three-weight message-passing algorithm takes the form of a decomposition-coordination procedure in which the solutions to small local subproblems are coordinated to find a solution to a large global problem. The algorithmic structure leads naturally to concurrent processing at multiple levels of execution (e.g., MapReduce, multi-core, GPU).
As discussed above, the exemplary optimization algorithm may be interpreted as an iteration loop that consists of two phases: (1) minimizing each cost function locally, and (2) concurring on the result of the local computation. The messages of the three-weight message-passing algorithm operate on the edges of the corresponding graphical model as opposed to directly on variables which raises the dimensionality of the problem, thereby allowing the exemplary three-weight message-passing algorithm to more effectively search the variable space and avoid invalid solutions.
To summarize, the minimization phase takes an input pair of message and weight for each edge and produces, for each edge, a corresponding output pair by using Equation 14. The set of variable values jointly minimize the sum of the local cost function while remaining close to the incoming message set as balanced by each edge's weight. Proper use of the three weight classes lead to dramatic performance gains (as will be explained below) and is crucial for integration with higher-level knowledge. The logic that implements this minimization step may itself be a general optimization algorithm but also be customized to each cost function in which the custom minimization almost always leads to improvements in algorithm performance and is typically the bulk of the implementation effort. The concur phase combines incoming messages about each variable from all associated local cost functions and computes a single assignment value using a fixed logic routine (e.g., a weighted average). After each variable node has concurred, it is possible to extract this set of values as the present solution state.
Each message incorporates both an assignment value and an accumulation over previous disagreements between the value computed by a local cost function and the concurred variable value. Due to this disagreement variable term, each edge, even those connected to the same variable node, often communicates a different message (i.e., each edge has a different view of the variable as informed by an aggregation over local iteration history). The three-weight message-passing algorithm converges when outgoing messages from all variable nodes do not change significantly between two subsequent iterations.
All message passing data is local to each edge within the graphical model and it is trivial to coarsely parallelize the two main phases of the algorithm (i.e., all factors may be minimized in parallel and then all variables may be concurred upon in parallel). For complex cost functions, fine grained parallelization within the minimization routine may lead to additional performance gains.
The following description relates to the utility of non-standard weights (i.e., infinite weight and zero weight) on non-convex problems. In a first non-convex problem, a constraint satisfaction puzzle such as Sudoku is illustrated. The optimization algorithm (i.e., message-passing algorithm used with the ADMM algorithm) is allowed to propagate certain information in which infinite weights are useful. In a second non-convex problem, a circle packing scenario is illustrated. The optimization algorithm is allowed to ignore irrelevant constraints in which zero weights are useful.
FIG. 3A shows an unsolved puzzle 300 while FIG. 3B shows a solved puzzle 350 using the device 100 of FIG. 1 according to an exemplary embodiment. Specifically, the optimization functionality executed by the processor 105 may be used for a Sudoku puzzle in order to successfully solve the unsolved puzzle 300. Again, the infinite weights provide invaluable information for such a non-convex problem.
A Sudoku puzzle is a partially completed row-column grid of cells partitioned into N regions, each of size N cells to be filled in using a prescribed set of N distinct symbols such that each row, column, and region contains exactly one of each element of the set. A well-formed Sudoku puzzle has exactly one solution. The Sudoku puzzle is an example of an exact-cover constraint-satisfaction problem and is NP-complete when generalized to N×N grids. Generally, the Sudoku puzzle is a 9×9 grid containing nine 3×3 regions such as the unsolved puzzle 300. However, as will be explained in further detail below, the optimization algorithm of the exemplary embodiments is designed to generate the solution regardless of the size or configuration of the puzzle.
To represent an N×N square-in-square Sudoku puzzle as an optimization problem, O(N3) binary indicator variables and O(N2) hard constraints may be used. For all open cells (those that have not been supplied as “clues”), a binary indicator variable may be used, designated as v(row, column, digit), to represent each possible digit assignment. For example, the variables v(1, 3, 1), v(1, 3, 2), . . . v(1, 3, 9) represent that the cell in row 1, column 3 can take values 1 through 9. A hard “one-on” constraint may be applied to enforce digit distinctiveness: a—one-on constraint requires that a single variable is “on” (e.g., 1.0) and any remaining are “off” (e.g., 0.0). One-on constraints may be applied to four classes of variable sets:

    
    
        (1) ∀r∀c{v(row, col, dig): row=r, col=c}
        
            one digit assignment per cell
        
        
        (2) ∀r∀d{v(row, col, dig): row=r,dig d}
        
            one of each digit assigned per row
        
        
        (3) •c•d{(row, col, dig): col=c,dig d}
        
            one of each digit assigned per column
        
        
        (4) ∀s∀d{v(row, col, dig):sq(row,col)=s,dig=d}
        
            one of each digit assigned per square
        
        
    
    


The exemplary optimization algorithm integrates an implementation of infinite weights within the one-on minimizers that serve to reduce the search space of the problem instance. The introduction of messages that are certain using the three-weight message-passing algorithm improves the time-to-solution, particularly in puzzles where constraint propagation is sufficient to logically deduce most or all of the puzzle solution with no search being required.
To interpret the Sudoku puzzle as nodes of the graphical model, each hard one-on constraint is implemented as a cost function on the left side of the graphical model. For this class of constraint, minimizing Equation 8 involves a linear scan: select the sole “on” edge as that which is certain and “on” for ijt=∞ and nijt=1.0 or, in absence of such an edge, that with the greatest incoming message value and a standard weight.
Outgoing weights {right arrow over (ρ)}ijt default to ρ0 with three exceptions. In a first exception, if a single edge is certain and “on” for ijt=∞ and nijt=1.0, all outgoing assignments are certain. In a second exception, if all but a single incoming edge is certain and “off” for ijt=∞ and nijt=0.0, all outgoing assignments are certain. In a third exception, incoming certainty for an edge is maintained in its outgoing weight for ijt=∞{right arrow over (ρ)}ijt=∞.
Using the exemplary optimization algorithm in which the three-weight message-passing algorithm is utilized with the ADMM algorithm, Sudoku puzzles that represent the problem may be solved with increased efficiency, specifically from the introduction of the infinite weight. For example, experimental evidence of the use of the optimization algorithm has shown that an overall median improvement from using the exemplary optimization algorithm results in a 4.12× reduction in iterations with a maximum improvement of 61.21× on a single puzzle.
FIG. 4 shows a packing problem and solution 400 using the device 100 of FIG. 1 according to an exemplary embodiment. Specifically, the optimization functionality executed by the processor 105 may be used for a circle packing problem in order to successfully determine the solution. In this instance, the zero weights provide invaluable information for such a non-convex problem.
Circle packing is a problem of positioning a given number of congruent circles in such a way that the circles fit fully in a square without overlapping. A large number of circles makes finding a solution difficult due in part to the coexistence of many different circle arrangements with similar density. For example, as shown in FIG. 4, the packing may be rotated across either or both axes and the free circle in the upper-right corner (a “free circle” or a “rattle”) may be moved without affecting the density of the configuration.
To represent a circle-packing instance with N objects as an optimization problem, O(N) continuous variables and O(N2) constraints may be used in a substantially similar manner as discussed above with the Sudoku puzzle. Each object may have two variables: one representing each of its coordinates (or more generally, d variables for packing spheres in d dimensions). For each object, a single box-intersection constraint may be created that enforces that the objects stays within the box. Furthermore, for each pair of objects, a pairwise-intersection constraint may be created that enforces that no two objects overlap. The description below illustrates how the zero-weight messages in the three-weight message-passing algorithm used with the ADMM algorithm provides an improvement compared with a ADMM/DC approach and may be applied more generally to a variety of problems where local constraints that are “inactive” may otherwise send messages that would slow progress towards convergence.
Initially, both types of intersection constraints are implemented as cost functions on the left side of the graphical model 200. Box-intersection functions send messages of standard weight ρ0 for circles sending n messages outside the box to have an x position at the nearest box boundary. For those circles that are not outside the box, a zero-weight message is sent to stay at the present position such that {right arrow over (ρ)}ijt=0 and xijt=nijt. Pairwise-intersection functions are analogous: intersecting circles are sent standard weight messages reflecting updated x positions obtained by moving each circle along the vector connecting them such that Equation 14 is satisfied (e.g., if both circles send equal weight messages, they are moved an equal amount; if one has a standard weight and one has a zero weight, only the circle sending a zero-weight message is moved), while non-intersecting circles are sent zero-weight messages to remain at their present locations (xijt=nijt).
A comparison between iterations-to-convergence for a small number of circles (e.g., N=1 . . . 24) in a square between the three-weight message-passing algorithm used with the ADMM algorithm and the conventional ADDM algorithm (using only standard weights for message passing) has illustrated an improvement. For example, with four random conditions per N, ρ0=1, and α=0.01, for N<20, the three-weight message-passing algorithm improved performance of 42% of experimental trials of double or more with a median improvement in iterations being around 1.77×. However, for N≧20, the three-weight message-passing algorithm showed improvements of double or more on 90% of experimental trials with a median improvement being more than 116×. The results have shown that only a logarithmic growth in iterations occurred as N increased using the exemplary three-weight message-passing algorithm such that fewer than 1,220 iterations were required after N=901. However, the ADDM algorithm using only standard weights resulted in a quadratic increase to require more than 187,000 iterations. In the standard ADMM/DC approach, each circle is effectively being sent a message to stay where it is by every other circle in the system that is not overlapping with it and this may tremendously slow down convergence, particularly when N grows larger. By allowing for zero-weight messages from constraints that do not care about the location of a circle, the three-weight message-passing algorithm used with the ADMM algorithm becomes focused on those constraints that are relevant.
Although the examples described above relate to how using a particular weighted message provides an advantage in determining an optimization solution for a specific problem, those skilled in the art will understand how utilizing all three weighted messages in the three-weight message-passing algorithm used in the ADMM algorithm provide a combined improvement to solve convex and non-convex problems much faster than using DC algorithms and have important advantages over more widely used BP algorithms. Specifically, messages and beliefs in BP algorithms are probability distributions while ADMM uses a single value representing a current best guess for a variable; ADMM-based algorithms easily handle constraints that would be awkward to deal with in BP algorithms such as the constraint that circles cannot overlap in the circle-packing problem; and ADMM-based algorithms only converge to fixed points that satisfy all hard constraints in the problem whereas BP algorithms potentially converge to fixed points that fail to satisfy all the hard constraints. An application of the three-weight message-passing algorithm used in the ADMM algorithm will be discussed in further detail below in which the optimization algorithm may be used to solve a multi-agent trajectory planning problem.
The optimization algorithm using the three-weight message-passing algorithm applied to the ADMM algorithm may consider a variety of factors to further improve generating an optimization solution. Specifically, an approach to generating the optimization solution for quickly and flexibly developing hybrid cognitive capabilities that are efficient, scalable, and exploit knowledge may improve the speed and quality at which the optimization solution is found. The integration of knowledge with the exemplary optimization algorithm may result in improving expressiveness, efficiency, and scaling. The adjustment to the optimization algorithm by integrating knowledge will again be described with regard to a Sudoku puzzle and a circle packing problem.
The following description relates to an intermediate approach between highly specialized algorithms and uniform computation such that a hybrid architectural substrate seeks to leverage optimization over factor graphs or graphical models using the three-weight message-passing algorithm used with the ADMM algorithm as a general platform upon which to rapidly and flexibly develop diverse cognitive-processing modules. The integration of higher-level knowledge into the operation of the three-weight message-passing algorithm may utilize a variety of approaches. Two exemplary approaches include reasoner hierarchy and graph dynamics. These techniques are general and, when specialized for a particular problem, may lead to improved algorithm efficiency (iterations and wall-clock time), scaling, and expressiveness of constraints and heuristics. Accordingly, the processor 105 of the device 100 may further include a reasoner engine (not shown) as well as a graph dynamics engine (not shown), both engines providing data to the optimization engine 145 to generate the optimization solution.
With regard to reasoner hierarchy, the three-weight message-passing algorithm may be augmented to introduce a two-level hierarchy of local and global reasoners. Local reasoners may be implemented as a special class of factor within the program graph. Local reasoners may also be able to send/receive messages like other factors but incoming message values always reflect the concurred upon variable value. The default operation of the local reasoners is to send zero-weight messages (having no impact on the problem) but in particular circumstances, the local reasoners may also affect change through non-zero weight messages. Local reasoners may also have a Reason method that supports arbitrary logic. This class of reasoners may also be considered local due to its local view of the problem via connected variables (to thereby be executed concurrently). However, the added reasoning step affords communication with global reasoners. Global reasoners may be implemented as code modules via a single Reason method and are not connected to the problem graph. Instead, the global reasoners have a global view of the problem via access to the concurred values of any variable, as well as any information transmitted via local reasoners. Global reasoners may affect the problem via three main methods: (1) requesting that local reasoners send non-zero weighted messages, (2) detecting problem-specific termination conditions and halting iteration, and (3) modifying the problem graph.
A major reason for a hierarchy of reasoners is to exploit parallelism in order to better scale to large problems. Thus, where possible, local reasoners serve as a filtering step such that global reasoners need not inspect/operate on the full variable set. For example, as will be described in further detail below, with the Sudoku puzzle, the hierarchy yields an event-based discrimination network where the local reasoners pass along information about changes to possible cell states and a global reasoner implements relational logic that would be difficult and inefficient to represent within the graphical model of the problem. Local reasoners may filter for certain information (e.g., when weights are infinite), propagate implications to global reasoning, and express certain results in outgoing messages. Again with the Sudoku puzzle, propagating certainty via the reasoner hierarchy maintains real-time reactivity for very large problem instances by pruning unnecessary constraints.
Accordingly, the two-level reasoner hierarchy improves the exemplary three-weight message-passing algorithm along several dimensions. In a first example, global reasoners may implement arbitrary logic, including complex indexing structures, search algorithms, etc. Information may be extracted from the problem, optionally filtered through local reasoners; processed via the global reasoner; and then re-integrated through the interface of message, weight pairs in the local reasoners. Certainty weighting (i.e., infinite weights) allows for fluid integration with processing mechanisms that operate over symbolic representations. In a second example, the two-level hierarchy supports relational reasoning over the inherently prepositional graphical model representation without incurring combinatorial explosion. The global reasoner may incorporate richer representations such as rules, higher-order logics, explicit notions of uncertainty, perceptual primitives, etc. In a third example, operations of the local reasoners may be parallelized in the same way as factor minimization and variable-value concurrence in the problem graph. Furthermore, effective use of local filtering may greatly reduce the set of variables considered by global reasoners, thereby supporting scaling to large, complex problems.
With regard to graph dynamics, integrating knowledge may support four classes of graph dynamics by global reasoners: (1) adding/removing edges (variables nodes may be automatically pruned if there are no remaining incoming edges); (2) adding/removing factor nodes (removing a factor node has the effect of removing any incoming edges); (3) adding new variables; and (4) re-parameterizing factors. These actions enable adapting the representation of a single problem instance over time as well as reusing a single mechanism for multiple problem instances.
Removing graph edges and factors has two main effects: (1) potentially changing variable assignments and (2) improving runtime performance. First, if an edge is disconnected from a factor/variable node, the outgoing variable assignments are now no longer dependent upon that input, and therefore the objective cost may yield a different outcome. Second, while removing an edge is analogous to sending a zero-weight message, the system is no longer required to expend computation time and thus, wall-clock time, per iteration, may improve. That is, processing requirements may be decreased in an overall sense.
The ability to add and remove edges allows the graph to represent and reason about dynamically sized sets of variables. For example, in the Sudoku puzzle, the three-weight message-passing algorithm may consider a set of possible symbols for each cell and reduce its option set over time as logically certain assignments are made within the row/column/square. That is, infinite weights may be applied in which again, the processing requirements may be decreased.
Supporting factor re-parameterization may support numerous capabilities. First, it is possible to reflect incremental environmental changes without having to reconstruct the graph. It is also possible to reflect changes to the objective function which may come about due to environmental change, task transitions, dynamic agent preference/goals, etc. Additionally, repurposing existing factors helps keep memory costs stable that support scaling to large, complex problems.
Accordingly, graph dynamics improves the three-weight message-passing algorithm along several dimensions. In a first example, global reasoners may dynamically restructure and reconfigure the problem graph to reflect changes in the state of the environment, task structure, as well as agent preferences, goals, and knowledge. In a second example, changing edge connectivity supports dynamic sets of variables without the necessity of enumerating all possibilities. In a third example, performance of the three-weight message-passing algorithm iteration loop depends upon the size of the graph that may be dynamically maintained such as to represent only those factor nodes, variables, and edges that are necessary.
As discussed above, the knowledge integration in the three-weight message-passing algorithm may be applied for a variety of problems, convex and non-convex. That is, the knowledge integration manner may be used for determining a solution for a general optimization problem. To illustrate how the knowledge integration is used, exemplary problems such as Sudoku and circle-packing show how this modification to the three-weight message-passing algorithm further improves convergence and further reduces processing requirements to generate the optimization solution.
With regard to a Sudoku puzzle, a local reasoner may be added per cell of the puzzle with the express purpose of maintaing the set of possible digit assignments (i.e., the values for which there was not a certain message for the “off” binary indicator variable discussed above) and, when a possibility was removed, communicate this information to a global reasoner. The global reasoner may utilize this information to improve solver efficiency and performance by utilizing graph dynamics to reduce the problem search space. Specifically, when a possibility is no longer viable logically, the global reasoner may remove four edges (those that connected the binary indicator variable to the constraints of the cell, row, column, and square) as well as any factors that are left without an option in the process. The graph dynamics reduce the problem-graph size by removing constraints as they became unnecessary and thereby improved iteration time of the second phase of numerical solving. Experimental data has shown that the integration of knowledge into the three-weight message-passing algorithm further improves upon the already improved manner of utilizing the three-weight message-passing algorithm in the ADMM algorithm for solving Sudoku puzzles.
Accordingly, with regard to Sudoku puzzles, local reasoners were used to create an efficient and scalable discrimination network of possibility-set changes for consumption by a global reasoner. A global reasoner responded to changes in possibility-sets in order to implement graph dynamics. The graph dynamics pruned the problem set size (more than an order of magnitude for large puzzles) and had the result of improving efficiency and scaling, including reactive iterations for the largest puzzles.
With regard to circle packing, the knowledge integration manner may be used for two distinct purposes: (1) drastically improving scaling via R-tree integration and (2) macro movements via human assistance. As discussed above, using zero-weight messages for inactive constraints improves iterations-to-convergence for circle packing. However, O(N2) intersection constraints are still required and thus, iteration-loop time (msec/iteration) still prevents scaling to large numbers of circles. However, R-trees provide an efficient detection of intersections of spatial objects. Since the three-weight message-passing algorithm tends to move circles in a local fashion, a graph dynamics global reasoner that integrates an R-tree to add/remove intersection constraints for each iteration may be used. Initially, the R-tree may be updated with the locations of all circles via concurred variable values. Axis-aligned bounding boxes (AABBs) may describe object locations and add a 5% “neighborhood” buffer to each box. All object intersections may be queried and categorize them as either existing or new based upon a current graphical model. The algorithm may iterate over factors in the problem graph such as (1) for existing intersections with corresponding factors, perform nothing and (2) for those factors that are no longer intersecting (or in the neighborhood buffer), these factors are removed from the graph and added to a “pool” of constraints. For all new intersections, a constraint may be added to the graphical model to reconfigure past factors if available in the pool. Therefore, only those objects that are in the graph need to be updated in the R-tree since those are the only coordinates the three-weight message-passing algorithm have altered. The knowledge integrated three-weight message-passing algorithm scales with the number of active factors and intersections which is small compared to the full O(N2) intersection constraint set. Old factors may be pooled to re-parameterize them dynamically in order to bound memory consumption. Experimental data has shown that the integration of knowledge into the three-weight message-passing algorithm further improves upon the already improved manner of utilizing the three-weight message-passing algorithm in the ADMM algorithm for solving circle packing problems.
It is noted that because the iteration dynamics of circle packing includes the three-weight message-passing algorithm to make local changes to circle positions, it may be difficult for the knowledge integrated three-weight message-passing algorithm to make larger, coordinated movements to escape local minima. Accordingly, this algorithm may incorporate two sets of reasoners for human assistance (via the I/O device 120) into the low-level circle-packing optimization problem. The first implementation requires only a single local reasoner where a user clicks on and drags a circle in an interface. The coordinates of the resulting position may be sent as weighted messages to the circle's coordinate variables (such that graph dynamics are used to connect the local-reasoner factor on demand). The inclusion of such “knowledge” significantly improves upon the already improved three-weight message-passing algorithm. For example, the global reasoner may extract meta-data from the R-tree intersections to identify the circle with the largest current overlap in which the user may identify a free or open space. The local reasoner may transport (via weighted messages) the “distressed” circle to the “vacant” region. Accordingly, this may decrease iterations-to-convergence.
Accordingly, with regard to circle packing problems, a local reasoner may integrate manual inputs from a user. A global reasoner may integrate an R-tree with the main iteration loop for improved efficiency and problem-size scaling as well as to inform user-assisted packing via an integration of intersection metadata. Graph dynamics may maintain the set of constraints as informed by the R-tree. These constraints may re-parameterize to limit memory usage.
As discussed above, one application for using the three-weight message-passing algorithm with the ADMM algorithm relates to planning trajectories for a plurality of agents such as cameras, robots in an amusement park, etc. That is, the exemplary optimization algorithm may compute collision-free global trajectories for p agents with specified initial and final configurations. In a particular embodiment, a specialization of the exemplary optimization algorithm may be used for local motion planning by solving the problem of joint optimization in velocity space. As will be described in further detail below, the three-weight message-passing algorithm improves upon an efficiency in generating an optimization solution for this problem. For example, the advantages described above in the examples regarding Sudoku and circle packing may be incorporated into generating the optimization solution for a multi-agent trajectory planning problem.
Robot navigation relies on at least three subtasks. The subtasks may include localization, mapping, and motion planning. Motion planning may be described as an optimization problem in which the lowest cost path (i.e., trajectory) is computed between an initial and final configuration. Trajectory planning for multiple agents is an important problem in robotics, computer animation, and crowd simulation. Conventional methods of generating the trajectory planning solution prioritize agents and compute trajectories sequentially which ultimately lead to suboptimal solutions. In contrast, the exemplary optimization algorithm plans for all agents simultaneously by taking advantage of the parallelization feature.
The exemplary optimization algorithm may perform the task of determining a global trajectory planning in the following manner. Given p agents of different radii {ri}i−1p with a given desired initial position {x(0)}i=1p and final position {xi(T)}i=1p, along with a cost function over trajectories, collision-free trajectories may be computed for all agents that minimize the cost function. That is, a set of intermediate points {xi(t)}i=1p may be found, where t∈(0,T), that satisfy the “hard” collision-free constraints of ∥xi(t)−xj(t)∥>ri+rj for all i, j, and t and that minimizes the cost function.
The manner of utilizing the exemplary optimization algorithm searches for a solution within the space of piece-wise linear trajectories where the trajectory of an agent is completely specified by a set of positions at a fixed set of time instants {ts}s=0a. These time instants may be defined as break-points which are the same for all agents to simplify the description below. All other intermediate points of the trajectories may be computed by assuming that each agent moves with constant velocity in between break-points. That is, if t1 and t2>t1 are consecutive break-points, then




  
    
      x
      i
    
    ⁡
    
      (
      t
      )
    
  
  =
  
    
      
        1
        
          
            t
            2
          
          -
          
            t
            1
          
        
      
      ⁢
      
        (
        
          
            
              (
              
                
                  t
                  2
                
                -
                t
              
              )
            
            ⁢
            
              
                x
                i
              
              ⁡
              
                (
                
                  t
                  1
                
                )
              
            
          
          +
          
            
              (
              
                t
                -
                
                  t
                  1
                
              
              )
            
            ⁢
            
              
                x
                i
              
              ⁡
              
                (
                
                  t
                  2
                
                )
              
            
          
        
        )
      
      ⁢
      
        
      
      ⁢
      for
      ⁢
      
        
      
      ⁢
      t
    
    ∈
    
      
        [
        
          
            t
            1
          
          ,
          
            t
            2
          
        
        ]
      
      .
    
  




Along with the set of initial and final configurations, the number of interior break-points (n−1) may be used as an input. Increasing n yields trajectories that are more flexible and smooth with possibly higher quality, but at the cost of more computation.

As will be described in further detail below, the exemplary optimization algorithm provides various features in generating the optimization solution. The global-planning task may be formulated as a decomposable optimization problem. The resulting sub-problems may be solved exactly and efficiently despite non-convexity and exemplary optimization algorithm coordinate the solutions of the sub-problems using message-passing. As described above, the three-weight message-passing algorithm is easily parallelized, does not require parameter tuning, and has good scalability with p. Within the decomposable framework, different minimizations of sub-problems, defined as “minimizers,” ensure the trajectories satisfy a separate criterion. The exemplary optimization algorithm is flexible and may consider different combinations of minimizers. A particularly crucial minimizer ensures that there are no inter-agent collisions. Other minimizers allow for finding trajectories for the agents with minimal total energy, avoiding static obstacles, and/or imposing dynamic constraints such as maximum/minimum agent velocity. The method may also be specialized to perform local planning by solving the problem of joint optimization in the velocity space.
The exemplary optimization algorithm may relate to a variety of different minimizers for the problem of multi-agent trajectory planning. Once a formulation for an unconstrained optimization without static obstacles, a formulation for unconstrained optimization with static obstacles, and a message-passing formulation has been defined, an agent-agent collision minimizer, an agent-obstacle collision minimizer, and a minimum energy and maximum (or minimum) velocity minimizer may be solved for the problem.
FIG. 5A shows a small part of graphical model 500 as applied to a multi-agent trajectory planning problem according to an exemplary embodiment. FIG. 5B shows input and output variables (I/O) 550 for each minimizer block of the graphical model of FIG. 5A according to an exemplary embodiment. As illustrated, the graphical model 500 may be a bipartite graph with one minimizer-node (g1, g2, and g3) on the left for each function making up the overall objective function and one equality-node (E1, E2, and E3) on the right per variable in the problem. It is noted that most embodiments would typically have many more minimizer nodes on the left and equality nodes on the right (e.g. hundreds of nodes on each side); the exemplary embodiment in FIG. 5A only shows three on each side for illustrative purposes. The I/O variables 550 illustrate the input and output variables for each minimizer block g1, g2, and g3. Specifically, based upon the graphical model 500, the minimizer block g1 includes input variables n1,1,1,1 and n1,3,1,3 with output variables x1,1,{right arrow over (ρ)}1,1 and x1,3,{right arrow over (ρ)}1,3; the minimizer block g2 includes input variables n2,1,2,1, n2,2, 2,2, and n2,3,2,3 with output variables x2,1,{right arrow over (ρ)}2,1, x2,2,{right arrow over (ρ)}2,2 and x2,3,{right arrow over (ρ)}2,3; and the minimizer block g3 includes input variable n3,3,3,3 with output variable x3,3,{right arrow over (ρ)}3,3.
The optimization problem for the multi-agent trajectory planning illustrated in FIG. 5A would be minx∈R3f(x)=min{x1,x2,x3} ƒ1(x1,x3)+ƒ2(x1, x2, x3)+ƒ3(x3), where ƒi(·)∈R∪{∞). The functions may represent soft costs or hard equality or inequality constraints. The three-weight message-passing algorithm may be used to solve this optimization problem iteratively by passing messages on the bipartite graph of graphical model 500 in which one minimizer node is represented per function ƒb, one equality node is represented per variable xj and an edge (b,j), connecting b and j, if ƒb depends on xj.
The way that minimizer nodes solve sub-problems and coordinate for robot trajectory optimization is similar to the description already made for circle packing. The minimizer node g1 may solve a small optimization problem over its local variables x1 and x3. The estimates x1,1 and x1,3 may then combine with running sums of the differences between the minimizer estimates and the equality-node consensus estimates to obtain messages m1,1 and m1,3 on each neighboring edge that are sent to the neighboring equality-nodes along with corresponding certainty weights {right arrow over (ρ)}1,1 and {right arrow over (ρ)}1,3. All other minimizers (e.g., g2 and g3) may also be solved accordingly.
The equality nodes E1, E2, and E3 may then receive these local messages and weights and produce consensus estimates for all variables by computing an average of the incoming messages, weighted by the incoming certainty weights {right arrow over (ρ)}. From these consensus estimates, correcting messages may be computed and communicated back to the minimizers to help them reach a consensus. A certainty weight for the correcting messages  may also be communicated back to the minimizers. For example, the minimizer g1 may receive correcting messages n1,1 and n1,3 with corresponding certainty weights 1,1 and 1,3.
As in the previous description of the three-weight message-passing optimization algorithm applied to circle packing, when producing new local estimates, the bth minimizer node may compute its local estimates {xj} by choosing a point that minimizes the sum of the local function ƒb and weighted squared distance from the incoming messages where ties are broken randomly with the following:




  
    
      
        
          
            {
            
              x
              
                b
                ,
                j
              
            
            }
          
          j
        
        =
        
          
            
              g
              b
            
            ⁡
            
              (
              
                
                  
                    {
                    
                      n
                      
                        b
                        ,
                        j
                      
                    
                    }
                  
                  j
                
                ,
                
                  
                    {
                    
                      
                        ρ
                        ←
                      
                      
                        b
                        ,
                        j
                      
                      k
                    
                    }
                  
                  j
                
              
              )
            
          
          ≡
          
            
              
                arg
                ⁢
                
                  
                
                ⁢
                min
              
              
                
                  {
                  
                    x
                    j
                  
                  }
                
                j
              
            
            [
            
              
                
                  f
                  b
                
                ⁡
                
                  (
                  
                    
                      {
                      
                        x
                        j
                      
                      }
                    
                    j
                  
                  )
                
              
              +
              
                
                  1
                  2
                
                ⁢
                
                  
                    ∑
                    j
                  
                  ⁢
                  
                    
                      
                        
                          ρ
                          ←
                        
                        
                          b
                          ,
                          j
                        
                      
                      ⁡
                      
                        (
                        
                          
                            x
                            j
                          
                          -
                          
                            n
                            
                              b
                              ,
                              j
                            
                          
                        
                        )
                      
                    
                    2
                  
                
              
            
            ]
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          A
        
        )
      
    
  




where { }j and Σj run over all equality nodes connected to b. In the three-weight message-passing algorithm, the certainty weights {{right arrow over (ρ)}b,i} that this minimizer outputs must be 0 (uncertain); ∞ (certain); or ρ0, set to some fixed weight. The logic for setting weights from minimizer nodes depends on the problem and with trajectory planning problems, only 0 or ρ0 weights may be used. The zero weight of the three-weight message-passing algorithm enables equality nodes to ignore inactive constraints and traverse the search space much faster. It is again noted that all minimizers may operate simultaneously and the same is true for the consensus calculation performed by each equality node, thereby being easy to parallelize.

To begin the decomposition specific to the global trajectory planning optimization problem, the variables to be defined in the optimization problem are defined. The variables are defined as positions {xi(s)}i∈[p] where the trajectories are indexed by i and break-points are indexed by a discrete variable s taking values between 1 and n−1. Accordingly, the points of the trajectories are not tracked by a continuous-time variable taking values in [0,T]. It is noted that {xi(0)}i∈[p] and {xi(n)}i∈[p] are the initial and final configuration, sets of fixed values, not variables to optimize. That is, these values are received as an initial input for the optimization problem and are not altered during the computation.
In a first optimization, a formulation for unconstrained optimization without static obstacles may be computed. In terms of these variables the non-collision constraints may be defined as:

∥(αxi(s+1)+(1−α)xi(s))−(αxj(s+1)+(1−α)xj(s))∥≧ri+rj  (Equation B)

for all i,j∈[p], s∈{0, . . . , m−1} and α∈[0,1]. The parameter α is used to trace out the constant-velocity trajectories of agents i and j between break-points s+1 and s. The parameter s has no units and is a normalized time rather than an absolute time. If t1 is the absolute time of the break-point with integer index s and t2 is the absolute time of the break-point with integer index s+1 and t parameterizes the trajectories in absolute time, then α=(t−t1)/(t2−t1). It is noted that in the above formulation, absolute time does not appear and any solution is simply a set of paths that, when traveled by each agent at constant velocity between break-points, leads to no collisions. When converting this solution into trajectories parameterized by absolute time, the break-points do not need to be chosen uniformly spaced in absolute time.

The constraints represented in Equation B may be formally incorporated into an unconstrained optimization problem with the following form:




  
    
      
        
          
            min
            
              
                {
                
                  
                    x
                    i
                  
                  ⁡
                  
                    (
                    s
                    )
                  
                
                }
              
              
                i
                ,
                s
              
            
          
          ⁢
          
            
              f
              cost
            
            ⁡
            
              (
              
                
                  {
                  
                    
                      x
                      i
                    
                    ⁡
                    
                      (
                      s
                      )
                    
                  
                  }
                
                
                  i
                  ,
                  s
                
              
              )
            
          
        
        +
        
          
            ∑
            
              s
              =
              0
            
            
              n
              -
              1
            
          
          ⁢
          
            
              ∑
              
                i
                >
                j
              
            
            ⁢
            
              
                f
                
                  
                    r
                    i
                  
                  ,
                  
                    r
                    j
                  
                
                coll
              
              ⁡
              
                (
                
                  
                    
                      x
                      i
                    
                    ⁡
                    
                      (
                      s
                      )
                    
                  
                  ,
                  
                    
                      x
                      i
                    
                    ⁡
                    
                      (
                      
                        s
                        +
                        1
                      
                      )
                    
                  
                  ,
                  
                    
                      x
                      j
                    
                    ⁡
                    
                      (
                      s
                      )
                    
                  
                  ,
                  
                    
                      x
                      j
                    
                    ⁡
                    
                      (
                      
                        s
                        +
                        1
                      
                      )
                    
                  
                
                )
              
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          C
        
        )
      
    
  




where {xi(0)}p and {xi(n)}p are constants rather than optimization variables and where the function ƒcost is a function that represents some cost to be minimized (e.g., the integrated kinetic energy or the maximum velocity over all the agents) and the function ƒr,r′call is defined as:

ƒr,r′call(x,x,x′,x′,x′)=J(∥α(x−x′)+(1−α)(x−x′)∥≧r+r′∀α∈[0,1])  (Equation D)

In this section, x and x represent the position of an arbitrary agent of radius r at two consecutive break-points and x′ and x′ represent the position of a second arbitrary agent of radius r′ at the same break-points. The expression J(·) equals 0 whenever its argument, a clause, is true and takes equals +∞ otherwise. Therefore, an infinite cost is paid in ƒr,r′call whenever there is a collision and zero is paid otherwise.

In Equation C, the cost function ƒcost(·) may be used to enforce a preference for trajectories satisfying specific properties. For example, there may be a preference for trajectories for which the total kinetic energy spent by the set of agents is small. Accordingly, the cost function may be defined as ƒCcost(x, x)=C∥x−x∥2 to result in the following:




  
    
      
        
          
            f
            cost
          
          ⁡
          
            (
            
              
                {
                
                  
                    x
                    i
                  
                  ⁡
                  
                    (
                    s
                    )
                  
                
                }
              
              
                i
                ,
                s
              
            
            )
          
        
        =
        
          
            1
            pn
          
          ⁢
          
            
              ∑
              
                i
                =
                1
              
              p
            
            ⁢
            
              
                ∑
                
                  s
                  =
                  0
                
                
                  n
                  -
                  1
                
              
              ⁢
              
                
                  f
                  
                    C
                    
                      i
                      ,
                      s
                    
                  
                  cost
                
                ⁡
                
                  (
                  
                    
                      
                        x
                        i
                      
                      ⁡
                      
                        (
                        s
                        )
                      
                    
                    ,
                    
                      
                        x
                        i
                      
                      ⁡
                      
                        (
                        
                          s
                          +
                          1
                        
                        )
                      
                    
                  
                  )
                
              
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          E
        
        )
      
    
  




where the coefficients {Ci,s} account for agents with different masses, different absolute-time intervals between break-points, or different preferences regarding which agents are to be less active and which agents are allowed to move faster.

To further simplify the computation, trajectories in which agents move faster than a certain amount may be excluded but without distinguishing among all remaining trajectories. Accordingly, the following may be used:

ƒCcost(x,x)=J(∥x−x∥≦C)  (Equation F)

In this case, associating each break-point to a time instant, the coefficients {Ci,s} in Equation E may represent different limits on the velocity of different agents between different sections of the trajectory. If all agents are to be forced to have a minimum velocity, the inequality of Equation F may be reversed.

In a second optimization, a formulation for unconstrained optimization with static obstacles may be computed. In many scenarios, agents should also avoid collisions with static obstacles. Given two points in space xL and xR, the optimization problem may include parameters to forbid all agents from crossing a line segment from xL to xR by adding the following term to Equation C: Σi=1pΣs=0n−1ƒXL,xR,riwall(xi(s),xi(s+1)). Since ri is the radius of agent i, the following may be used:

ƒxL,xR,rwall(x,x)=J∥(αx+(1−α)x)−(βxR+)(1−β)xL)∥≧r  (Equation G)

for all α,β∈[0,1]. It is noted that ƒcall may be expressed using ƒwall. Specifically, the following may result:

ƒr,r′call(x,x,x′,x′)=ƒ0,0,r+r′wall(x′−x,x′−x′)   (Equation H)

This modification may be used subsequently to express the minimizer associated with agent-agent collisions using the minimizer associated with agent-obstacle collisions, as will be described in further detail below.

When agents move in a plane (e.g., xi(s)∈R2 for all s+1∈[n+1]), being able to avoid collisions with a general static line segment allows to automatically avoid collisions with multiple static obstacles of arbitrary polygonal shape. However, for simplicity, the numerical experiments only consider agents in the plane and so, only the minimizer block for wall collision for a 2-dimensional world is described. However, those skilled in the art will understand that the computations described herein may also be modified to incorporate higher dimensions.
In a third optimization, a formulation for message-passing may be computed. To solve Equation C using the three-weight message-passing algorithm, the topology of the bipartite graphical model 500 associated with the unconstrained Equation C is specified and the operation performed by every minimizer (i.e., the {right arrow over (ρ)} weight update logic and the x variable update equations). It is noted that the choice of initial message values and weights and internal parameters are also included and will be described in further detail below.
To be concrete, the cost function is assumed to have the form of Equation E. The unconstrained Equation C then indicates that the global objective function is the sum of np(p+1)/2 terms: np(p−1)/2 functions ƒcall and np functions ƒCcall. These functions involve a total of (n+1)p variables out of which only (n−1)p are free (since the initial and final configurations are fixed). Correspondingly, the graphical model 500 along which messages are passed has np(p+1)/2 minimizer nodes (left side of graphical model 500) that connect to the (n+1)p equality nodes (right side of graphical model 500). Specifically, the equality node associated with the break-point variable xi(s),n>s>0 is connected to 2(p−1) different gcall minimizer nodes and two different gCcost minimizer nodes. If s=0 or s=n, the equality node only connects to half as many gcall nodes and gCcost nodes.
The following provides a description of the different minimizers g1, g2, and g3. Specifically, the minimizers include an agent-agent collision minimizer, an agent-obstacle collision minimizer, and a minimum energy and maximum (minimum) velocity minimizer. It is noted that every minimizer is a special case of Equation A.
In a first minimizer, the agent-agent collision minimizer may be a minimizer associated with the functions ƒcall that may be denoted by gcall. This minimizer receives as parameters the radius r and r′ of the two agents whose collision it is avoiding. The minimizer takes as input a set of incoming n messages {n,n,n′,n′} and associated  weights, {,,,}. The minimizer outputs a set of updated x variables according to the following:




  
    
      
        
          
            g
            coll
          
          ⁡
          
            (
            
              
                n
                _
              
              ,
              
                n
                _
              
              ,
              ,
              ,
              
                
                  ρ
                  _
                
                ←
              
              ,
              
                
                  ρ
                  _
                
                ←
              
              ,
              ,
              ,
              r
              ,
            
            )
          
        
        =
        
          
            ⁢
            ⁢
            
              (
              
                
                  x
                  _
                
                ,
                
                  x
                  _
                
                ,
                ,
                
                  
                    x
                    ′
                  
                  _
                
              
              )
            
          
          +
          
            
              
                
                  ρ
                  ←
                
                _
              
              2
            
            ⁢
            
              
                
                
                  
                    x
                    _
                  
                  -
                  
                    n
                    _
                  
                
                
              
              2
            
          
          +
          
            
              
                
                  ρ
                  _
                
                ←
              
              2
            
            ⁢
            
              
                
                
                  
                    x
                    _
                  
                  -
                  
                    n
                    _
                  
                
                
              
              2
            
          
          +
          
            
              2
            
            ⁢
            
              
                
                
                  -
                
                
              
              2
            
          
          +
          
            
              2
            
            ⁢
            
              
                
                
                  -
                
                
              
              2
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          I
        
        )
      
    
  




Messages n and n come from the two equality nodes associated with the positions of one of the agents at two consecutive break points and n′ and n′ from the corresponding equality nodes for the other agent.

The update logic for the weights {right arrow over (ρ)} for this minimizer may be determined with the following. If the trajectory from n to n for an agent of radius r does not collide with the trajectory from n′ to n′ for an agent of radius r′, then all the outgoing weights {right arrow over (ρ)} may be set to zero. Otherwise, the outgoing weights {right arrow over (ρ)} may be set to ρ0. The outgoing zero weights indicate to the receiving equality nodes in the graphical model 500 that the collision constraint for this pair of agents is inactive and that the values it receives from this minimizer node should be ignored when computing the consensus values of the receiving equality nodes.
In a second minimizer, the agent-obstacle collision minimizer may be the minimizer for the function ƒwall that may be denoted as gwall. This minimizer may be parameterized by the obstacle position {xL,xR} as well as the radius of the agent that needs to avoid the obstacle. This minimizer may receive two n messages {n,n} and corresponding weights {,} from the equality nodes associated with two consecutive positions of an agent that needs to avoid the obstacle. The minimizer may output the x-variables as defined by the following:




  
    
      
        
          
            g
            wall
          
          ⁡
          
            (
            
              
                n
                _
              
              ,
              
                n
                _
              
              ,
              r
              ,
              
                x
                L
              
              ,
              
                x
                R
              
              ,
              
                
                  ρ
                  _
                
                ←
              
              ,
              
                
                  ρ
                  _
                
                ←
              
            
            )
          
        
        =
        
          
            
              
                arg
                ⁢
                
                  
                
                ⁢
                min
              
              
                {
                
                  
                    x
                    _
                  
                  ,
                  
                    x
                    _
                  
                
                }
              
            
            ⁢
            
              
                f
                
                  
                    x
                    L
                  
                  ,
                  
                    X
                    R
                  
                  ,
                  r
                
                wall
              
              ⁡
              
                (
                
                  
                    x
                    _
                  
                  ,
                  
                    x
                    _
                  
                
                )
              
            
          
          +
          
            
              
                
                  ρ
                  ←
                
                _
              
              2
            
            ⁢
            
              
                
                
                  
                    x
                    _
                  
                  -
                  
                    n
                    _
                  
                
                
              
              2
            
          
          +
          
            
              
                
                  ρ
                  _
                
                ←
              
              2
            
            ⁢
            
              
                
                
                  
                    x
                    _
                  
                  -
                  
                    n
                    _
                  
                
                
              
              2
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          J
        
        )
      
    
  




When agents move in the plane (as discussed above, a two-dimensional plane), this minimizer may be solved by reformulating the optimization in Equation 10 as a mechanical problem involving a system of springs that are solved exactly and efficiently.

The update logic for the {right arrow over (ρ)} weights is similar to that of the gcall minimizer. Specifically, if an agent of radius r going from n and n does not collide with the line segment from xL to xR, then all outgoing weights may be set to zero because the constraint is inactive. Otherwise, all the outgoing weights {right arrow over (ρ)} may be set to ρ0. It is noted that from Equation H, the agent-agent minimizer gcall may be expressed using gwall. Specifically, the following may result: gcall(n, n, n′, n′, , , , , r, r′)=M2gwall (M1·{n, n, n′, n′, , , , , r, r′}) for a constant rectangular matrix M1 and a matrix M2 that depend on {n, n, n′, n′, , , , }.
In a third minimizer, the minimum energy and maximum velocity (or minimum velocity) minimizer may be the minimizer may be associated with the functions ƒcost that is denoted by gcost and receives as input, two n messages {n,n} and corresponding weights {,} when ƒcost is decomposed as in Equation E. The messages may come from two equality nodes associated with two consecutive positions of an agent. The minimizer may also be parameterized by a cost factor c. Accordingly, this minimizer may output a set of updated x messages defined with the following:




  
    
      
        
          
            g
            cost
          
          ⁡
          
            (
            
              
                n
                _
              
              ,
              
                n
                _
              
              ,
              
                
                  ρ
                  _
                
                ←
              
              ,
              
                
                  ρ
                  _
                
                ←
              
              ,
              c
            
            )
          
        
        =
        
          
            
              
                arg
                ⁢
                
                  
                
                ⁢
                min
              
              
                {
                
                  
                    x
                    _
                  
                  ,
                  
                    x
                    _
                  
                
                }
              
            
            ⁢
            
              
                f
                c
                cost
              
              ⁡
              
                (
                
                  
                    x
                    _
                  
                  ,
                  
                    x
                    _
                  
                
                )
              
            
          
          +
          
            
              
                
                  ρ
                  ←
                
                _
              
              2
            
            ⁢
            
              
                
                
                  
                    x
                    _
                  
                  -
                  
                    n
                    _
                  
                
                
              
              2
            
          
          +
          
            
              
                
                  ρ
                  _
                
                ←
              
              2
            
            ⁢
            
              
                
                
                  
                    x
                    _
                  
                  -
                  
                    n
                    _
                  
                
                
              
              2
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          K
        
        )
      
    
  




The update logic for the {right arrow over (ρ)} weights of the minimum energy minimizer may always set all outgoing weights {right arrow over (ρ)} to ρ0. For the update logic of the {right arrow over (ρ)} weights of the maximum velocity minimizer, the following rule may be used: if ∥n−n∥≦c, then all outgoing weights may be set to zero, while otherwise, the outgoing weights may be set to ρ0. For the update logic of the {right arrow over (ρ)} weights of the minimum velocity minimizer, the following rule may be used: if ∥n−n∥≧c, then all the outgoing weights may be set to zero, while otherwise, the outgoing weights may be set to ρ0.
The following description provides how the exemplary optimization algorithm in which a three-weight message-passing algorithm used with the ADMM algorithm is configured to generate an optimization solution for local trajectory planning based on velocity obstacles. Given the current position {xi(0)}i∈[p] and radius {ri} of all agents that are velocity obstacles, a new velocity command may be computed jointly for all agents minimizing the distance to their preferred velocity {viref}i∈[p]. This new velocity command is used to guarantee that the trajectories of all agents remain collision free for at least a time horizon τ. New collision free velocities are computed every ατ seconds where α<1 until all agents reach their final configuration. Assuming an obstacle-free environment and first order dynamics, the collision-free velocities may be given by




  
    
      
        minimize
        
          
            {
            
              v
              i
            
            }
          
          
            i
            ∈
            
              
              p
              
            
          
        
      
      ⁢
      
        
          ∑
          
            i
            ∈
            
              [
              p
              ]
            
          
        
        ⁢
        
          
            C
            i
          
          ⁢
          
            
              
              
                
                  v
                  i
                
                -
                
                  v
                  i
                  ref
                
              
              
            
            2
          
          ⁢
          
            s
            .
            t
            .
            
              
            
            ⁢
            
              
              
                
                  (
                  
                    
                      
                        x
                        i
                      
                      ⁡
                      
                        (
                        0
                        )
                      
                    
                    +
                    
                      
                        v
                        i
                      
                      ⁢
                      t
                    
                  
                  )
                
                -
                
                  (
                  
                    
                      
                        x
                        j
                      
                      ⁡
                      
                        (
                        0
                        )
                      
                    
                    +
                    
                      
                        v
                        j
                      
                      ⁢
                      t
                    
                  
                  )
                
              
              
            
          
        
      
    
    ≥
    
      
        r
        i
      
      +
      
        
          r
          j
        
        ⁢
        
          ∀
          
            i
            ∈
            
              [
              p
              ]
            
          
        
      
    
  
  ,
  
    t
    ∈
    
      [
      
        0
        ,
        τ
      
      ]
    
  




Since the velocities {vi}i∈[p] are related linearly to the final position of each object after τ seconds, a simple change of variables as {xi(τ)}i∈[p] enables the above problem to be reformulated as the following:





  
    
      
        
          
            
              minimize
              
                
                  {
                  
                    v
                    i
                  
                  }
                
                
                  i
                  ∈
                  
                    
                    p
                    
                  
                
              
            
            ⁢
            
              
                ∑
                
                  i
                  ∈
                  
                    [
                    p
                    ]
                  
                
              
              ⁢
              
                ⁢
                
                  
                    
                    
                      
                        x
                        i
                      
                      -
                      
                        x
                        i
                        ref
                      
                    
                    
                  
                  2
                
                ⁢
                
                  s
                  .
                  t
                  .
                  
                    
                  
                  ⁢
                  
                    
                    
                      
                        
                          (
                          
                            1
                            -
                            α
                          
                          )
                        
                        ⁢
                        
                          (
                          
                            
                              
                                x
                                i
                              
                              ⁡
                              
                                (
                                0
                                )
                              
                            
                            -
                            
                              
                                x
                                j
                              
                              ⁡
                              
                                (
                                0
                                )
                              
                            
                          
                          )
                        
                      
                      +
                      
                        α
                        ⁡
                        
                          (
                          
                            
                              x
                              i
                            
                            ⁢
                            
                              x
                              j
                            
                          
                          )
                        
                      
                    
                    
                  
                
              
            
          
          ≥
          
            
              r
              i
            
            +
            
              
                r
                j
              
              ⁢
              
                ∀
                
                  
                    j
                    >
                    i
                  
                  ∈
                  
                    [
                    p
                    ]
                  
                
              
            
          
        
        ,
        
          α
          ⁡
          
            [
            
              0
              ,
              1
            
            ]
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          L
        
        )
      
    
  




where Ci′=Ct/τ2,xiref=xi(0)+virefτ and the τ in xi(τ) is dropped. This may be extended to account for collisions with the static line segments {xRk, xLk}k, thereby formulating an unconstrained form using the functions ƒcost, ƒcall, and ƒwall as the following:





  
    
      
        
          ⁢
          
            
              ∑
              
                i
                ∈
                
                  [
                  p
                  ]
                
              
            
            ⁢
            
              ⁢
              
                (
                
                  
                    x
                    i
                  
                  ,
                  
                    x
                    i
                    ref
                  
                
                )
              
            
          
        
        +
        
          
            ∑
            
              i
              >
              j
            
          
          ⁢
          
            
              f
              
                
                  r
                  i
                
                ,
                
                  r
                  j
                
              
              coll
            
            ⁡
            
              (
              
                
                  
                    x
                    i
                  
                  ⁡
                  
                    (
                    0
                    )
                  
                
                ,
                
                  x
                  i
                
                ,
                
                  
                    x
                    j
                  
                  ⁡
                  
                    (
                    0
                    )
                  
                
                ,
                
                  x
                  j
                
              
              )
            
          
        
        +
        
          
            ∑
            
              i
              ∈
              
                [
                p
                ]
              
            
          
          ⁢
          
            
              ∑
              k
            
            ⁢
            
              
                f
                
                  
                    x
                    Rk
                  
                  ,
                  
                    x
                    Lk
                  
                  ,
                  
                    r
                    i
                  
                
                wall
              
              ⁡
              
                (
                
                  
                    
                      x
                      i
                    
                    ⁡
                    
                      (
                      0
                      )
                    
                  
                  ,
                  
                    x
                    i
                  
                
                )
              
            
          
        
      
    
    
      
        (
        
          Equation
          ⁢
          
            
          
          ⁢
          M
        
        )
      
    
  




It is noted that {xi(0)}i and {xiref}i are constants, not variables being optimized.

Given this formulation, the three-weight message-passing algorithm may be used to solve the optimization problem. All corresponding minimizers are special cases of minimizers derived in the manner discussed above for global trajectory planning. Experimental testing has shown that using the exemplary optimization algorithm provides improvements to conventional manners such as a Mixed Integer Quadratic Programming (MIQP) approach. Specifically, the time for convergence (e.g., when there is a consensus) to determine the optimization solution for the problem is improved for each number of agents p as well as for when the number of agents p increases. Once the optimization solution has been determined, a signal may be generated based upon the optimization solution and transmitted to the agents.
The exemplary optimization algorithm provides a manner in which global and local planning of the trajectory of multiple distinct agents may be determined. The optimization solution to be determined is based on solving a non-convex optimization problem using a modified ADMM algorithm in which a three-weight message-passing algorithm is incorporated. The optimization algorithm may flexibly account for obstacles and different cost functions. For agents in a plane or space, explicit expressions that account for static obstacles, moving obstacles, and dynamic constraints on the velocity and energy may be derived such as in the manner described above.
FIG. 6 shows a method 600 of generating an optimization solution according to an exemplary embodiment. In step 605, the device 100 receives the optimization problem. For example, via the I/O device 120, a user may provide the optimization problem. Within the optimization problem are a plurality of cost functions and a plurality of variables or constraints. Thus, in step 610 and 615, the device 100 (via the variable engine 125 and the cost engine 130) determine the cost functions and the variables of the optimization problem. It should be noted that the user may simply provide the cost functions and the variables. Thus, the device 100 may receive the cost functions and the variables.
It should be noted that the device 100 may utilize a factor graph such as the graphical models 200, 500. Thus, the method 600 may include steps that generate these factor graphs. For example, the factor graph engine 135 may generate the graphical models via the constraint engine 140 providing data as to the relationships between the cost functions and the variables. Specifically, the relationships may be those shown in the graphical models 200, 500 for the nodes on the left and the nodes on the right.
In step 620, the device 100 generates messages to be passed from the cost functions to the variables (e.g., the messages to be passed from the left nodes to the right nodes on the factor graph). In step 625, the device 100 also determine a respective weight to be applied to each message. As discussed above, the exemplary optimization algorithm utilizes a three-weight message-passing algorithm. Accordingly, each message that is generated also includes a corresponding weight of zero, infinite, and standard. In step 630, the message and corresponding weight are transmitted from the cost functions to the variables.
In step 635, based upon the messages and weights that were received from the cost functions to the variables, the device 100 generates further messages to be passed from the variables to the cost functions. In step 640, the device 100 determines a respective weight to be applied to each of these further messages. In step 645, the further message and corresponding weight are transmitted from the variables to the cost functions.
In step 650, the device determines whether there is a consensus between the messages being transmitted between the cost functions and the variables. If a consensus has not been reached, the method 600 returns to step 620 where steps 620-645 are repeated. This may repeat as many times as required for the consensus to be reached. Once the consensus is reached, the method 600 continues to step 655. In step 655, the optimization solution is generated as a function of the messages transmitted that were determined to form the consensus. That is, the optimization solution relates to a minimum being reached for each of the cost functions while still satisfying all the requirements of the variables.
It should be noted that the method 600 may include further steps. For example, the method 600 may include a criterion in which the consensus is only formed when a predetermined number of consecutive messages that are passed between the cost functions and the variables holds true. Thus, after step 650, the method 600 may include a further determination such that the satisfaction of the predetermined number is determined. In another example, the method 600 may include steps after step 625 and step 640 in which a determination is made whether nodes or input/outputs to the nodes may be ignored or removed from further computation. For example, when a zero weight or infinite weight is determined for a particular message, the graphical model may be modified to incorporate this knowledge to further improve efficiency in determining the optimization solution.
The exemplary embodiments provide a method and device for generating an optimization solution for an optimization problem. Specifically, the optimization solution is a modified version of the ADMM algorithm in which the messages in the message-passing algorithm includes various weights. The weights may indicate a “certainty,” a “no opinion,” or otherwise a “standard.” Through incorporating the weights, the optimization solution may be determined in a more efficient manner. The optimization solution may also be modified through incorporation of knowledge that is determined during the course of the message-passing. Whether a reasoner hierarchy or graph dynamics are used, the knowledge may further improve the efficiency at which the optimization solution is determined.
Those skilled in the art will understand that the above-described exemplary embodiments may be implemented in any suitable software or hardware configuration or combination thereof. An exemplary hardware platform for implementing the exemplary embodiments may include, for example, an Intel x86 based platform with compatible operating system, a Mac platform and MAC OS, etc. In a further example, the exemplary embodiments of the above described method may be embodied as a program containing lines of code stored on a non-transitory computer readable storage medium that, when compiled, may be executed on a processor or microprocessor.
It will be apparent to those skilled in the art that various modifications may be made in the present invention, without departing from the spirit or scope of the invention. Thus, it is intended that the present invention cover the modifications and variations of this invention provided they come within the scope of the appended claims and their equivalents.

## Claims

1. A method, comprising:
(a) receiving an optimization problem, the optimization problem including a plurality of cost functions and a plurality of variables, the cost functions representing possible costs for values of the variables in the optimization problem, each of the cost functions having a predetermined relationship with select ones of the variables,
wherein the optimization problem comprises a trajectory planning problem;

(b) generating a first message for each of the cost functions for each corresponding variable based upon the respective predetermined relationship, the first message indicating a first belief that the corresponding variable has a first value when the optimization problem is solved, the first message having a respective first weight indicating a certainty of the first message;
(c) generating a second message for each of the variables for each corresponding cost function based upon the respective predetermined relationship, the second message indicating a second belief that the corresponding variable has a second value when the optimization problem is solved, the second message having a respective second weight indicating a certainty of the second message, wherein each of the first and second weights is one of a zero weight, an infinite weight, and a standard weight;
(d) generating a disagreement variable for each corresponding pair of variables and cost functions measuring a disagreement value between the first and second beliefs;
(e) repeating steps (b), (c), and (d) until a consensus is formed between the first and second messages based on the cost functions and the variables, a subsequent first message being modified based upon the second message, its corresponding second weight, and the corresponding disagreement variable, and a subsequent second message being modified based upon the first message, its corresponding first weight, and the corresponding disagreement variable; and
(f) determining an optimization solution based upon the consensus,
wherein the optimization solution comprises a solution to the trajectory planning problem.

2. The method of claim 1, wherein the zero weight indicates that the corresponding message is not to be used in determining the optimization solution.

3. The method of claim 1, wherein the infinite weight indicates that the corresponding message specifies an unchangeable value for a variable in determining the optimization solution.

4. The method of claim 1, wherein the standard weight is given a predetermined value between zero and infinity.

5. The method of claim 1, wherein the first messages for the cost functions are computed in parallel to determine the optimization solution.

6. The method of claim 1, further comprising:
determining knowledge data based upon the first and second messages and the first and second weights, respectively;
incorporating the knowledge data in the subsequent first and second messages,
wherein the knowledge data indicates whether at least one of at least one of the cost functions and at least one of the variables are to be removed in determining the optimization solution.

7. The method of claim 6, wherein the knowledge data is at least one of a reasoner hierarchy between local and global aspects of the optimization problem and a graph dynamic between the predetermined relationship of the cost functions and the variables.

8. The method of claim 1, wherein the optimization problem is one of convex and non-convex.

9. The method of claim 1, wherein the optimization problem is a multi-agent trajectory planning problem.

10. A device, comprising:
a processor coupled to a memory, wherein the processor is programmed to determine an optimization solution to an optimization problem by:
(a) receiving the optimization problem, the optimization problem includes a plurality of cost functions and a plurality of variables, the cost functions representing possible costs for values of the variables in the optimization problem, each of the cost functions having a predetermined relationship with select ones of the variables,
wherein the optimization problem comprises a trajectory planning problem;

(b) generating a first message for each of the cost functions for each corresponding variable based upon the respective predetermined relationship, the first message indicating a first belief that the corresponding variable has a first value when the optimization problem is solved, the first message having a respective first weight indicating a certainty of the first message;
(c) generating a second message for each of the variables for each corresponding cost function based upon the respective predetermined relationship, the second message indicating a second belief that the corresponding variable has a second value when the optimization problem is solved, the second message having a respective second weight indicating a certainty of the second message, wherein each of the first and second weights is one of a zero weight, an infinite weight, and a standard weight;
(d) generating a disagreement variable for each corresponding pair of variables and cost functions measuring a disagreement value between the first and second beliefs;
(e) repeating steps (b), (c), and (d) until a consensus is formed between the first and second messages based on the cost functions and the variables, a subsequent first message being modified based upon the second message, its corresponding second weight, and the corresponding disagreement variable, and a subsequent second message being modified based upon the first message, its corresponding first weight, and the corresponding disagreement variable; and
(f) determining the optimization solution based upon the consensus,
wherein the optimization solution comprises a solution to the trajectory planning problem.

11. The method of claim 10, wherein the zero weight indicates that the corresponding message is not to be used in determining the optimization solution.

12. The method of claim 10, wherein the infinite weight indicates that the corresponding message specifies an unchangeable value for a variable in determining the optimization solution.

13. The method of claim 10, wherein the standard weight is given a predetermined value between zero and infinity.

14. The device of claim 10, wherein the first messages for the cost functions are computed in parallel to determine the optimization solution.

15. The device of claim 10, wherein the processor further determines knowledge data based upon the first and second messages and the first and second weights, respectively, and incorporates the knowledge data in the subsequent first and second messages, wherein the knowledge data indicates whether at least one of at least one of the cost functions and at least one of the variables are to be removed in determining the optimization solution.

16. The device of claim 15, wherein the knowledge data is at least one of a reasoner hierarchy between local and global aspects of the optimization problem and a graph dynamic between the predetermined relationship of the cost functions and the variables.

17. The device of claim 10, wherein the optimization problem is one of convex and non-convex.

18. A non-transitory computer readable storage medium with an executable program stored thereon, wherein the program instructs a microprocessor to perform operations comprising:
(a) receiving an optimization problem, the optimization problem includes a plurality of cost functions and a plurality of variables, the cost functions representing possible costs for values of the variables in the optimization problem, each of the cost functions having a predetermined relationship with select ones of the variables,
wherein the optimization problem comprises a trajectory planning problem;

(b) generating a first message for each of the cost functions for each corresponding variable based upon the respective predetermined relationship, the first message indicating a first belief that the corresponding variables has a first value when the optimization problem is solved, the first message having a respective first weight indicating a certainty of the first message;
(c) generating a second message for each of the variables for each corresponding cost function based upon the respective predetermined relationship, the second message indicating a second belief that the corresponding variable has a second value when the optimization is solved, the second message having a respective second weight indicating a certainty of the second message, wherein each of the first and second weights is one of a zero weight, an infinite weight, and a standard weight;
(d) generating a disagreement variable for each corresponding pair of variables and cost functions measuring a disagreement value between the first and second beliefs;
(e) repeating steps (b), (c), and (d) until a consensus is formed between the first and second messages based on the cost functions and the variables, a subsequent first message being modified based upon the second message, its corresponding second weight, and the corresponding disagreement variable, and a subsequent second message being modified based upon the first message, its corresponding first weight, and the corresponding disagreement variable; and
(f) determining an optimization solution based upon the consensus,
wherein the optimization solution comprises a solution to the trajectory planning problem.