Optimal Power and Bandwidth Allocation in a Gaussian Channel ============================================================ by Robert Gowers, Roger Hill, Sami Al-Izzi, Timothy Pollington and Keith Briggs from Boyd and Vandenberghe, Convex Optimization, exercise 4.62 page 210 Consider a system in which a central node transmits messages to :math:n receivers. Each receiver channel :math:i \in \{1,...,n\} has a transmit power :math:P_i and bandwidth :math:W_i. A fraction of the total power and bandwidth is allocated to each channel, such that :math:\sum_{i=1}^{n}P_i = P_{tot} and :math:\sum_{i=1}^{n}W_i = W_{tot}. Given some utility function of the bit rate of each channel, :math:u_i(R_i), the objective is to maximise the total utility :math:U = \sum_{i=1}^{n}u_i(R_i). Assuming that each channel is corrupted by Gaussian white noise, the signal to noise ratio is given by :math:\beta_i P_i/W_i. This means that the bit rate is given by: :math:R_i = \alpha_i W_i \log_2(1+\beta_iP_i/W_i) where :math:\alpha_i and :math:\beta_i are known positive constants. One of the simplest utility functions is the data rate itself, which also gives a convex objective function. The optimisation problem can be thus be formulated as: minimise :math:\sum_{i=1}^{n}-\alpha_i W_i \log_2(1+\beta_iP_i/W_i) subject to :math:\sum_{i=1}^{n}P_i = P_{tot} \quad \sum_{i=1}^{n}W_i = W_{tot} \quad P \succeq 0 \quad W \succeq 0 Although this is a convex optimisation problem, it must be rewritten in DCP form since :math:P_i and :math:W_i are variables and DCP prohibits dividing one variable by another directly. In order to rewrite the problem in DCP format, we utilise the :math:\texttt{kl_div} function in CVXPY, which calculates the Kullback-Leibler divergence. :math:\text{kl_div}(x,y) = x\log(x/y)-x+y :math:-R_i = \text{kl_div}(\alpha_i W_i, \alpha_i(W_i+\beta_iP_i)) - \alpha_i\beta_iP_i Now that the objective function is in DCP form, the problem can be solved using CVXPY. .. code:: python #!/usr/bin/env python3 # @author: R. Gowers, S. Al-Izzi, T. Pollington, R. Hill & K. Briggs import numpy as np import cvxpy as cp .. code:: python def optimal_power(n, a_val, b_val, P_tot=1.0, W_tot=1.0): # Input parameters: α and β are constants from R_i equation n = len(a_val) if n != len(b_val): print('alpha and beta vectors must have same length!') return 'failed', np.nan, np.nan, np.nan P = cp.Variable(shape=n) W = cp.Variable(shape=n) alpha = cp.Parameter(shape=n) beta = cp.Parameter(shape=n) alpha.value = np.array(a_val) beta.value = np.array(b_val) # This function will be used as the objective so must be DCP; # i.e. elementwise multiplication must occur inside kl_div, # not outside otherwise the solver does not know if it is DCP... R = cp.kl_div(cp.multiply(alpha, W), cp.multiply(alpha, W + cp.multiply(beta, P))) - \ cp.multiply(alpha, cp.multiply(beta, P)) objective = cp.Minimize(cp.sum(R)) constraints = [P>=0.0, W>=0.0, cp.sum(P)-P_tot==0.0, cp.sum(W)-W_tot==0.0] prob = cp.Problem(objective, constraints) prob.solve() return prob.status, -prob.value, P.value, W.value Example ------- Consider the case where there are 5 channels, :math:n=5, :math:\alpha = \beta = (2.0,2.2,2.4,2.6,2.8), :math:P_{\text{tot}} = 0.5 and :math:W_{\text{tot}}=1. .. code:: python np.set_printoptions(precision=3) n = 5 # number of receivers in the system a_val = np.arange(10,n+10)/(1.0*n) # α b_val = np.arange(10,n+10)/(1.0*n) # β P_tot = 0.5 W_tot = 1.0 status, utility, power, bandwidth = optimal_power(n, a_val, b_val, P_tot, W_tot) print('Status: {}'.format(status)) print('Optimal utility value = {:.4g}'.format(utility)) print('Optimal power level:\n{}'.format(power)) print('Optimal bandwidth:\n{}'.format(bandwidth)) .. parsed-literal:: Status: optimal Optimal utility value = 2.451 Optimal power level: [1.151e-09 1.708e-09 2.756e-09 5.788e-09 5.000e-01] Optimal bandwidth: [3.091e-09 3.955e-09 5.908e-09 1.193e-08 1.000e+00]