有纪律的参数化编程

注意：DPP需要CVXPY版本1.1.0或更高。

Parameters 是 常量的符号表示。使用参数可以让你修改常量的值而无需重新构建整个问题。当你的参数化问题根据Disciplined Parametrized Programming (DPP)构建时，对于不同的参数值重复求解可以比重复求解一个新问题快得多。

如果您打算多次解决DCP或DGP问题，针对不同的数值数据，或者如果您想通过DCP或DGP问题的解决方案图进行微分，您应该阅读本教程。

什么是DPP？

DPP 是一套规则集，用于生成参数化的 DCP 或 DGP 兼容问题，CVXPY 可以非常快速地重新规范化这些问题。第一次解决符合 DPP 的问题时，CVXPY 会编译它并缓存从参数到问题数据的映射。因此，后续重写 DPP 问题可以显著加快。CVXPY 允许您解决不符合 DPP 的参数化问题，但在这样做时您不会看到速度的提升。

DPP规则集

DPP 对参数如何进入 DCP 和 DGP 问题中的表达式施加了轻微的限制。首先，我们描述了 DCP 问题的 DPP 规则集。然后，我们描述了 DGP 问题的 DPP 规则集。

DCP问题。 在DPP中，如果一个表达式不涉及变量并且在其参数中是仿射的，则称其为参数仿射的；如果它没有参数，则称其为无参数的。DPP对DCP引入了两个限制：

  1. 在DPP下，所有参数都被分类为仿射，就像变量一样。

  2. 在DPP下，当至少一个表达式是常数时，或者当一个表达式是参数仿射而另一个是无参数时，两个表达式的乘积是仿射的。

如果一个表达式在遵循这两个限制的情况下是DCP兼容的，那么它就是DPP兼容的。你可以通过调用is_dcp方法并设置关键字参数dpp=True来检查一个表达式或问题是否是DPP兼容的（默认情况下，这个关键字参数是False）。例如，

import cvxpy as cp


m, n = 3, 2
x = cp.Variable((n, 1))
F = cp.Parameter((m, n))
G = cp.Parameter((m, n))
g = cp.Parameter((m, 1))
gamma = cp.Parameter(nonneg=True)

objective = cp.norm((F + G) @ x - g) + gamma * cp.norm(x)
print(objective.is_dcp(dpp=True))

打印 True。我们可以通过DPP分析来理解为什么 objective 符合DPP。乘积 (F + G) @ x 在DPP下是仿射的，因为 F + G 是参数仿射的，而 x 是无参数的。差值 (F + G) @ x - g 是仿射的，因为加法原子是仿射的，且 (F + G) @ x 和 - g 都是仿射的。同样地，gamma * cp.norm(x) 在DPP下是仿射的，因为 gamma 是参数仿射的，而 cp.norm(x) 是无参数的。最终的目标在DPP下是仿射的，因为加法是仿射的。

一些表达式符合DCP但不符合DPP。例如，DPP禁止对两个参数化表达式进行乘积操作：

import cvxpy as cp


x = cp.Variable()
gamma = cp.Parameter(nonneg=True)
problem = cp.Problem(cp.Minimize(gamma * gamma * x), [x >= 1])
print("Is DPP? ", problem.is_dcp(dpp=True))
print("Is DCP? ", problem.is_dcp(dpp=False))

这段代码片段打印

Is DPP? False
Is DCP? True

正如可以将非DCP问题重写为符合DCP的方式一样，也可以将非DPP问题重新表达为符合DPP的方式。例如，上述问题可以等价地写为

import cvxpy as cp


x = cp.Variable()
y = cp.Variable()
gamma = cp.Parameter(nonneg=True)
problem = cp.Problem(cp.Minimize(gamma * y), [y == gamma * x])
print("Is DPP? ", problem.is_dcp(dpp=True))
print("Is DCP? ", problem.is_dcp(dpp=False))

这段代码打印

Is DPP? True
Is DCP? True

在其他情况下，您可以通过在DSL之外执行参数的非DPP转换来表示它们，例如在NumPy中。例如，如果P是一个参数，x是一个变量，cp.quad_form(x, P)不是DPP。您可以像这样表示一个参数化的二次形式：

import cvxpy as cp
import numpy as np
import scipy.linalg


