# How to represent the objective function of the Weapon Target Assignment problem in CVXPY?

I am trying to use CVXPY to analyse a problem and the objective function for this problem involves calculating a product and a sum as per the problem description below (taken from a draft paper I am writing):

However, I cannot figure out how to specify such an objective function in CVXPY given that it uses a product to combine all of the probability values for a particular target. In order to replicate the objective function, I need to multiply all of the elements of each row in a matrix together (i.e. multiplication along an axis). I am hoping to achieve something like the function shown below, which doesn't work because of the use of the numpy function.

objective = cp.Maximize(cp.sum(np.prod(1-cp.exp(cp.multiply(x, cp.log(q))), axis=1)))


Is there any way that I can perform multiplication along an axis in CVXPY? Or, is there an alternative way to represent my objective function?

• Independence of destruction (kill)( across weapons is not realistic in most cases. But I believe you could introduce a reasonable shot to shot correlation model of one type or another, while the objective would still be a posynomial, and DGP could be applied as peer my answer. Mar 9, 2023 at 18:14
• I think you can use the idea docs.mosek.com/modeling-cookbook/… to reformulate the objective. Mar 9, 2023 at 18:17
• Thank you both for your help, your answers helped me to figure out how to resolve my issue. Mar 16, 2023 at 5:00

TLDR: This can be formulated and solved in CXVPY with Mosek as solver, as a Mixed-Integer (generalized) Geometric Programming problem, using CVXPY's Disciplined Geometric Programming (DGP) capability. It can also be formulated and solved in CVX using its gp mode.

Details: This satisfies all of CVXPY's Disciplined Geometric Programming (DGP) rules. The objective function is a sum of monomials, which is a posynomial, and is log-log convex. So this can be formulated as a (generalized) Mixed-Integer Geometric Programming problem. Mosek is the best (only?) solver to use, because it has native exponential cone capability, which is invoked by CVXPY to solve the generalized Geometric Programming problem, and can solve such problems subject to integer constraints on some or all of the variables.

Now $$\begin{array}{rcl} t_j & \geq & \prod_{i=1}^m (1-p_{ij})^{x_{ij}} \\ & = & exp(\ln(\prod_{i=1}^m (1-p_{ij})^{x_{ij}})) \\ & = & exp(\sum_{i=1}^m x_{ij} \ln(1-p_{ij})) \\ \end{array}$$ This is the same as saying $$(t_j;1;\sum_{i=1}^m x_{ij} \ln(1-p_{ij}))$$ has to belong to the primal exponential cone.

So you can easily rewrite your problem as a conic optimization problem involving the exponential cone.

Mosek for sure can solve those problems even if some variables are constrained to be integers. I guess Cvxpy+Mosek should be handle that formulation. Not sure how hard the problem will be given the integer constraints.

• CVXPY, CVX, and CVXR will automatically do this reformulation for you if DGP rules (mode) are followed, per my answer. Mar 9, 2023 at 18:39
• True. However, since the trick is so simple I think might be valuable to see what it is. For instance, if you want to do this in a non cvx* context. Mar 10, 2023 at 5:45
• True. I gave my answer the way I gave it because the question was how to do it in CVXPY. So both our answers provide value. Mar 10, 2023 at 12:00

There is no need for DGP rules - the problem can be formulated directly as a mixed integer DCP, and you were pretty close to it. I know you said you figured it out, but I'm adding this for future visits to this question (since neither give a CVXPY implementation).

Somewhat in line with your code, I'm setting this up initially as:

$$\min \sum_{j=1}^n \prod_{i=1}^m q_{ij}^{x_{ij}}$$

s.t. $$\sum_{j=1}^n x_{ij} \leq W_i \quad \forall i$$

$$x_{ij} \in \mathbb{Z}_+$$

where $$q_{ij} = 1-p_{ij}$$

The objective function becomes $$\sum_{j=1}^n \prod_{i=1}^m e^{ln(q_{ij})x_{ij}}$$

and then $$\sum_{j=1}^n e^{\sum_{i=1}^m ln(q_{ij})x_{ij}}$$.

Since $$x$$ is the variable, we have the sum of exponentials of a linear function of $$x$$ -- it's a mixed integer DCP.

Here's a quick script for it in CVXPY:

import cvxpy as cp
import numpy as np

# Define the problem data
n = 5 # number of targets
m = 3 # number of weapon types
np.random.seed(1)
q = np.random.rand(n,m)*.9 + .1 # Survival probability
W = np.random.randint(1,10,m) # Number of weapons of each type

# Define the CVXPY problem.
x = cp.Variable((n,m), integer=True)
weighted_weapons = cp.multiply(x, np.log(q)) # (n,m)
survival_probs = cp.exp(cp.sum(weighted_weapons, axis=1)) # (n,)
objective = cp.Minimize(cp.sum(survival_probs))
cons = [cp.sum(x, axis=0) <= W, x >= 0]

# Solve
prob = cp.Problem(objective, cons)
prob.solve(verbose=True)
prob.status # Optimal
prob.value # 0.0785
np.set_printoptions(precision=3, suppress=True)
x.value
# array([[0., 0., 2.],
#        [0., 3., 0.],
#        [3., 0., 0.],
#        [1., 4., 0.],
#        [4., 0., 0.]])
$$$$
`