Skip to main content

Divisor must be a constant

Answered

Comments

7 comments

  • Eli Towle
    Gurobi Staff Gurobi Staff

    As the error suggests, Gurobi doesn't support dividing by variables. For more information, see What types of models can Gurobi solve?.

    You could model this problem by introducing a few auxiliary variables and bilinear constraints. However, your problem looks very similar to problem of maximizing the Sharpe ratio (though you use variance instead of standard deviation), for which there is a nice convex reformulation. Let \( \mu \in \mathbb{R}^n \) be the vector of expected returns, and let \( Q \in \mathbb{S}_{+}^{n \times n} \) be the (symmetric, positive semidefinite) covariance matrix. I'll assume you are trying to maximize the Sharpe ratio \( \mu^\top x / \sqrt{x^\top Q x} \); I'm not aware of a common risk-adjusted return measure that uses variance instead of standard deviation.

    Your problem, which we refer to as \( {\bf(1)} \), can be written as

    $$\begin{alignat*}{2} \max_x\ && \frac{\mu^\top x}{\sqrt{x^\top Q x}} & \\ \textrm{s.t.}\ && \sum_{i=1}^n x_i ={} &1 \\ && x \geq {} &0.\end{alignat*}$$

    Consider the following problem, which we refer to as \( {\bf(2)} \):

    $$\begin{alignat*}{2} \max_y\ && \frac{1}{\sqrt{y^\top Q y}} & \\ \textrm{s.t.}\ && \mu^\top y = {} &1 \\ && y \geq {} &0.\end{alignat*}$$

    We'll show that \( {\bf(1)} \) and \( {\bf(2)} \) are "equivalent", in the sense that given a solution to either problem, we can construct a solution to the other of equal or better objective value. We assume there exists some \( \hat{x} \) for which \( \mu^\top \hat{x} > 0 \), which should hold in any reasonable portfolio optimization problem. It follows from this assumption that the optimal solution to \( {\bf(1)} \) has a positive expected return.

    First, let \( \bar{x} \) be a solution to \( {\bf(1)} \). We construct a solution to \( {\bf(2)} \) of equal objective value. For all \( i = 1, \ldots, n \), define \( \bar{y}_i := \bar{x}_i / \mu^\top \bar{x} \). It holds that \( \bar{y} \geq 0 \), because \( \bar{x} \) is nonnegative and \( \mu^\top \bar{x} > 0 \) by our assumption. Additionally, \( \mu^\top \bar{y} = \mu^\top \bar{x} / \mu^\top \bar{x} = 1 \). Thus, \( \bar{y} \) satisfies the constraints of \( {\bf(2)} \). Finally, the objective of \( {\bf(2)} \) with respect to the solution \( \bar{y} \) is \( 1 / \sqrt{\bar{y}^\top Q \bar{y}} = 1 / \sqrt{ \bar{x}^\top Q \bar{x} / (\mu^\top \bar{x})^2} = \mu^\top \bar{x} / \sqrt{\bar{x}^\top Q \bar{x}} \).

    Next, let \( \bar{y} \) be a solution to \( {\bf(2)} \). We construct a solution to \( {\bf(1)} \) of equal objective value. Let \( \bar{x}_i := \bar{y}_i / \sum_{j=1}^n \bar{y}_j \) for \( i = 1, \ldots, n \). Using the same ideas as above, it's easy to verify that \( \bar{x} \) satisfies the constraints of \( {\bf(1)} \) and has objective value equal to \( 1 / \sqrt{\bar{y}^\top Q \bar{y}} \).

    So, instead of solving your original problem \( {\bf(1)} \), we can just solve \( {\bf(2)} \) and map the optimal solution \( \bar{y} \) back to the original problem using the translation \( \bar{x}_i := \bar{y}_i / \sum_{j=1}^n \bar{y}_j \) for \( i = 1, \ldots, n \). However, because \( Q \) is positive semidefinite, the solutions of \( {\bf(2)} \) are the same as the solutions to

    $$\begin{alignat*}{2} \min_y\ && y^\top Q y \enspace & \\ \textrm{s.t.}\ && \mu^\top y = {} &1 \\ && y \geq {} &0.\end{alignat*}$$

    This is a convex QP that Gurobi can solve directly. A solution \( \bar{y} \) to the above problem can be mapped to the original problem using the same transformation \( \bar{x}_i := \bar{y}_i / \sum_{j=1}^n \bar{y}_j \) for \( i = 1, \ldots, n \).

    If you have other side constraints to \( {\bf(1)} \), see Erwin Kalvelagen's blog post on maximizing the Sharpe ratio. Basically, additional constraints \( Ax \leq b \) in \( {\bf(1)} \) translate to constraints \( Ay \leq b \kappa \) in \( {\bf(2)} \), where \( \kappa = \sum_{j=1}^n y_j \).

    2
  • Marek K.
    Gurobi-versary
    First Question
    First Comment

    My goodness Eli, this was a fantastic explanation. Thank you for taking the time to work through this and write this detailed of a reply. I am happy to say that I implemented this into my code without issue and got the exact result I was expecting. And yes, this is exactly the problem of maximizing the Sharpe ratio, but I had left out the square root due to other issues in an attempt to solve it in piecemeal fashion.

    Thanks again!!

    Marek

    0
  • Moritz Spatz

    Hi together, I have a similar problem, I understand Elis answer but am unable to implement the solution. How do I implement the described approach?

    I have:

    for tk in TK:
        model.addConstr(R[tk]*A[tk]==E[tk]) 

    R, A and E are all variables in the model.

     

    I would be glad if somebody could help here.

    Kind regards

    Moritz

    0
  • Eli Towle
    Gurobi Staff Gurobi Staff

    Which part of the implementation do you need help with? Your code looks okay, depending on how you defined \( \texttt{TK} \), \( \texttt{R} \), \( \texttt{A} \), and \( \texttt{E} \). For example, the following code runs without error:

    import gurobipy as gp

    TK = [1, 2, 3]

    model = gp.Model()

    R = model.addVars(TK, name='R')
    A = model.addVars(TK, name='A')
    E = model.addVars(TK, name='E')

    for tk in TK:
        model.addConstr(R[tk] * A[tk] == E[tk])

    If you encounter an error, it would be helpful to see a minimal, self-contained code snippet that reproduces the error.

    0
  • Moritz Spatz

    Hi again,

    below I send the snippet. It took a while to reduce the code to the problem, sorry for coming back late. I marked the line which makes trouble with a comment. Its about R[t]. Can you please support here? I am not aware of a a reasonable solution.

    Thanks in advance and kind regards.

    Moritz

     

    from gurobipy import *

    # =============================================================================
    model = Model("Smart Factory Strategie")
    model.modelSense = GRB.MAXIMIZE

    # =============================================================================
    A=['1','2','3']
    D={'1':1.5,'2':1.8,'3':1.5}
    KWS={'1':[-10842,-1355,-1355,-1355,-1355,-1355,-1355,-1355,-1355,-1355,-1355,-1355,0,0,0,0],'2':[-5421,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],'3':[-10842,-1355,-1355,-1355,-1355,-1355,-1355,-1355,-1355,-1355,-1355,-1355,0,0,0,0]}
    KWE={'1':[35090,35090,35090,35090,35090,35090,35090,35090,35090,35090,35090,35090,0,0,0,0],'2':[22562,22562,22562,22562,22562,22562,22562,22562,22562,22562,22562,22562,0,0,0,0],'3':[18492,18492,18492,18492,18492,18492,18492,18492,18492,18492,18492,18492,0,0,0,0]}
    Tmax=   10
    T=      [i for i in range(0,Tmax+1)] 


    # =============================================================================

    Y={}
    for a in A:
        for t in T:
            Y[a,t]=model.addVar(name="Y_%s_%s" % (a,t), vtype=GRB.BINARY)         

    Au={}
    for t in T:
        Au[t]=model.addVar(name="Au_%s" % (t), vtype=GRB.CONTINUOUS, lb=-1*10**30)
            
    BW={}
    for t in T:
        BW[t]=model.addVar(name="BW_%s" % (t), vtype =GRB.CONTINUOUS, lb=-1*10**30)
        
    BWX= model.addVar(vtype =GRB.CONTINUOUS , obj = 1 , name='BWX') 
                
    Er={}
    for t in T:
        Er[t]=model.addVar(name="Er_%s" % (t), vtype=GRB.CONTINUOUS, lb=-1*10**30)  
        
    R={}
    for t in T:
        R[t]=model.addVar(name="R_%s" % (t), vtype =GRB.CONTINUOUS, lb=-1*10**30)


    # =============================================================================

        
    model.addConstr(BWX==BW[Tmax])

    for t in T:
        model.addConstr(Au[t] == quicksum(KWS[a][t] * Y[a,t] for a in A))
        model.addConstr(Er[t] == quicksum(KWE[a][t] * Y[a,t] for a in A))
    for t in T:
        model.addConstr(BW[t]==Er[t]+Au[t]) 
        
    for t in T:
        model.addConstr(R[t]*Au[t]==Er[t]) #This line of code makes trouble

    # =============================================================================
    model.update()
    model.optimize()
    # =============================================================================

    0
  • Eli Towle
    Gurobi Staff Gurobi Staff

    It looks like you just need to set the NonConvex parameter to 2 to direct Gurobi to solve the model (which is made non-convex by the quadratic equality constraints):

    model.Params.NonConvex = 2

    Note that solving non-convex quadratic problems is only supported in Gurobi 9.0 onward.

    0
  • Moritz Spatz

    Thanks a lot. It works

    0

Please sign in to leave a comment.