n = 4
L = np.random.randn(n, n)
P = L.T @ L
P_sqrt = cp.Parameter((n, n))
x = cp.Variable((n, 1))
quad_form = cp.sum_squares(P_sqrt @ x)
P_sqrt.value = scipy.linalg.sqrtm(P)
assert quad_form.is_dcp(dpp=True)

再举一个例子，当p是参数时，商expr / p不符合DPP规范，但这可以重写为expr * p_tilde，其中p_tilde是表示1/p的参数。

DGP问题。 正如DGP是DCP的对数-对数类比，DGP的DPP是DCP的DPP的对数-对数类比。DPP为DGP引入了两个限制：

  1. 在DPP下，所有正参数都被分类为对数-对数-仿射，就像正变量一样。

  2. 在DPP下，幂原子 x**p（以 x 为底，p 为指数）是log-log仿射的，只要 x 和 p 不同时被参数化。

请注意，对于幂运算，指数 p 必须是一个数值常量或参数；尝试构造一个指数为复合表达式的幂原子，例如 x**(p + p)，其中 p 是一个参数，将会导致 ValueError。

如果一个参数在DGP问题中作为指数出现，它可以有任何符号。如果参数在DGP问题的其他地方出现，它必须是正的，即必须使用cp.Parameter(pos=True)来构造。

你可以通过调用is_dgp方法并设置关键字参数dpp=True（默认情况下，此关键字参数为False）来检查表达式或问题是否符合DPP。例如，

import cvxpy as cp


x = cp.Variable(pos=True)
y = cp.Variable(pos=True)
a = cp.Parameter()
b = cp.Parameter()
c = cp.Parameter(pos=True)

monomial = c * x**a * y**b
print(monomial.is_dgp(dpp=True))

打印 True。表达式 x**a 和 y**b 是对数-对数仿射的，因为 x 和 y 不包含参数。参数 c 是对数-对数仿射的 因为它是正的，而单项式表达式是对数-对数仿射的，因为 对数-对数仿射表达式的乘积也是对数-对数仿射的。

一些表达式符合DGP但不符合DPP。例如，DPP禁止将参数化表达式提升为幂：

import cvxpy as cp


x = cp.Variable(pos=True)
a = cp.Parameter()

monomial = (x**a)**a
print("Is DPP? ", monomial.is_dgp(dpp=True))
print("Is DGP? ", monomial.is_dgp(dpp=False))

这段代码片段打印

Is DPP? False
Is DGP? True

你可以通过在CVXPY外部进行参数的转换来表示非DPP转换，例如在NumPy中。例如，你可以将上述程序重写为以下符合DPP的程序

import cvxpy as cp


a = 2.0
x = cp.Variable(pos=True)
b = cp.Parameter(value=a**2)

monomial = x**b

重复解决DPP问题

以下示例演示了参数如何加速重复解决符合DPP的DCP问题。

import cvxpy as cp
import numpy
import matplotlib.pyplot as plt
import time

n = 15
m = 10
numpy.random.seed(1)
A = numpy.random.randn(n, m)
b = numpy.random.randn(n)
# gamma must be nonnegative due to DCP rules.
gamma = cp.Parameter(nonneg=True)

x = cp.Variable(m)
error = cp.sum_squares(A @ x - b)
obj = cp.Minimize(error + gamma*cp.norm(x, 1))
problem = cp.Problem(obj)
assert problem.is_dcp(dpp=True)

gamma_vals = numpy.logspace(-4, 1)
times = []
new_problem_times = []
for val in gamma_vals:
    gamma.value = val
    start = time.time()
    problem.solve()
    end = time.time()
    times.append(end - start)
    new_problem = cp.Problem(obj)
    start = time.time()
    new_problem.solve()
    end = time.time()
    new_problem_times.append(end - start)

plt.rc('text', usetex=True)
plt.rc('font', family='serif')
plt.figure(figsize=(6, 6))
plt.plot(gamma_vals, times, label='Re-solving a DPP problem')
plt.plot(gamma_vals, new_problem_times, label='Solving a new problem')
plt.xlabel(r'$\gamma$', fontsize=16)
plt.ylabel(r'time (s)', fontsize=16)
plt.legend()
../../_images/re_solving_dpp.png

对于DGP问题，可以获得类似的速度提升。

敏感性分析和梯度

