## 4 Solution Approach  
The solution approach to the problem formulation is based on Lagrangian Relaxation. Constraints (IP 1), (IP 3) and (IP 5) are relaxed and multiplied by nonnegative Lagrangian multipliers, \(\mu^{i}_{*}\), \(\dot{\mu^{i}_{i}}\), and \(\mu^{j}_{*n}\), respectively. They are added to the objective functions as follows. The constraints are relaxed in such a way that the corresponding Lagrangian multipliers \(\mu^{i}_{*}\), \(\dot{\mu^{i}_{i}}\), and \(\mu^{j}_{*n}\) can be relaxed to objective function. This approach is not guaranteed to have a feasible solution due to accessing a big value of \(\gamma_{*}\) into the network which will result in unsatisfactory QoS requirement to the O-D pair \(w\). \[ \begin{array}{l}\hline\\ Z_{zz}(\mu^{i}_{*},\mu^{i}_{t},\mu^{i}_{*n})=\\ \min\sum_{i\in L}D_{i}(g_{t})g_{t}\\ +\sum_{i\in L}\mu^{i}_{*}(\sum_{j\in L}(g_{t})\gamma_{*n}-D_{*})\\ +\sum_{i\in L}\mu^{i}_{*n}(\sum_{j\in L}x_{*}\beta_{*n}\gamma_{*n}-g_{t})\\ +\sum_{i\in L}\mu^{i}_{*n}(\sum_{j\in L}x_{*}\beta_{*n}-\gamma_{*n})\\ \end{array} \] (LR) \[ \begin{array}{l}\\ \sum_{n=l_{*}}x_{*}=1\\ 0\leq g_{t}\leq\frac{C_{t}}{I(l)}\\ x_{*}=0\,\mathrm{or}\,1\\ y_{*n}=0\,\mathrm{or}\,1\\ \end{array} \] (IP 2) \[ \begin{array}{l}\\ \forall\,\in\,L\\ \forall\,\in\,P_{*},w\,\in\,W\\ \forall\,\in\,W,\forall\,\in\,L\,.\\ \end{array} \] (IP 6)  
We use the Lagrangian Relaxation Method to relax constraints in the formulation and decompose (LR) into two independent sub-problems. In the first Sub-Problem (SP 1), the nonnegative weight is calculated as \((\mu^{i}_{*i}+\mu^{i}_{t}\gamma_{*})\) on each link. The actual meaning of the multipliers \(\mu^{i}_{*n}\) and \(\mu^{i}_{t}\) are the link mean delay and the derivative of queue length, respectively; both can be derived from the solution of the problem (LR). By the problem (LR), the multiplier \(\mu_{\omega}^{j}\) is related to the relaxation of link selection for each O-D pair. When one unit on the decision variable of link selection is added (_i.e._, \(y_{wl}\)), it generates one unit of the corresponding mean delay as formulated in QoS constraint, but when the problem (LR) is solved, the decision variable \(y_{wl}\) is either 1 or 0 not only one unit, thus, the multiplier \(\mu_{\omega}^{j}\) is equivalent to the mean delay on the chosen link. Besides, the link selection is based on the sum of link mean delay \(D_{\lambda}(g_{\lambda})\) of each O-D pair which is weighted by \(\mu_{\omega}^{j}\), that is, how much the end-to-end delay will impact the objective function (IP). Moreover, in the problem (SP 2.1'), the first step for solving it is also related to the multiplier \(\mu_{\omega}^{j}\) and chooses the links for each O-D pair according to \(y_{\omega}^{\cdot}(g_{\lambda})=\mu_{\omega}^{j}D_{\lambda}(g_{\lambda})-\mu_ {\omega}^{j}=0\). It means that whenever a link is selected for an O-D pair, the link will produce corresponding link mean delay \(\mu_{\omega}^{j}\) and affect the end-to-end delay constraint which is weighted by multiplier \(\mu_{\omega}^{j}\) in the problem (LR). On the other hand, the multiplier \(\mu_{i}^{j}\) is related to the relaxation of aggregate flow over each link. The physical meaning of \(\mu_{i}^{j}\) would be that when \(\mu_{i}^{j}\) can be relaxed on the link flow, it will affect the objective function (IP), which means \(\mu_{i}^{j}\) would be the sensitivity of the average number of packets over the link with respect to the aggregate flow. A brief illustration is shown in Table 4, \(\mu_{\omega}^{j}\) implies the link mean delay of each link, and \(\mu_{i}^{2}\) indicates the derivative of queue length. The term \(\gamma_{\omega}\) in the arc weight form represents the weighting factor between \(\mu_{\omega}^{j}\) and \(\mu_{i}^{2}\). Thus, this metric consists of the link mean delay and the derivative of queue length, considers both system perspective and user perspective for our distributed routing protocol, and uses the traffic requirement of each O-D pair as the weighting factor to combine these two parameters.