**Adjustable Robust Optimization**This example solves an inventory problem (Ben-Tal et. al. 2004) that involves multi-stage decision-making. Consider a inventor system comprised of a warehouse and \(I\) factories in a horizon with \(T\) time steps. The uncertainty appears in the demand for the product at each time step \(t\), denoted by \(d_t\). An uncertainty set \(\mathcal{D}\) expressed as the following equation is used to define the upper and lower bounds of random demand.

\begin{align}

\mathcal{D}=\left\{\pmb{d}:~d_t\in\left[d_t^*-\delta d_t^*, d_t^*+\delta d_t^*\right],~\forall t=1,...,T\right\}

\end{align} where \(d_t^*\) is the expected demand at time \(t\), and \(\delta\) is a constant implying the deviation of random demand. In this example, \(\delta\) is set to be 0.2, and \(d_t^*\) follows the equation below.

\begin{align}

d_t^*=1000\left(1+\frac{1}{2}\sin\left(\frac{\pi(t-1)}{12}\right)\right), ~~~~~\forall t=1,...,T

\end{align} The amount of product in the warehouse at the beginning of time \(t\) under uncertainty outcome \(\pmb{d}\) is denoted by \(v_t(\pmb{d})\), and the amount of product produced during time \(t\) by factor \(i\) under uncertainty \(\pmb{d}\) is represented by \(p_{it}(\pmb{d})\). Let \(P_{it}\) be the corresponding maximum production capacity, and \(c_{it}\) be the cost of producing one unit of product, expressed by the following function.

\begin{align}

c_{it}=\alpha_i\left(1+\frac{1}{2}\left(\frac{\pi(t-1)}{12}\right)\right), ~~~\forall i=1,2,3,..I,~\forall t=1,2,3,...T

\end{align} The robust inventory problem can be formulated as follows.

\begin{align}

\min~&\max\limits_{\pmb{d}\in\mathcal{D}}\sum\limits_{t=1}^T\sum\limits_{i=1}^Tc_{it}p_{it}(\pmb{d})\\

\text{s.t.}~&0\leq p_{it}(\pmb{d})\leq P_{it}, ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\forall \pmb{d}\in\mathcal{D}, \forall i=1, 2, 3, ..., I,~\forall t\in 1, 2, 3, ..., T \\

&\sum\limits_{t=1}^Tp_{it}\leq Q_i, ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\forall \pmb{d}\in\mathcal{D}, \forall i=1, 2, 3, ..., I \\

&v_{t+1}(\pmb{d})=v_t(\pmb{d})+\sum\limits_{i}^Ip_{it}(\pmb{d})-d_t,, ~~~~~~\forall \pmb{d}\in\mathcal{D}, \forall t=1, 2, 3, ..., T \\

& V_{min}\leq v_{t}(\pmb{d})\leq V_{max}, , ~~~~~~~~~~~~~~~~~~~~~~~~~~\forall \pmb{d}\in\mathcal{D}, \forall t=1, 2, 3, ..., T \\

\end{align} where the constant \(I=3\), and the number of time steps is \(T=24\). The production capacity \(P_{it}=567\), and the total production cannot exceed the constant \(Q_i=13600\). The minimum and maximum capacity limitations of the warehouse are \(V_{max}=2000\) and \(V_{min}=500\), respectively.

The adjustable decisions \(\pmb{p}(\pmb{d})\) are approximated by linear decision rule function below, assuming "standard information basis" in paper (Ben-Tal et. al. 2004).

\begin{align}

p_{it}(\pmb{d})=\pi_{it}^0+\sum\limits_{\tau\in I_t}\pi_{it}^{\tau}d_t, ~~~~~\forall i=1,2,..,I, ~\forall t=1,2,3,..T

\end{align} Details of defining the decision rule are presented in the following code.

T=24; % number of time steps I=3; % number of factories; d0=1000*(1+1/2*sin(pi*((tt-1)/12))); % vector of expected demands c =[1;1.5;2]*(1+1/2*sin(pi*(tt-1)/12)); % matrix of production cost P=567; % production capacity Q=13600; % limitation of total products number Vmin=500; % mimimum capacity of the warehouse Vmax=2000; % maximum capacity of the warehouse V0 =1500; % initial storage of the warehouse delta=0.2; % constant indicating the deviations model=xprog('inventory'); % create a model, named 'inventory' d=model.random(1,T); % define random demand model.uncertain(d<=(1+delta)*d0); % upper bound of random demand model.uncertain(d>=(1-delta)*d0); % lower bound of random demand p=model.recourse(I,T); % define recourse decision as decision rule for t=2:T % define dependency in a for loop I_t=1:t-1; % standard information basis assumption p(:,t).depend(d,I_t); % dependency based on standard information basis end % the t-th column of p depends on the entries of d indexed by I_t model.min(sum(sum(c.*p))); % define objective function model.add(p>=0); % lower bound of p model.add(p<=P); % upper bound of p model.add(sum(p,2)<=Q) % total production capacity Q v=V0; % initial inventory for t=1:T v=v+(sum(p(:,t))-d(t)); % update inventory v_{t+1} model.add(v<=Vmax); model.add(v>=Vmin); end model.solve; % solve the problem

The objective value is 4206.9. It can be seen that the recourse decisions are approximated by the linear decision rule functions, and the dependency of these decision rules upon random variables can be conveniently defined by using the function "depend" in a "for" loop.