在版本1.1.0中新增：此功能需要CVXPY版本1.1.0或更高版本

优化问题可以被视为将参数映射到解的函数。这个解映射有时是可微分的。CVXPY 内置支持计算问题的最优变量值相对于参数的小扰动的导数（即，问题中出现的 Parameter 实例）。

问题类暴露了两个与计算导数相关的方法。 derivative 方法评估给定参数扰动时的导数。这使您能够计算问题解决方案在参数发生小变化时的变化，而无需重新解决问题。 backward 方法评估导数的伴随，计算解决方案相对于参数的梯度。当与自动微分软件结合使用时，这可能非常有用。

导数和反向传播方法仅在问题包含参数时才有意义。为了使问题可微分，它必须符合DPP。CVXPY可以计算任何符合DPP的DCP或DGP问题的导数。在不可微分的点，CVXPY计算一个启发式量。

示例。

作为第一个例子，我们解决一个具有解析解的简单问题，以说明backward和derivative函数的使用。在下面的代码块中，我们构建了一个具有标量变量x和标量参数p的问题。问题是最小化二次函数(x - 2*p)**2。

import cvxpy as cp

x = cp.Variable()
p = cp.Parameter()
quadratic = cp.square(x - 2 * p)
problem = cp.Problem(cp.Minimize(quadratic))

接下来，我们解决特定值 p == 3 的问题。请注意， 在解决问题时，我们向 solve 方法提供了关键字参数 requires_grad=True。

p.value = 3.
problem.solve(requires_grad=True)

解决了requires_grad=True的问题后，我们现在可以使用backward和derivative来对问题进行微分。首先，我们通过调用backward()方法来计算解相对于其参数的梯度。作为副作用，backward()方法会填充所有参数的gradient属性，该属性表示解相对于该参数的梯度。

problem.backward()
print("The gradient is {0:0.1f}.".format(p.gradient))

在这种情况下，问题有一个简单的解析解 2*p，因此梯度就是2。所以，正如预期的那样，上面的代码会打印

The gradient is 2.0.

接下来，我们使用derivative方法来查看p的微小变化如何影响解x。我们将通过设置p.delta = 1e-5来扰动p，并调用derivative方法，该方法将用一阶近似预测的x的变化（即dx/dp * p.delta）填充x的delta属性。

p.delta = 1e-5
problem.derivative()
print("x.delta is {0:2.1g}.".format(x.delta))

在这种情况下，解决方案是微不足道的，其导数只是 2*p，所以我们知道 x 中的增量应该是 2e-5。正如预期的那样，输出是

x.delta is 2e-05.

我们强调这个例子是简单的，因为它有一个简单的解析解，且导数也是简单的。backward() 和 forward() 方法之所以有用，是因为绝大多数凸优化问题没有解析解：在这些情况下，CVXPY 可以计算解及其导数，尽管手工推导这些解是不可能的。

注意。 在这个简单的例子中，变量 x 是一个标量，因此 backward 方法计算了 x 相对于 p 的梯度。 当有多个标量变量时，默认情况下，backward 计算最优变量值的总和相对于参数的梯度。

更一般地，backward 方法可用于计算标量值函数 f 相对于参数的最优变量的梯度。如果 x(p) 表示变量（可能是向量或矩阵）在参数 p 的特定值下的最优值，并且 f(x(p)) 是一个标量，那么可以使用 backward 来计算 f 相对于 p 的梯度。设 x* = x(p)，并假设 f 相对于 x* 的导数是 dx。要计算 f 相对于 p 的导数，在调用 problem.backward() 之前，只需设置 x.gradient = dx。

当与自动微分软件结合使用时，backward 方法可以非常强大。我们推荐使用软件包 CVXPY Layers，它为CVXPY问题提供了可微分的PyTorch和TensorFlow封装。

反向传播还是求导？ 当你需要解的（标量值函数）相对于参数的梯度时，应该使用 backward 方法。如果你只想进行敏感性分析，也就是说，如果你只关心当一个或多个参数改变时解会如何变化，你应该使用 derivative 方法。当存在多个变量时，使用求导方法计算敏感性比计算整个雅可比矩阵（可以通过多次调用反向传播方法，每次针对一个标准基向量）要高效得多。

下一步。 查看关于导数的入门笔记本。