7.8 Library of basic functions

This section contains a library of small models of basic functions frequently appearing in optimization models. It is essentially an implementation of the mathematical models from the MOSEK Modeling Cookbook using Fusion API for .NET. These short code snippets can be seen as illustrative examples, can be copy-pasted to other code, and can even be directly called when assembling optimization models as we show in Sec. 7.8.6 (Model assembly example) (although this may be more suitable for prototyping; also note that additional variables and constraints will be introduced and there is no error checking).

7.8.1 Variable and constraint management

Variable duplication

\(x=y\)

Listing 7.16 Duplicate variables. Click here to download.
    // x = y
    public static void dup(Model M, Variable x, Variable y) {
      M.Constraint(Expr.Sub(x,y), Domain.EqualsTo(0.0));
    }

7.8.2 Linear operations

Absolute value

\(t\geq |x|\)

Listing 7.17 Absolute value. Click here to download.
    // t >= |x|, where t, x have the same shape
    public static void abs(Model M, Variable t, Variable x) {
      M.Constraint(Expr.Add(t,x), Domain.GreaterThan(0.0));
      M.Constraint(Expr.Sub(t,x), Domain.GreaterThan(0.0));
    }

1-norm

\(t\geq \sum_i |x_i|\)

Listing 7.18 \(1\)-norm. Click here to download.
    // t >= sum( |x_i| ), x is a vector Variable
    public static void norm1(Model M, Variable t, Variable x) {
      Variable u = M.Variable(x.GetShape(), Domain.Unbounded());
      abs(M, u, x);
      M.Constraint(Expr.Sub(t, Expr.Sum(u)), Domain.GreaterThan(0.0));
    }

7.8.3 Quadratic and power operations

Square

\(t\geq x^2\)

Listing 7.19 Square. Click here to download.
    // t >= x^2
    public static void sq(Model M, Variable t, Variable x) {
      M.Constraint(Expr.Hstack(0.5, t, x), Domain.InRotatedQCone());
    }

2-norm

\(t\geq \sqrt{\sum_i x_i^2}\)

Listing 7.20 \(2\)-norm. Click here to download.
    // t >= sqrt(x_1^2 + ... + x_n^2) where x is a vector
    public static void norm2(Model M, Variable t, Variable x) {
      M.Constraint(Expr.Vstack(t, x), Domain.InQCone());
    }

Powers

\(t\geq |x|^p\), \(p>1\)

Listing 7.21 Power. Click here to download.
    // t >= |x|^p (where p>1)
    public static void pow(Model M, Variable t, Variable x, double p) {
      M.Constraint(Expr.Hstack(t, 1, x), Domain.InPPowerCone(1.0/p));
    }

\(t\geq 1/x^p,\ x>0\), \(p>0\)

Listing 7.22 Power reciprocal. Click here to download.
    // t >= 1/|x|^p, x>0 (where p>0)
    public static void pow_inv(Model M, Variable t, Variable x, double p) {
      M.Constraint(Expr.Hstack(t, x, 1), Domain.InPPowerCone(1.0/(1.0+p)));
    }

p-norm

\(t\geq (\sum_i |x_i|^p)^{1/p}\), \(p>1\)

Listing 7.23 \(p\)-norm. Click here to download.
    // t >= \|x\|_p (where p>1), x is a vector Variable
    public static void pnorm(Model M, Variable t, Variable x, double p) {
      int n = (int) x.GetSize();
      Variable r = M.Variable(n);
      M.Constraint(Expr.Sub(t, Expr.Sum(r)), Domain.EqualsTo(0.0));
      M.Constraint(Expr.Hstack(Var.Repeat(t,n), r, x), Domain.InPPowerCone(1.0-1.0/p));
    }

Geometric mean

\(t\leq (x_1\cdot\cdots\cdot x_n)^{1/n}\), \(x_i>0\)

Listing 7.24 Geometric mean. Click here to download.
    // |t| <= (x_1...x_n)^(1/n), x_i>=0, x is a vector Variable of length >= 1
    public static void geo_mean(Model M, Variable t, Variable x) {
      int n = (int) x.GetSize();
      if (n==1) {
          abs(M, x, t);
      }
      else {
          Variable t2 = M.Variable();
          M.Constraint(Expr.Hstack(t2, x.Index(n-1), t), Domain.InPPowerCone(1.0-1.0/n));
          geo_mean(M, t2, x.Slice(0,n-1));
      }
    }

7.8.4 Exponentials and logarithms

log

\(t\leq \log{x},\ x>0\)

Listing 7.25 Logarithm. Click here to download.
    // t <= log(x), x>=0
    public static void log(Model M, Variable t, Variable x) {
      M.Constraint(Expr.Hstack(x, 1, t), Domain.InPExpCone());
    }

exp

\(t\geq e^{x}\)

Listing 7.26 Exponential. Click here to download.
    // t >= exp(x)
    public static void exp(Model M, Variable t, Variable x) {
      M.Constraint(Expr.Hstack(t, 1, x), Domain.InPExpCone());
    }

Entropy

\(t\geq x\log{x},\ x>0\)

