Continuous water-filling optimization problem

Question

Disclaimer: this question has been previously posted on Math StackExchange. I reposted it here since I did not receive any satisfactory answer there and a user suggested to re-post it here.

Let $x\in\mathbb{R}^n$ be an optimization variable and $\alpha\in\mathbb{R}^n$ be a given $n$-dimensional vector. The standard water-filling problem is formulated as

\begin{equation} \begin{array}{ll} \underset{x}{\operatorname{minimize}} & -\displaystyle\sum_{i=1}^{n} \log \left(\alpha_{i}+x_{i}\right) \\ \text { subject to } & x \succeq 0, \quad \mathbf{1}^{T} x=1, \end{array} \end{equation}

and it has a well-known solution (see Boyd & Vandenberghe page 245).

I was thinking about the case in which we have continuous communication channel slots. Intuitively this may be thought of as the case when the sizes of the communication channels approach zero. How is it possible to generalize this optimization problem to this case?

I believe that $\alpha$ and $x$, which from now on we will denote by $\alpha(x)$ and $f(x)$ respectively, would in this case be continuous real function and the problem would be:

\begin{equation} \begin{array}{ll} \underset{f \in \mathcal{F}}{\operatorname{minimize}} & -\displaystyle\int_{x} \log \left(\alpha(x)+f(x)\right)dx \\ \text { subject to } & f(x) \geq 0\; \forall x, \quad \int_x f(x)dx=1 \end{array} \end{equation}

where $\mathcal{F}$ is a given class of functions (e.g. Hölder continuous, Lipschitz, etc). We can also assume that the domain of the functions $f(x)$ and $\alpha(x)$ is compact, i.e. $x\in\mathcal{X}\subseteq \mathbb{R}$, with $\mathcal{X}$ a compact subspace of $\mathbb{R}$.

Am I following the right path?

Do you have any idea on how to solve these types of problems? It looks it is related to the calculus of variations, but I have never seen these types of problems and I have no idea how to solve them.

Answer 1

Using Calculus of Variations as an inspiration we have the lagrangian

$$ \mathcal{L}(f,\lambda) = \int_0^{x_{max}}\left(-\ln\left(\alpha(x)+f(x)\right)+\lambda\left(f(x)-\frac{1}{x_{max}}\right)\right)dx\ \ \ \text{s. t.}\ \ \ f(x)\ge 0 $$

The variation regarding $f$ gives

$$ -\frac{1}{\alpha(x)+f(x)}+\lambda = 0\Rightarrow f(x) = \frac{1}{\lambda}-\alpha(x)\ge 0 $$

Proposing then

$$ f_{\lambda}(x) = \max\left(0,\frac{1}{\lambda}-\alpha(x)\right) $$

now choosing $\lambda^*$ such that $\int_0^{x_{max}}f_{\lambda^*}(x)dx = 1$ we follow with the solution

$$ f_{\lambda^*}(x) = \max\left(0,\frac{1}{\lambda^*}-\alpha(x)\right) $$

Attached a MATHEMATICA script to handle the necessary computations.

Clear[g, f, binsearch]
alpha[t_] := (Sin[2 t] + Cos[3 t - 1/2] + 3)/5
xmax = 3 Pi;
f[x_, lambda_] := Max[0, 1/lambda - alpha[x]]
g[lambda_] := NIntegrate[f[x, lambda], {x, 0, xmax}] - 1
binsearch[xi0_, xs0_] := Module[{error = 0.0001, xi = xi0, xs = xs0, gi, gs, xm, gm},
  gi = g[xi];
  gs = g[xs];
  If[Abs[gi - gs] > error,
    xm = 0.5 (xi + xs);
    gm = g[xm];
    If[gm > 0, xi = xm, xs = xm];
    binsearch[xi, xs], Return[0.5 (xi + xs)]
  ]
]

lambda0 = binsearch[1, 3]
Plot[{alpha[x], f[x, lambda0]}, {x, 0, xmax}]
Plot[{alpha[x], alpha[x] + f[x, lambda0]}, {x, 0, xmax}]

Answer 2

I can't tell if your path is right (because it depends on that you want to do), but it is sensible.

Your problem structure is quite similar to one occurring in the field of optimal control. I'm not aware of any closed form solution of your approach, so you will have to pick a numerical approach to get started.

As you are probably aware, you can't express arbitrary functions over the real numbers in a computer with finite memory, so you will have to approximate your function using finite amount of information. I see two approaches being applicable here: collocation and function vector space composition. I will detail the second approach as results in simpler constraints.

We choose some basis of (preferably orthogonal) functions $G$ of length $m$ with sufficient properties and choose a set of real weights $h$ such that $f(x) = \sum_{i=1}^m h_iG_i(x)$.

$\forall x: f(x) \geq 0$ can be expressed as $\forall i\in \{1, ..., m\}: h_i \geq 0 $ and $\forall i\in \{1, ..., m\}, \forall x:G_i(x) \geq 0$.

$\quad \int_x f(x)dx=1$ can be expressed as a linear constraint $\sum_{i=1}^m h_i\int_x G_i(x)dx = 1$ where the integral can be evaluated ahead of time.

Now we are left with this constrained optimization problem:

\begin{equation} \begin{array}{ll} \underset{h}{\operatorname{minimize}} & -\displaystyle\int_{x} \log \left(\alpha(x)+\sum_{i=1}^m h_iG_i(x)\right)dx \\ \text { subject to } & \forall i\in \{1, ..., m\}: h_i \geq 0; \\ & \forall i\in \{1, ..., m\}, \forall x:G_i(x) \geq 0; \\ & \sum_{i=1}^m h_i\int_x G_i(x)dx = 1 \end{array} \end{equation}

This type of problem for a given basis can be solved using optimal control methods. For smooth problems i in general would recommend giving the optimizer of your choosing access to the derivatives of your objective in regards to $h$. One way to accomplish that is using sensitives of the ODE solutions using this observation that $$\frac{dS_j}{dt}=\frac{d}{dt}\frac{\partial u}{\partial p_{j}}=\frac{\partial f}{\partial y}\frac{\partial y}{\partial p_{j}}+\frac{\partial f}{\partial p_{j}}=J\cdot S_{j}+F_{j}$$ where $S_j(t)$ is the local sensitivity at point $t$. For more details also see this software manual.

By calculating the sensitivities you turned you Quadrature problem into an ODE problem but since the sensitivities are the gradients of your objective you can use a wider variety of constrained optimization algorithms.

Sadly implementing a practical and robust solution to this requires some engineering but i'm positive that within 5 years a plug-and-play solution will exist for this kind of problem.