DGP基础

本文档将介绍**标准几何规划**（DGP）的基本概念，DGP允许您在CVXPY中制定和解决*log-log凸规划*（LLCP）。

LLCP是在变量、目标函数和约束函数被其对数替代后变为凸性的问题，我们称之为对数对数变换。LLCP是几何规划的推广。

import cvxpy as cp

1. 对数对数凸性

与CVXPY中的每个表达式都具有曲率（常数、仿射、凸、凹或未知）一样，每个表达式也都有对数对数曲率。

如果函数 \(f : D \subseteq \mathbf{R}^{n}_{++} \to \mathbf{R}\) 是**对数对数凸**的，那么函数 \(F(u)=\log f(e^u)\) 是一个凸函数，其定义域为 \(\{u \in \mathbf{R}^n : e^u \in D\}`（其中 :math:\)mathbf{R}^{n}_{++}` 表示正实数的集合，对数和指数运算是逐元素进行的）；函数 \(F\) 被称为 \(f\) 的对数对数变换。如果 \(F\) 是凹的，则函数 \(f\) 是**对数对数凹**的；如果 \(F\) 是仿射的，则函数 \(f\) 是**对数对数仿射**的。

注意，如果一个函数具有对数对数凸性，则其定义域和值域只能包括正数。

可以通过 .log_log_curvature 属性访问``Expression``的对数对数曲率。要使表达式具有已知的对数对数曲率，它引用的所有``Constant``、``Variable``和``Parameter``必须逐元素为正。.. code:: python

f = cp.sum(cp.log1p(x))

f(x) = log(1 + x_1) + log(1 + x_2) + ldots + log(1 + x_n)

The following functions are special cases of log-log affine functions:

  1. Exponential function

    f = cp.exp(x)

    f(x) = exp(x_1) + exp(x_2) + ldots + exp(x_n)

  2. Power function

    f = cp.power(x, p)

    f(x) = x_1^p_1 cdot x_2^p_2 ldots x_n^p_n

where \(p\) is a positive constant parameter.

\[f(x) = cx_1^{a_1}x_2^{a_2} \ldots x_n^{a_n},\]

其中 \(c > 0\)，而 \(a_i\) 是实数。在几何规划的背景下，这样的函数被称为单项式。

x = cp.Variable(shape=(3,), pos=True, name="x")
c = 2.0
a = [0.5, 2.0, 1.8]

单项式 = c * x[0] ** a[0] * x[1] ** a[1] * x[2] ** a[2]
# 单项式不是凸函数。
assert not 单项式.is_convex()

# 然而它们是对数-对数仿射的。
print(单项式, ":", 单项式.log_log_curvature)
assert 单项式.is_log_log_affine()

2.0 * power(x[0], 1/2) * power(x[1], 2) * power(x[2], 9/5) : 对数-对数仿射

多个单项式函数的和是对数-对数凸的；在几何规划的背景下，这样的函数被称为多项式。有一些函数不是多项式但仍然是对数-对数凸的。

x = cp.Variable(pos=True, name="x")
y = cp.Variable(pos=True, name="y")

常数 = cp.Constant(2.0)
单项式 = 常数 * x * y
多项式 = 单项式 + (x ** 1.5) * (y ** -1)
倒数 = 多项式 ** -1
未知 = 倒数 + 多项式

print(常数, ":", 常数.log_log_curvature)
print(单项式, ":", 单项式.log_log_curvature)
print(多项式, ":", 多项式.log_log_curvature)
print(倒数, ":", 倒数.log_log_curvature)
print(未知, ":", 未知.log_log_curvature).. parsed-literal::

2.0: 对数-对数常数
2.0 * x * y: 对数-对数仿射
2.0 * x * y + power(x, 3/2) * power(y, -1): 对数-对数凸
power(2.0 * x * y + power(x, 3/2) * power(y, -1), -1): 对数-对数凹
power(2.0 * x * y + power(x, 3/2) * power(y, -1), -1) + 2.0 * x * y + power(x, 3/2) * power(y, -1): 未知

2. 对数-对数曲率规则集

CVXP 具有具有已知对数-对数曲率和单调性的原子函数库。它使用这些信息对每个“Expression”（即原子函数的组合）进行标记，并分配对数-对数曲率。具体而言，

一个函数 \(f(expr_1,expr_2,...,expr_n)\) 是对数-对数凸的， 如果 \(f\) 是一个对数-对数凸函数，并且对于每个 \(expr_i\) 至少满足下列条件之一：

  • \(f\) 关于第 \(i\) 个参数递增，且 \(expr_i\) 是对数-对数凸的。

  • \(f\) 关于第 \(i\) 个参数递减，且 \(expr_i\) 是对数-对数凹的。

  • \(expr_i\) 是对数-对数仿射的。

一个函数 \(f(expr_1,expr_2,...,expr_n)\) 是对数-对数凹的， 如果 \(f\) 是一个对数-对数凹函数，并且对于每个 \(expr_i\) 至少满足下列条件之一：

  • \(f\) 关于第 \(i\) 个参数递增，且 \(expr_i\) 是对数-对数凹的。

  • \(f\) 关于第 \(i\) 个参数递减，且 \(expr_i\) 是对数-对数凸的。

  • \(expr_i\) 是对数-对数仿射的。

一个函数 \(f(expr_1,expr_2,...,expr_n)\) 是对数-对数仿射的， 如果 \(f\) 是一个对数-对数仿射函数，并且每个 \(expr_i\) 都是对数-对数仿射的。

如果以上三个规则都不适用，表达式 \(f(expr_1,expr_2,...,expr_n)\) 将被标记为具有未知曲率。

如果一个Expression满足组合规则，我们通常说该Expression“是DGP”。你可以通过调用方法``is_dgp()``来检查一个Expression是否是DGP。

x = cp.Variable(pos=True, name="x")
y = cp.Variable(pos=True, name="y")

单项式 = 2.0 * x * y
多项式 = 单项式 + (x ** 1.5) * (y ** -1)

输出(单项式, "是 DGP 吗?", 单项式.is_dgp())
输出(多项式, "是 DGP 吗?", 多项式.is_dgp())

2.0 * x * y is dgp? True
2.0 * x * y + power(x, 3/2) * power(y, -1) is dgp? True

3. DGP 问题

LLCP 是以下形式的优化问题

\[\begin{split}\begin{equation} \begin{array}{ll} \mbox{最小化} & f_0(x) \\ \mbox{满足} & f_i(x) \leq \tilde{f_i}, \quad i=1, \ldots, m\\ & g_i(x) = \tilde{g_i}, \quad i=1, \ldots, p, \end{array} \end{equation}\end{split}\]

其中函数 \(f_i\) 是对数对数凸函数，\(\tilde{f_i}\) 是对数对数凹函数，函数 \(g_i\) 和 \(\tilde{g_i}\) 是对数对数仿射函数。如果一个优化问题的约束满足上述形式， 并且目标是最大化一个对数对数凹函数，则它也是一个 LLCP 问题。

如果所有的函数都是 DGP，则问题是 DGP 问题。可以通过调用 CVXPY Problem 的 .is_dgp() 方法来检查一个问题是否是 DGP 问题。

x = cp.Variable(pos=True, name="x")
y = cp.Variable(pos=True, name="y")
z = cp.Variable(pos=True, name="z")

objective_fn = x * y * z
constraints = [
  4 * x * y * z + 2 * x * z <= 10, x <= 2*y, y <= 2*x, z >= 1]
assert objective_fn.is_log_log_concave()
assert all(constraint.is_dgp() for constraint in constraints)
problem = cp.Problem(cp.Maximize(objective_fn), constraints)

print(problem)
print("Is this problem DGP?", problem.is_dgp())

最大化 x * y * z
约束条件 4.0 * x * y * z + 2.0 * x * z <= 10.0
           x <= 2.0 * y
           y <= 2.0 * x
           1.0 <= z
这个问题是否是DGP？ True

解决DGP问题

您可以通过调用``gp=True``的``solve``方法解决DGP“问题”。

problem.solve(gp=True)
print("最优值:", problem.value)
print(x, ":", x.value)
print(y, ":", y.value)
print(z, ":", z.value)
print("对偶值: ", list(c.dual_value for c in constraints))

最优值: 1.9999999926890524
x : 0.9999999989968756
y : 1.9999999529045318
z : 1.000000020895385
对偶值:  [1.1111111199586956, 1.94877846244994e-09, 0.1111111217156332, 0.11111112214962586]

如果您忘记提供``gp=True``，将会引发错误。.. code:: python

try:

problem.solve()

except cp.DCPError as e:

print(e)

问题不符合DCP规则。然而，问题确实符合DGP规则。考虑使用`gp=True`调用此函数。

4. 下一步

原子

CVXPY有一个大型的对数凸函数库，包括常见的函数，如 \(\exp\)， \(\log\)，和两个数之间的差。请查看我们网站上的教程以获取完整的原子函数列表：https://www.cvxpy.org/tutorial/dgp/index.html

参考

有关DGP的参考，请参阅以下论文： https://web.stanford.edu/~boyd/papers/dgp.html