Listing 7.27 Entropy. Click here to download.
    // t >= x * log(x), x>=0
    public static void ent(Model M, Variable t, Variable x) {
      M.Constraint(Expr.Hstack(1, x, Expr.Neg(t)), Domain.InPExpCone());
    }

Relative entropy

\(t\geq x\log{x/y},\ x,y>0\)

Listing 7.28 Relative entropy. Click here to download.
    // t >= x * log(x/y), x,y>=0
    public static void relent(Model M, Variable t, Variable x, Variable y) {
      M.Constraint(Expr.Hstack(y, x, Expr.Neg(t)), Domain.InPExpCone());
    }

Log-sum-exp

\(\log{\sum_i e^{x_i}}\leq t\)

Listing 7.29 Log-sum-exp. Click here to download.
    // log( sum_i(exp(x_i)) ) <= t, where x is a vector
    public static void logsumexp(Model M, Variable t, Variable x) {
      int n = (int) x.GetSize();
      Variable u = M.Variable(n);
      M.Constraint(Expr.Hstack(u, Expr.ConstTerm(n, 1.0), Expr.Sub(x, Var.Repeat(t, n))), Domain.InPExpCone());
      M.Constraint(Expr.Sum(u), Domain.LessThan(1.0));
    }

7.8.5 Integer Modeling

Semicontinuous variable

\(x\in\{0\}\cup[a,b]\), \(b>a>0\)

Listing 7.30 Semicontinuous variable. Click here to download.
    // x = 0 or a <= x <= b
    public static void semicontinuous(Model M, Variable x, double a, double b) {
      Variable u = M.Variable(x.GetShape(), Domain.Binary());
      M.Constraint(Expr.Sub(x, Expr.Mul(a, u)), Domain.GreaterThan(0.0));
      M.Constraint(Expr.Sub(x, Expr.Mul(b, u)), Domain.LessThan(0.0));
    }

Indicator variable

\(x\neq 0 \implies t=1\). We assume \(x\) is a priori normalized so \(|x_i|\leq 1\).

Listing 7.31 Indicator variable. Click here to download.
    // x!=0 implies t=1. Assumes that |x|<=1 in advance.
    public static void indicator(Model M, Variable t, Variable x) {
      M.Constraint(t, Domain.InRange(0,1));
      t.MakeInteger();
      abs(M, t, x);
    }

Logical OR

At least one of the conditions is true.

Listing 7.32 Logical OR. Click here to download.
    // x OR y, where x, y are binary
    public static void logic_or(Model M, Variable x, Variable y) {
      M.Constraint(Expr.Add(x, y), Domain.GreaterThan(1.0));
    }
    // x_1 OR ... OR x_n, where x is a binary vector
    public static void logic_or_vect(Model M, Variable x) {
      M.Constraint(Expr.Sum(x), Domain.GreaterThan(1.0));
    }

Logical NAND

At most one of the conditions is true (also known as SOS1).

Listing 7.33 Logical NAND. Click here to download.
    // at most one of x_1,...,x_n, where x is a binary vector (SOS1 Constraint)
    public static void logic_sos1(Model M, Variable x) {
      M.Constraint(Expr.Sum(x), Domain.LessThan(1.0));
    }
    // NOT(x AND y), where x, y are binary
    public static void logic_nand(Model M, Variable x, Variable y) {
      M.Constraint(Expr.Add(x, y), Domain.LessThan(1.0));
    }

Cardinality bound

At most \(k\) of the continuous variables are nonzero. We assume \(x\) is a priori normalized so \(|x_i|\leq 1\).

Listing 7.34 Cardinality bound. Click here to download.
    // At most k of entries in x are nonzero, assuming in advance |x_i|<=1.
    public static void card(Model M, Variable x, int k) {
      Variable t = M.Variable(x.GetShape(), Domain.Binary());
      abs(M, t, x);
      M.Constraint(Expr.Sum(t), Domain.LessThan(k));
    }

7.8.6 Model assembly example

We now demonstrate how to quickly build a simple optimization model for the problem

(7.19)\[\begin{split}\begin{array}{rl} \maximize & -\sqrt{x^2 + y^2} + \log{y} - x^{1.5}, \\ \st & x\geq y+3, \end{array}\end{split}\]

or equivalently

\[\begin{split}\begin{array}{rl} \maximize & -t_0+t_1-t_2, \\ \st & x\geq y+3, \\ & t_0\geq\sqrt{x^2+y^2}, \\ & t_1\leq\log{y}, \\ & t_2\geq x^{1.5}. \end{array}\end{split}\]
Listing 7.35 Modeling (7.19). Click here to download.
    public static void testExample() {
      Model M = new Model();
      Variable x = M.Variable();
      Variable y = M.Variable();
      Variable t = M.Variable(3);

      M.Constraint(Expr.Sub(x, y), Domain.GreaterThan(3.0));
      norm2(M, t.Index(0), Var.Vstack(x,y));
      log  (M, t.Index(1), y);
      pow  (M, t.Index(2), x, 1.5);

      M.Objective(ObjectiveSense.Maximize, Expr.Dot(t, new double[]{-1,1,-1}));