001/*******************************************************************************
002 * Copyright (c) 2016 Pablo Pavon Mariņo.
003 * All rights reserved. This program and the accompanying materials
004 * are made available under the terms of the GNU Lesser Public License v2.1
005 * which accompanies this distribution, and is available at
006 * http://www.gnu.org/licenses/lgpl.html
007 ******************************************************************************/
008package com.net2plan.examples.ocnbook.onlineSim;
009
010
011
012
013import java.io.File;
014import java.util.HashMap;
015import java.util.List;
016import java.util.Map;
017import java.util.Random;
018
019import cern.colt.matrix.tdouble.DoubleFactory1D;
020import cern.colt.matrix.tdouble.DoubleFactory2D;
021import cern.colt.matrix.tdouble.DoubleMatrix1D;
022import cern.colt.matrix.tdouble.DoubleMatrix2D;
023import cern.jet.math.tdouble.DoubleFunctions;
024
025import com.net2plan.examples.ocnbook.offline.Offline_cba_congControLinkBwSplitTwolQoS;
026import com.net2plan.interfaces.networkDesign.Demand;
027import com.net2plan.interfaces.networkDesign.Link;
028import com.net2plan.interfaces.networkDesign.Net2PlanException;
029import com.net2plan.interfaces.networkDesign.NetPlan;
030import com.net2plan.interfaces.networkDesign.Node;
031import com.net2plan.interfaces.networkDesign.Route;
032import com.net2plan.interfaces.simulation.IEventProcessor;
033import com.net2plan.interfaces.simulation.SimEvent;
034import com.net2plan.utils.Constants.RoutingType;
035import com.net2plan.utils.GradientProjectionUtils;
036import com.net2plan.utils.InputParameter;
037import com.net2plan.utils.Pair;
038import com.net2plan.utils.TimeTrace;
039import com.net2plan.utils.Triple;
040
041/** 
042 * This module implements a distributed primal-decomposition-based gradient algorithm, for a coordinated adjustment of the congestion control of two types of demands (with different utility functions), and the fraction of each link capacity to grant to the traffic of each type, to maximize the network utility enforcing a fair allocation of the resources.
043 *
044 * Ths event processor is adapted to permit observing the algorithm performances under user-defined conditions, 
045 * including asynchronous distributed executions, where signaling can be affected by losses and/or delays, and/or measurement errors. 
046 * The time evolution of different metrics can be stored in output files, for later processing. 
047 * As an example, see the <a href="../../../../../../graphGeneratorFiles/fig_sec11_3_congControlAndQoSLinkCap_primalDecomp.m">{@code fig_sec11_3_congControlAndQoSLinkCap_primalDecomp.m}</a> MATLAB file used for generating the graph/s of the case study in the 
048 * <a href="http://eu.wiley.com/WileyCDA/WileyTitle/productCd-1119013356.html">book</a> using this algorithm.
049 * 
050 * To simulate a network with this module, use the {@code Online_evGen_doNothing} generator.
051 * 
052 * @net2plan.keywords Bandwidth assignment (BA), Capacity assignment (CA), Distributed algorithm, Primal decomposition 
053 * @net2plan.ocnbooksections Section 11.3
054 * @net2plan.inputParameters 
055 * @author Pablo Pavon-Marino
056 */
057public class Online_evProc_congControlAndQoSTwoClassesPrimalDecomp extends IEventProcessor 
058{
059        private double PRECISIONFACTOR;
060        private Random rng;
061
062        private int N,E,D1,D2,D;
063        private Map<String,String> algorithmParameters , net2PlanParameters;
064        
065        private static final int MAC_UPDATE_WAKEUPTOUPDATE = 202;
066        private static final int CC_SIGNALING_WAKEUPTOSENDMESSAGE = 400;
067        private static final int CC_SIGNALING_RECEIVEDMESSAGE = 401;
068        private static final int CC_UPDATE_WAKEUPTOUPDATE = 402;
069
070        private InputParameter mac_update_isSynchronous = new InputParameter ("mac_update_isSynchronous", false , "true if all the distributed agents involved wake up synchronousely to update its state");
071        private InputParameter mac_update_averageInterUpdateTime = new InputParameter ("mac_update_averageInterUpdateTime", 50.0 , "Average time between two updates of an agent" , 0 , false , Double.MAX_VALUE , true);
072        private InputParameter mac_update_maxFluctuationInterUpdateTime = new InputParameter ("mac_update_maxFluctuationInterUpdateTime", 10.0 , "Max fluctuation in time in the update interval of an agent, in absolute time values. The update intervals are sampled from a uniform distribution within the given interval" , 0 , true , Double.MAX_VALUE , true);
073        private InputParameter mac_gradient_gammaStep = new InputParameter ("mac_gradient_gammaStep", 5.0 , "Gamma step in the gradient algorithm" , 0 , false , Double.MAX_VALUE , true);
074        private InputParameter mac_gradient_maxGradientCoordinateChange = new InputParameter ("mac_gradient_maxGradientCoordinateChange", 1000.0 , "Maximum change in an iteration of a gradient coordinate" , 0 , false , Double.MAX_VALUE , true);
075
076        private InputParameter mac_simulation_maxNumberOfUpdateIntervals = new InputParameter ("mac_simulation_maxNumberOfUpdateIntervals", 600.0 , "Maximum number of update intervals in average per agent" , 0 , false , Double.MAX_VALUE , true);
077        private InputParameter simulation_randomSeed = new InputParameter ("simulation_randomSeed", (long) 1 , "Seed of the random number generator");
078        private InputParameter simulation_outFileNameRoot = new InputParameter ("simulation_outFileNameRoot", "crossLayerCongControlTwoQoSPrimalDecomp" , "Root of the file name to be used in the output files. If blank, no output");
079        
080        private InputParameter cc_signaling_isSynchronous = new InputParameter ("cc_signaling_isSynchronous", false , "true if all the distributed agents involved wake up synchronously to send the signaling messages");
081        private InputParameter cc_signaling_averageInterMessageTime = new InputParameter ("cc_signaling_averageInterMessageTime", 1.0 , "Average time between two signaling messages sent by an agent" , 0 , false , Double.MAX_VALUE , true);
082        private InputParameter cc_signaling_maxFluctuationInterMessageTime = new InputParameter ("cc_signaling_maxFluctuationInterMessageTime", 0.5 , "Max fluctuation in time between two signaling messages sent by an agent" , 0 , true , Double.MAX_VALUE , true);
083        private InputParameter cc_signaling_averageDelay = new InputParameter ("cc_signaling_averageDelay", 3.0 , "Average time between signaling message transmission by an agent and its reception by other or others" , 0 , true , Double.MAX_VALUE , true);
084        private InputParameter cc_signaling_maxFluctuationInDelay = new InputParameter ("cc_signaling_maxFluctuationInDelay", 0.5 , "Max fluctuation in time in the signaling delay, in absolute time values. The signaling delays are sampled from a uniform distribution within the given interval" , 0 , true , Double.MAX_VALUE , true);
085        private InputParameter cc_signaling_signalingLossProbability = new InputParameter ("cc_signaling_signalingLossProbability", 0.05 , "Probability that a signaling message transmitted is lost (not received by other or others involved agents)" , 0 , true , Double.MAX_VALUE , true);
086
087        private InputParameter cc_update_isSynchronous = new InputParameter ("cc_update_isSynchronous", false , "true if all the distributed agents involved wake up synchronousely to update its state");
088        private InputParameter cc_update_averageInterUpdateTime = new InputParameter ("cc_update_averageInterUpdateTime", 1.0 , "Average time between two updates of an agent" , 0 , false , Double.MAX_VALUE , true);
089        private InputParameter cc_update_maxFluctuationInterUpdateTime = new InputParameter ("cc_update_maxFluctuationInterUpdateTime", 0.5 , "Max fluctuation in time in the update interval of an agent, in absolute time values. The update intervals are sampled from a uniform distribution within the given interval" , 0 , true , Double.MAX_VALUE , true);
090
091        private InputParameter cc_gradient_gammaStep = new InputParameter ("cc_gradient_gammaStep", 0.001 , "Gamma step in the gradient algorithm" , 0 , false , Double.MAX_VALUE , true);
092        private InputParameter cc_gradient_maxGradientAbsoluteNoise = new InputParameter ("cc_gradient_maxGradientAbsoluteNoise", 0.0 , "Max value of the added noise to the gradient coordinate in absolute values" , 0 , true , Double.MAX_VALUE , true);
093
094        private InputParameter cc_simulation_maxNumberOfUpdateIntervals = new InputParameter ("cc_simulation_maxNumberOfUpdateIntervals", 600.0 , "Maximum number of update intervals in average per agent" , 0 , false , Double.MAX_VALUE , true);
095        private InputParameter cc_control_minHd = new InputParameter ("cc_control_minHd", 0.1 , "Minimum traffic assigned to each demand" , 0 , true , Double.MAX_VALUE , true);
096        private InputParameter cc_control_maxHd = new InputParameter ("cc_control_maxHd", 1e6 , "Maximum traffic assigned to each demand" , 0 , true , Double.MAX_VALUE , true);
097        private InputParameter cc_control_initialiLinkPrices = new InputParameter ("cc_control_initialiLinkPrices", 1.0 , "Link prices in the algorithm initialization" , 0 , true , Double.MAX_VALUE , true);
098
099        private InputParameter cc_control_fairnessFactor_1 = new InputParameter ("cc_control_fairnessFactor_1", 1.0 , "Fairness factor in utility function of congestion control for demands of class 1" , 0 , true , Double.MAX_VALUE , true);
100        private InputParameter cc_control_fairnessFactor_2 = new InputParameter ("cc_control_fairnessFactor_2", 1.0 , "Fairness factor in utility function of congestion control for demands of class 2" , 0 , true , Double.MAX_VALUE , true);
101        private InputParameter cc_control_weightFairness_1 = new InputParameter ("cc_control_weightFairness_1", 1.0 , "Weight factor in utility function demands type 1" , 0 , true , Double.MAX_VALUE , true);
102        private InputParameter cc_control_weightFairness_2 = new InputParameter ("cc_control_weightFairness_2", 2.0 , "Weight factor in utility function demands type 2" , 0 , true , Double.MAX_VALUE , true);
103        private InputParameter mac_minCapacity_1 = new InputParameter ("mac_minCapacity_1", 0.1 , "Minimum capacity in each link, allocated to traffic of type 1" , 0 , true , Double.MAX_VALUE , true);
104        private InputParameter mac_minCapacity_2 = new InputParameter ("mac_minCapacity_2", 0.1 , "Minimum capacity in each link, allocated to traffic of type 2" , 0 , true , Double.MAX_VALUE , true);
105
106        private DoubleMatrix1D demandType;
107        private DoubleMatrix1D congControl_price1_e;
108        private DoubleMatrix1D congControl_price2_e;
109        private DoubleMatrix2D control_mostUpdatedLinkPriceKnownDemand_de;
110        private DoubleMatrix1D mac_u1; 
111        private DoubleMatrix1D mac_u2; 
112        
113        private TimeTrace traceOf_pi1_e;
114        private TimeTrace traceOf_pi2_e;
115        private TimeTrace traceOf_u1_e;
116        private TimeTrace traceOf_u2_e;
117        private TimeTrace traceOf_y1_e;
118        private TimeTrace traceOf_y2_e;
119        private TimeTrace traceOf_h_d1;
120        private TimeTrace traceOf_h_d2;
121        private TimeTrace traceOf_objFunction;
122
123        private NetPlan currentNetPlan , copyInitialNetPlan;
124        
125        
126        @Override
127        public String getDescription()
128        {
129                return "This module implements a distributed primal-decomposition-based gradient algorithm, for a coordinated adjustment of the congestion control of two types of demands (with different utility functions), and the fraction of each link capacity to grant to the traffic of each type, to maximize the network utility enforcing a fair allocation of the resources.";
130        }
131
132        @Override
133        public List<Triple<String, String, String>> getParameters()
134        {
135                /* Returns the parameter information for all the InputParameter objects defined in this object (uses Java reflection) */
136                return InputParameter.getInformationAllInputParameterFieldsOfObject(this);
137        }
138
139        @Override
140        public void initialize(NetPlan currentNp, Map<String, String> algorithmParameters, Map<String, String> simulationParameters, Map<String, String> net2planParameters)
141        {
142                /* Initialize all InputParameter objects defined in this object (this uses Java reflection) */
143                InputParameter.initializeAllInputParameterFieldsOfObject(this, algorithmParameters);
144
145                this.currentNetPlan = currentNp;
146                this.algorithmParameters = algorithmParameters;
147                this.net2PlanParameters = net2planParameters;
148                if (currentNetPlan.getNumberOfLayers() != 1) throw new Net2PlanException ("This algorithm works in single layer networks");
149
150                /* INITIALIZE VARIABLES COMMON FOR BOTH LAYERS */
151                this.copyInitialNetPlan = currentNp.copy();
152                if (currentNetPlan.getVectorLinkCapacity().getMaxLocation() [0]< mac_minCapacity_1.getDouble() + mac_minCapacity_2.getDouble()) throw new Net2PlanException ("Minimum link capacities for each class are too high for the available link capacities");
153                
154                this.PRECISIONFACTOR = Double.parseDouble(net2planParameters.get("precisionFactor"));
155                
156                /* READ INPUT PARAMETERS BOTH LAYERS */
157                this.rng = new Random (simulation_randomSeed.getLong());
158                this.N = this.currentNetPlan.getNumberOfNodes();
159                this.E = this.currentNetPlan.getNumberOfLinks();
160                if (E == 0) throw new Net2PlanException ("The input design should have links");
161                this.D1 = N*(N-1);
162                this.D2 = N*(N-1);
163                this.D = D1+D2;
164                
165                /* Initialize the demands and routers per demand */
166                this.currentNetPlan.removeAllDemands();
167                this.currentNetPlan.setRoutingType(RoutingType.SOURCE_ROUTING);
168                this.demandType = DoubleFactory1D.dense.make (D1+D2);
169                for (Node n1 : this.currentNetPlan.getNodes())
170                        for (Node n2 : this.currentNetPlan.getNodes())
171                                if (n1 != n2) 
172                                { 
173                                        final Demand d1 = this.currentNetPlan.addDemand(n1, n2, 0.0, null); d1.setAttribute("type" , "1"); demandType.set(d1.getIndex (), 1); 
174                                        final Demand d2 = this.currentNetPlan.addDemand(n1, n2, 0.0, null); d1.setAttribute("type" , "2"); demandType.set(d2.getIndex (), 2); 
175                                }
176                
177                /* Remove all routes, and create one with the shortest path in km for each demand */
178                currentNetPlan.removeAllUnicastRoutingInformation();
179                currentNetPlan.setRoutingType(RoutingType.SOURCE_ROUTING);
180                currentNetPlan.addRoutesFromCandidatePathList(currentNetPlan.getVectorLinkLengthInKm().toArray()  , "K" , "1");
181                for (Route r : currentNp.getRoutes ()) r.setCarriedTraffic(cc_control_minHd.getDouble() , cc_control_minHd.getDouble());
182                
183                /* INITIALIZE control information  */
184                this.congControl_price1_e = DoubleFactory1D.dense.make (E , cc_control_initialiLinkPrices.getDouble());
185                this.congControl_price2_e = DoubleFactory1D.dense.make (E , cc_control_initialiLinkPrices.getDouble());
186                
187                /* Initialize the information each demand knows of the prices of all the links */
188                this.control_mostUpdatedLinkPriceKnownDemand_de = DoubleFactory2D.dense.make (D,E,cc_control_initialiLinkPrices.getDouble());
189
190                this.mac_u1 = DoubleFactory1D.dense.make (E , mac_minCapacity_1.getDouble());
191                this.mac_u2 = currentNetPlan.getVectorLinkCapacity().copy ().assign(DoubleFunctions.minus (mac_minCapacity_1.getDouble()));
192                
193                /* Each demand adjusts its rate to the pi_e */
194                /* Compute the traffic each demand injects, set it in net2plan, and tell routing layer there was an update in this */
195                for (Demand d : currentNetPlan.getDemands())
196                {
197                        final double h_d = this.computeHdFromPrices(d);
198                        d.setOfferedTraffic(h_d);
199                        d.getRoutes().iterator().next ().setCarriedTraffic(h_d , h_d);
200                }
201                
202                /* Initially all nodes receive a "wake up to transmit" event, aligned at time zero or y asynchr => randomly chosen */
203                for (Link e : currentNp.getLinks())
204                {
205                        final double updateTime = (mac_update_isSynchronous.getBoolean())? mac_update_averageInterUpdateTime.getDouble() : Math.max(0 , mac_update_averageInterUpdateTime.getDouble() + mac_update_maxFluctuationInterUpdateTime.getDouble() * (rng.nextDouble() - 0.5));
206                        this.scheduleEvent(new SimEvent (updateTime , SimEvent.DestinationModule.EVENT_PROCESSOR , MAC_UPDATE_WAKEUPTOUPDATE , e));
207                }
208                
209                /* INITIALIZATION CONGESTION CONTROL LAYER */
210                /* Initially all nodes receive a "wake up to transmit" event, aligned at time zero or y asynchr => randomly chosen */
211                for (Link e : currentNp.getLinks())
212                {
213                        final double signalingTime = (cc_signaling_isSynchronous.getBoolean())? cc_signaling_averageInterMessageTime.getDouble() : Math.max(0 , cc_signaling_averageInterMessageTime.getDouble() + cc_signaling_maxFluctuationInterMessageTime.getDouble() * (rng.nextDouble() - 0.5));
214                        this.scheduleEvent(new SimEvent (signalingTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_SIGNALING_WAKEUPTOSENDMESSAGE , e));
215                }
216                for (Demand d : currentNetPlan.getDemands())
217                {
218                        final double updateTime = (cc_update_isSynchronous.getBoolean())? cc_update_averageInterUpdateTime.getDouble() : Math.max(0 , cc_update_averageInterUpdateTime.getDouble() + cc_update_maxFluctuationInterUpdateTime.getDouble() * (rng.nextDouble() - 0.5));
219                        this.scheduleEvent(new SimEvent (updateTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_UPDATE_WAKEUPTOUPDATE , d));
220                }
221                
222                /* INITIALIZE STATISTIC VARIABLES FOR BOTH LAYERS */
223                this.traceOf_pi1_e = new TimeTrace  ();
224                this.traceOf_pi2_e = new TimeTrace  ();
225                this.traceOf_u1_e = new TimeTrace ();
226                this.traceOf_u2_e = new TimeTrace ();
227                this.traceOf_y1_e = new TimeTrace (); 
228                this.traceOf_y2_e = new TimeTrace (); 
229                this.traceOf_h_d1 = new TimeTrace (); 
230                this.traceOf_h_d2 = new TimeTrace (); 
231                this.traceOf_objFunction = new TimeTrace ();
232
233                
234                this.traceOf_pi1_e.add(0.0, this.congControl_price1_e.copy ());
235                this.traceOf_pi2_e.add(0.0, this.congControl_price2_e.copy ());
236                this.traceOf_u1_e.add(0.0, this.mac_u1.copy ());
237                this.traceOf_u2_e.add(0.0, this.mac_u2.copy ());
238                final Triple<Double,Pair<DoubleMatrix1D,DoubleMatrix1D>,Pair<DoubleMatrix1D,DoubleMatrix1D>> objFunctionInfo = computeObjectiveFunctionFromNetPlan(this.currentNetPlan);
239                this.traceOf_objFunction.add(0.0 , objFunctionInfo.getFirst());
240                this.traceOf_h_d1.add(0.0, objFunctionInfo.getSecond().getFirst().copy ());
241                this.traceOf_h_d2.add(0.0, objFunctionInfo.getSecond().getSecond().copy ());
242                this.traceOf_y1_e.add(0.0, objFunctionInfo.getThird().getFirst().copy ());
243                this.traceOf_y2_e.add(0.0, objFunctionInfo.getThird().getSecond().copy ());
244        }
245
246        @Override
247        public void processEvent(NetPlan currentNetPlan, SimEvent event)
248        {
249                final double t = event.getEventTime();
250                switch (event.getEventType())
251                {
252                case MAC_UPDATE_WAKEUPTOUPDATE: 
253                {
254                        final Link eMe = (Link) event.getEventObject();
255                        final double current_u1 = mac_u1.get(eMe.getIndex ());
256                        final double current_u2 = mac_u2.get(eMe.getIndex ());
257                        final double gradient_1 = congControl_price1_e.get(eMe.getIndex ());
258                        final double gradient_2 = congControl_price2_e.get(eMe.getIndex ());
259                        final double nextBeforeProjection_u1 = current_u1 + this.mac_gradient_gammaStep.getDouble() * gradient_1;
260                        final double nextBeforeProjection_u2 = current_u2 + this.mac_gradient_gammaStep.getDouble() * gradient_2;
261                        DoubleMatrix1D initialSolution = DoubleFactory1D.dense.make (new double [] { current_u1 , current_u2 });
262                        DoubleMatrix1D inOut = DoubleFactory1D.dense.make (new double [] { nextBeforeProjection_u1 , nextBeforeProjection_u2 });
263                        DoubleMatrix1D uMin = DoubleFactory1D.dense.make (new double [] { mac_minCapacity_1.getDouble() , mac_minCapacity_2.getDouble() });
264                        GradientProjectionUtils.euclideanProjection_sumInequality (inOut , null , uMin , eMe.getCapacity());
265                        if (mac_gradient_maxGradientCoordinateChange.getDouble() > 0)
266                                GradientProjectionUtils.scaleDown_maxAbsoluteCoordinateChange (initialSolution , inOut , null , mac_gradient_maxGradientCoordinateChange.getDouble());
267
268                        this.mac_u1.set(eMe.getIndex (), inOut.get(0));
269                        this.mac_u2.set(eMe.getIndex (), inOut.get(1));
270                        if (mac_u1.get(eMe.getIndex ()) + mac_u2.get(eMe.getIndex ()) > eMe.getCapacity() + PRECISIONFACTOR) throw new RuntimeException ("Bad");
271                        
272                        final double updateTime = mac_update_isSynchronous.getBoolean()? t + mac_update_averageInterUpdateTime.getDouble() : Math.max(t , t + mac_update_averageInterUpdateTime.getDouble() + mac_update_maxFluctuationInterUpdateTime.getDouble() * (rng.nextDouble() - 0.5));
273                        this.scheduleEvent(new SimEvent (updateTime , SimEvent.DestinationModule.EVENT_PROCESSOR , MAC_UPDATE_WAKEUPTOUPDATE,  eMe));
274
275                        this.traceOf_u1_e.add(t, this.mac_u1.copy ());
276                        this.traceOf_u2_e.add(t, this.mac_u2.copy ());
277
278                        if (t > this.mac_simulation_maxNumberOfUpdateIntervals.getDouble() * this.mac_update_averageInterUpdateTime.getDouble()) { this.endSimulation (); }
279
280                        break;
281                }
282
283                case CC_SIGNALING_RECEIVEDMESSAGE: // A node receives from an out neighbor the q_nt for any destination
284                {
285                        final Pair<Demand,Triple<Link,Double,Double>> signalInfo = (Pair<Demand,Triple<Link,Double,Double>>) event.getEventObject();
286                        final Demand d = signalInfo.getFirst();
287                        final int type = (int) demandType.get(d.getIndex ());
288                        final Link e = signalInfo.getSecond().getFirst();
289                        final double pi1 = signalInfo.getSecond().getSecond();
290                        final double pi2 = signalInfo.getSecond().getThird();
291                        this.control_mostUpdatedLinkPriceKnownDemand_de.set(d.getIndex () , e.getIndex () , (type == 1)? pi1 : pi2);  
292                        break;
293                }
294                
295                case CC_SIGNALING_WAKEUPTOSENDMESSAGE: // A node broadcasts signaling info to its 1 hop neighbors
296                {
297                        final Link eMe = (Link) event.getEventObject();
298
299                        /* Update the new price with the gradient approach */
300                        Triple<Double,Pair<DoubleMatrix1D,DoubleMatrix1D>,Pair<DoubleMatrix1D,DoubleMatrix1D>>  triple = computeObjectiveFunctionFromNetPlan (this.currentNetPlan);
301                        final double u1 = this.mac_u1.get(eMe.getIndex ());
302                        final double u2 = this.mac_u2.get(eMe.getIndex ());
303                        final double y1 = triple.getThird().getFirst().get(eMe.getIndex ());
304                        final double y2 = triple.getThird().getSecond().get(eMe.getIndex ());
305                        final double old_pie_1 = this.congControl_price1_e.get(eMe.getIndex ());
306                        final double old_pie_2 = this.congControl_price2_e.get(eMe.getIndex ());
307                        final double new_pie_1 = Math.max(0, old_pie_1 - this.cc_gradient_gammaStep.getDouble() * (u1 - y1) + 2*cc_gradient_maxGradientAbsoluteNoise.getDouble()*(rng.nextDouble()-0.5));
308                        final double new_pie_2 = Math.max(0, old_pie_2 - this.cc_gradient_gammaStep.getDouble() * (u2 - y2) + 2*cc_gradient_maxGradientAbsoluteNoise.getDouble()*(rng.nextDouble()-0.5));
309                        this.congControl_price1_e.set(eMe.getIndex (), new_pie_1);
310                        this.congControl_price2_e.set(eMe.getIndex (), new_pie_2);
311                        
312                        this.traceOf_pi1_e.add(t, this.congControl_price1_e.copy ());
313                        this.traceOf_pi2_e.add(t, this.congControl_price2_e.copy ());
314
315                        /* Create the info I will signal */
316                        Triple<Link,Double,Double> infoToSignal = Triple.of(eMe ,  new_pie_1 , new_pie_2);
317
318                        /* Send the events of the signaling information messages to all the nodes */
319                        for (Route route : eMe.getTraversingRoutes())
320                        {
321                                if (rng.nextDouble() < this.cc_signaling_signalingLossProbability.getDouble()) continue; // the signaling may be lost => lost to all demands
322                                final Demand d = route.getDemand ();
323                                final double signalingReceptionTime = t + Math.max(0 , cc_signaling_averageDelay.getDouble() + cc_signaling_maxFluctuationInDelay.getDouble() * (rng.nextDouble() - 0.5));
324                                this.scheduleEvent(new SimEvent (signalingReceptionTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_SIGNALING_RECEIVEDMESSAGE , Pair.of(d , infoToSignal)));
325                        }
326                        
327                        /* Re-schedule when to wake up again */
328                        final double signalingTime = cc_signaling_isSynchronous.getBoolean()? t + cc_signaling_averageInterMessageTime.getDouble() : Math.max(t , t + cc_signaling_averageInterMessageTime.getDouble() + cc_signaling_maxFluctuationInterMessageTime.getDouble() * (rng.nextDouble() - 0.5));
329                        this.scheduleEvent(new SimEvent (signalingTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_SIGNALING_WAKEUPTOSENDMESSAGE , eMe));
330                        break;
331                }
332
333                case CC_UPDATE_WAKEUPTOUPDATE: 
334                {
335                        final Demand dMe = (Demand) event.getEventObject();
336                        final double new_hd = computeHdFromPrices (dMe);
337                        
338                        dMe.setOfferedTraffic(new_hd);
339                        dMe.getRoutes ().iterator().next().setCarriedTraffic(new_hd , new_hd);
340                        
341                        final Triple<Double,Pair<DoubleMatrix1D,DoubleMatrix1D>,Pair<DoubleMatrix1D,DoubleMatrix1D>> objFunctionInfo = computeObjectiveFunctionFromNetPlan(this.currentNetPlan);
342                        this.traceOf_objFunction.add(t , objFunctionInfo.getFirst());
343                        this.traceOf_h_d1.add(t, objFunctionInfo.getSecond().getFirst().copy ());
344                        this.traceOf_h_d2.add(t, objFunctionInfo.getSecond().getSecond().copy ());
345                        this.traceOf_y1_e.add(t, objFunctionInfo.getThird().getFirst().copy ());
346                        this.traceOf_y2_e.add(t, objFunctionInfo.getThird().getSecond().copy ());
347
348                        final double updateTime = cc_update_isSynchronous.getBoolean()? t + cc_update_averageInterUpdateTime.getDouble() : Math.max(t , t + cc_update_averageInterUpdateTime.getDouble() + cc_update_maxFluctuationInterUpdateTime.getDouble() * (rng.nextDouble() - 0.5));
349                        this.scheduleEvent(new SimEvent (updateTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_UPDATE_WAKEUPTOUPDATE,  dMe));
350
351                        if (t > this.cc_simulation_maxNumberOfUpdateIntervals.getDouble() * this.cc_update_averageInterUpdateTime.getDouble()) { this.endSimulation (); }
352
353                        break;
354                }
355                
356                
357                default: throw new RuntimeException ("Unexpected received event");
358                }
359                
360                
361        }
362
363        public String finish (StringBuilder st , double simTime)
364        {
365                if (simulation_outFileNameRoot.getString().equals("")) return null;
366                traceOf_h_d1.printToFile(new File (simulation_outFileNameRoot.getString() + "_hd1.txt"));
367                traceOf_h_d2.printToFile(new File (simulation_outFileNameRoot.getString() + "_hd2.txt"));
368                traceOf_u1_e.printToFile(new File (simulation_outFileNameRoot.getString() + "_ue1.txt"));
369                traceOf_u2_e.printToFile(new File (simulation_outFileNameRoot.getString() + "_ue2.txt"));
370                traceOf_objFunction.printToFile(new File (simulation_outFileNameRoot.getString() + "_objFunc.txt"));
371
372                Map<String,String> param = new HashMap<String,String> (algorithmParameters);
373                param.put("solverName", "ipopt");
374                param.put("solverLibraryName", "");
375                param.put("maxSolverTimeInSeconds", "-1");
376                new Offline_cba_congControLinkBwSplitTwolQoS().executeAlgorithm(copyInitialNetPlan , param , this.net2PlanParameters);
377                DoubleMatrix1D jom_hd1 = DoubleFactory1D.dense.make (D1);
378                DoubleMatrix1D jom_hd2 = DoubleFactory1D.dense.make (D2);
379                DoubleMatrix1D jom_ue1 = DoubleFactory1D.dense.make (E);
380                DoubleMatrix1D jom_ue2 = DoubleFactory1D.dense.make (E);
381                int counter_1 = 0; int counter_2 = 0;
382                for (Demand d : copyInitialNetPlan.getDemands ())
383                        if (d.getAttribute("type").equals ("1")) jom_hd1.set (counter_1 ++ , d.getOfferedTraffic()); else jom_hd2.set (counter_2 ++ , d.getOfferedTraffic());
384                for (Link e : copyInitialNetPlan.getLinks ()) { jom_ue1.set (e.getIndex () , Double.parseDouble (e.getAttribute("u_1"))); jom_ue2.set (e.getIndex () , Double.parseDouble (e.getAttribute("u_2"))); }
385                final double jomObjcFunc = Double.parseDouble (copyInitialNetPlan.getAttribute("netUtility"));
386                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_hd1.txt"), jom_hd1);
387                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_hd2.txt"), jom_hd2);
388                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_ue1.txt"), jom_ue1);
389                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_ue2.txt"), jom_ue2);
390                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_objFunc.txt"), jomObjcFunc);
391                return null;
392        }
393        
394
395        private Triple<Double,Pair<DoubleMatrix1D,DoubleMatrix1D>,Pair<DoubleMatrix1D,DoubleMatrix1D>>  computeObjectiveFunctionFromNetPlan (NetPlan np)
396        {
397                DoubleMatrix1D y1 = DoubleFactory1D.dense.make (E);
398                DoubleMatrix1D y2 = DoubleFactory1D.dense.make (E);
399                DoubleMatrix1D h1 = DoubleFactory1D.dense.make (D1);
400                DoubleMatrix1D h2 = DoubleFactory1D.dense.make (D2);
401                double obj = 0;
402                int counter_d1 = 0; int counter_d2 = 0;
403                for (Demand d : np.getDemands())
404                {
405                        final double h_d = d.getCarriedTraffic();
406                        final int type = (int) this.demandType.get(d.getIndex());
407                        final List<Link> seqLinks = d.getRoutes ().iterator().next().getSeqLinksRealPath();
408                        for (Link e : seqLinks)
409                        {
410                                final int index_e = e.getIndex ();
411                                if (type == 1) y1.set(index_e , y1.get(index_e) + h_d); else y2.set(index_e , y2.get(index_e) + h_d); 
412                        }
413                        if (type == 1) 
414                        { 
415                                h1.set (counter_d1 , h_d);
416                                obj += (this.cc_control_fairnessFactor_1.getDouble() == 1)? this.cc_control_weightFairness_1.getDouble() * Math.log(h_d) : this.cc_control_weightFairness_1.getDouble() *  Math.pow(h_d, 1-this.cc_control_fairnessFactor_1.getDouble()) / (1-this.cc_control_fairnessFactor_1.getDouble());
417                                counter_d1 ++; 
418                        } else 
419                        { 
420                                h2.set (counter_d2 , h_d);
421                                obj += (this.cc_control_fairnessFactor_2.getDouble() == 1)? this.cc_control_weightFairness_2.getDouble() * Math.log(h_d) : this.cc_control_weightFairness_2.getDouble() *  Math.pow(h_d, 1-this.cc_control_fairnessFactor_2.getDouble()) / (1-this.cc_control_fairnessFactor_2.getDouble());
422                                counter_d2 ++; 
423                        } 
424                }
425                if ((counter_d1 != D1) || (counter_d2 != D2)) throw new RuntimeException ("Bad");
426                return Triple.of(obj, Pair.of(h1, h2),Pair.of(y1, y2));
427        }
428
429
430        private double computeHdFromPrices (Demand d)
431        {
432                DoubleMatrix1D infoIKnow_price_e = this.control_mostUpdatedLinkPriceKnownDemand_de.viewRow(d.getIndex ());
433
434                /* compute the demand price as weighted sum in the routes of route prices  */
435                double demandWeightedSumLinkPrices = 0;
436                double demandCarriedTraffic = 0; 
437                for (Route r : d.getRoutes ())
438                {
439                        final double h_r = r.getCarriedTraffic();
440                        demandCarriedTraffic += h_r;
441                        for (Link e : r.getSeqLinksRealPath())
442                                demandWeightedSumLinkPrices += h_r * infoIKnow_price_e.get(e.getIndex ());
443                }
444                //if (Math.abs(demandCarriedTraffic - this.currentNetPlan.getDemandCarriedTraffic(dIdMe)) > 1E-3) throw new RuntimeException ("Not all the traffic is carried");
445                demandWeightedSumLinkPrices /= demandCarriedTraffic;
446
447                /* compute the new h_d */
448                final double alpha = (demandType.get(d.getIndex ()) == 1)? this.cc_control_fairnessFactor_1.getDouble() : this.cc_control_fairnessFactor_2.getDouble();
449                final double weight = (demandType.get(d.getIndex ()) == 1)? this.cc_control_weightFairness_1.getDouble() : this.cc_control_weightFairness_2.getDouble();
450                final double new_hd = Math.max(this.cc_control_minHd.getDouble() , Math.min(this.cc_control_maxHd.getDouble(), Math.pow(demandWeightedSumLinkPrices/weight, -1/alpha)));
451
452                return new_hd;
453        }
454
455}