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_wirelessCongControlTransmissionPowerAssignment;
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.libraries.NetworkPerformanceMetrics;
035import com.net2plan.libraries.WirelessUtils;
036import com.net2plan.utils.Constants.RoutingType;
037import com.net2plan.utils.GradientProjectionUtils;
038import com.net2plan.utils.InputParameter;
039import com.net2plan.utils.Pair;
040import com.net2plan.utils.TimeTrace;
041import com.net2plan.utils.Triple;
042
043/** 
044 * This module implements a distributed dual-decomposition-based gradient algorithm, for a coordinated adjustment of the traffic to inject by each demand (congestion control), and the transmission power in each link of the underlying wireless network, to maximize the network utility enforcing a fair allocation of the resources.
045 *
046 * Ths event processor is adapted to permit observing the algorithm performances under user-defined conditions, 
047 * including asynchronous distributed executions, where signaling can be affected by losses and/or delays, and/or measurement errors. 
048 * The time evolution of different metrics can be stored in output files, for later processing. 
049 * As an example, see the <a href="../../../../../../graphGeneratorFiles/fig_sec11_4_congControlAndBackpressureRouting_dualDecomp.m">{@code fig_sec11_4_congControlAndBackpressureRouting_dualDecomp.m}</a> MATLAB file used for generating the graph/s of the case study in the 
050 * <a href="http://eu.wiley.com/WileyCDA/WileyTitle/productCd-1119013356.html">book</a> using this algorithm.
051 * 
052 * To simulate a network with this module, use the {@code Online_evGen_doNothing} generator.
053 * 
054 * @net2plan.keywords Bandwidth assignment (BA), Transmission power optimization, Wireless, Capacity assignment (CA), Distributed algorithm, Dual decomposition 
055 * @net2plan.ocnbooksections Section 11.4
056 * @net2plan.inputParameters 
057 * @author Pablo Pavon-Marino
058 */
059public class Online_evProc_congControlAndTransmissionPowerAssignmentDualDecomp extends IEventProcessor 
060{
061        private Random rng;
062
063        static final int MAC_SIGNALING_WAKEUPTOSENDMESSAGE = 200;
064        static final int MAC_SIGNALING_RECEIVEDMESSAGE = 201;
065        static final int MAC_UPDATE_WAKEUPTOUPDATE = 202;
066
067        static final int CC_SIGNALING_WAKEUPTOSENDMESSAGE = 400;
068        static final int CC_SIGNALING_RECEIVEDMESSAGE = 401;
069        static final int CC_UPDATE_WAKEUPTOUPDATE = 402;
070
071        
072        private InputParameter mac_signaling_isSynchronous = new InputParameter ("mac_signaling_isSynchronous", false , "true if all the distributed agents involved wake up synchronously to send the signaling messages");
073        private InputParameter mac_signaling_averageInterMessageTime = new InputParameter ("mac_signaling_averageInterMessageTime", 1.0 , "Average time between two signaling messages sent by an agent" , 0 , false , Double.MAX_VALUE , true);
074        private InputParameter mac_signaling_maxFluctuationInterMessageTime = new InputParameter ("mac_signaling_maxFluctuationInterMessageTime", 0.5 , "Max fluctuation in time between two signaling messages sent by an agent" , 0 , true , Double.MAX_VALUE , true);
075        private InputParameter mac_signaling_averageDelay = new InputParameter ("mac_signaling_averageDelay", 0.0 , "Average time between signaling message transmission by an agent and its reception by other or others" , 0 , true , Double.MAX_VALUE , true);
076        private InputParameter mac_signaling_maxFluctuationInDelay = new InputParameter ("mac_signaling_maxFluctuationInDelay", 0.0 , "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);
077        private InputParameter mac_signaling_signalingLossProbability = new InputParameter ("mac_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);
078
079        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");
080        private InputParameter mac_update_averageInterUpdateTime = new InputParameter ("mac_update_averageInterUpdateTime", 1.0 , "Average time between two updates of an agent" , 0 , false , Double.MAX_VALUE , true);
081        private InputParameter mac_update_maxFluctuationInterUpdateTime = new InputParameter ("mac_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);
082
083        private InputParameter mac_gradient_maxGradientAbsoluteNoise = new InputParameter ("mac_gradient_maxGradientAbsoluteNoise", 0.01 , "Max value of the added noise to the gradient coordinate in absolute values" , 0 , true , Double.MAX_VALUE , true);
084        private InputParameter mac_gradient_gammaStep = new InputParameter ("mac_gradient_gammaStep", 1.0 , "Gamma step in the gradient algorithm" , 0 , false , Double.MAX_VALUE , true);
085        private InputParameter mac_gradient_heavyBallBetaParameter = new InputParameter ("mac_gradient_heavyBallBetaParameter", 0.0 , "Beta parameter of heavy ball, between 0 and 1. Value 0 means no heavy ball" , 0 , true , 1.0 , true);
086        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);
087
088        private InputParameter mac_pathLossExponent = new InputParameter ("mac_pathLossExponent", 3.0 , "Exponent in the model for propagation losses" , 0 , true , Double.MAX_VALUE , true);
089        private InputParameter mac_worseRiseOverThermal_nu = new InputParameter ("mac_worseRiseOverThermal_nu", 10.0 , "Worse case ratio between interference power and thermal power at any receiver. Used to set the common thermal power in the receivers" , 0 , true , Double.MAX_VALUE , true);
090        private InputParameter mac_interferenceAttenuationFactor_nu = new InputParameter ("mac_interferenceAttenuationFactor_nu", 1.0e6 , "The interference power received in natural units is divided by this to reduce its effect" , 1 , true , Double.MAX_VALUE , true);
091        private InputParameter mac_maxTransmissionPower_logu = new InputParameter ("mac_maxTransmissionPower_logu", 3.0 , "The maximum link transmission power in logarithmic units (e.g. dBm)");
092        private InputParameter mac_minTransmissionPower_logu = new InputParameter ("mac_minTransmissionPower_logu", 0.0 , "The minimum link transmission power in logarithmic units (e.g. dBm)");
093        private InputParameter mac_simulation_maxNumberOfUpdateIntervals = new InputParameter ("mac_simulation_maxNumberOfUpdateIntervals", 2000.0 , "Maximum number of update intervals in average per agent" , 0 , false , Double.MAX_VALUE , true);
094
095        private InputParameter simulation_randomSeed = new InputParameter ("simulation_randomSeed", (long) 1 , "Seed of the random number generator");
096        private InputParameter simulation_outFileNameRoot = new InputParameter ("simulation_outFileNameRoot", "crossLayerCongControlPowerTransmissionAsignmentDualDecomp" , "Root of the file name to be used in the output files. If blank, no output");
097
098        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");
099        private InputParameter cc_signaling_averageInterMessageTime = new InputParameter ("cc_signaling_averageInterMessageTime", 10.0 , "Average time between two signaling messages sent by an agent" , 0 , false , Double.MAX_VALUE , true);
100        private InputParameter cc_signaling_maxFluctuationInterMessageTime = new InputParameter ("cc_signaling_maxFluctuationInterMessageTime", 5.0 , "Max fluctuation in time between two signaling messages sent by an agent" , 0 , true , Double.MAX_VALUE , true);
101        private InputParameter cc_signaling_averageDelay = new InputParameter ("cc_signaling_averageDelay", 30.0 , "Average time between signaling message transmission by an agent and its reception by other or others" , 0 , true , Double.MAX_VALUE , true);
102        private InputParameter cc_signaling_maxFluctuationInDelay = new InputParameter ("cc_signaling_maxFluctuationInDelay", 15.0 , "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);
103        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);
104        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");
105        private InputParameter cc_update_averageInterUpdateTime = new InputParameter ("cc_update_averageInterUpdateTime", 10.0 , "Average time between two updates of an agent" , 0 , false , Double.MAX_VALUE , true);
106        private InputParameter cc_update_maxFluctuationInterUpdateTime = new InputParameter ("cc_update_maxFluctuationInterUpdateTime", 5.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);
107
108        private InputParameter cc_gradient_gammaStep = new InputParameter ("cc_gradient_gammaStep", 0.001 , "Gamma step in the gradient algorithm" , 0 , false , Double.MAX_VALUE , true);
109        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);
110
111        private InputParameter cc_minHd = new InputParameter ("cc_minHd", 0.1 , "Minimum traffic assigned to each demand" , 0 , true , Double.MAX_VALUE , true);
112        private InputParameter cc_maxHd = new InputParameter ("cc_maxHd", 1.0E6 , "Maximum traffic assigned to each demand" , 0 , true , Double.MAX_VALUE , true);
113        private InputParameter cc_fairnessFactor = new InputParameter ("cc_fairnessFactor", 1.0 , "Fairness factor in utility function of congestion control" , 0 , true , Double.MAX_VALUE , true);
114        private InputParameter cc_initialLinkPrices = new InputParameter ("cc_initialLinkPrices", 1 , "Initial value of the link weights" , 0 , true , Double.MAX_VALUE , true);
115        private InputParameter cc_simulation_maxNumberOfUpdateIntervals = new InputParameter ("cc_simulation_maxNumberOfUpdateIntervals", 2000.0 , "Maximum number of update intervals in average per agent" , 0 , false , Double.MAX_VALUE , true);
116
117        private Map<String, String> algorithmParameters;
118        private Map<String, String> net2planParameters;
119
120        private DoubleMatrix1D cc_price_e;
121        private DoubleMatrix2D cc_mostUpdatedLinkPriceKnownByDemand_de; 
122        
123        private DoubleMatrix2D mac_g_nu_ee;
124        private DoubleMatrix1D mac_transmissionPower_logu_e; 
125        private DoubleMatrix1D mac_previousTransmissionPower_logu_e; 
126        private double mac_receptionThermalNoise_nu;
127        private DoubleMatrix2D mac_mostUpdatedMe2ValueKnownByLink1_e1e2; 
128                        
129        private TimeTrace traceOf_pi_e;
130        private TimeTrace traceOf_u_e;
131        private TimeTrace traceOf_y_e;
132        private TimeTrace traceOf_p_e;
133        private TimeTrace traceOf_h_d;
134        private TimeTrace traceOf_objFunction;
135
136        private NetPlan currentNetPlan;
137        private NetPlan copyInitialNetPlan;
138        private int E,N,D;
139        
140        @Override
141        public String getDescription()
142        {
143                return "This module implements a distributed dual-decomposition-based gradient algorithm, for a coordinated adjustment of the traffic to inject by each demand (congestion control), and the transmission power in each link of the underlying wireless network, to maximize the network utility enforcing a fair allocation of the resources.";
144        }
145
146        @Override
147        public List<Triple<String, String, String>> getParameters()
148        {
149                /* Returns the parameter information for all the InputParameter objects defined in this object (uses Java reflection) */
150                return InputParameter.getInformationAllInputParameterFieldsOfObject(this);
151        }
152
153        @Override
154        public void initialize(NetPlan currentNp, Map<String, String> algorithmParameters, Map<String, String> simulationParameters, Map<String, String> net2planParameters)
155        {
156                /* Initialize all InputParameter objects defined in this object (this uses Java reflection) */
157                InputParameter.initializeAllInputParameterFieldsOfObject(this, algorithmParameters);
158
159                this.currentNetPlan = currentNp;
160                this.copyInitialNetPlan = currentNp.copy();
161                this.algorithmParameters = algorithmParameters;
162                this.net2planParameters = net2planParameters;
163                this.N = this.currentNetPlan.getNumberOfNodes();
164                this.E = this.currentNetPlan.getNumberOfLinks();
165                if (currentNp.getNumberOfLayers() != 1) throw new Net2PlanException ("This algorithm works in single layer networks");
166                if (E == 0) throw new Net2PlanException ("The input design should have links");
167
168                this.rng = new Random (simulation_randomSeed.getLong ());
169                
170                /* INITIALIZE INFORMATION COMMON TO MAC AND CONGESTION CONTROL */
171                this.cc_price_e = DoubleFactory1D.dense.make (E , cc_initialLinkPrices.getDouble());
172                
173                /* Remove all routes, and create one with the shortest path in km for each demand */
174                this.currentNetPlan.removeAllDemands();
175                currentNetPlan.setRoutingType(RoutingType.SOURCE_ROUTING);
176                for (Node n1 : this.currentNetPlan.getNodes())
177                        for (Node n2 : this.currentNetPlan.getNodes())
178                                if (n1 != n2)
179                                        this.currentNetPlan.addDemand(n1, n2, cc_minHd.getDouble(), null);
180                currentNetPlan.addRoutesFromCandidatePathList(currentNetPlan.getVectorLinkLengthInKm().toArray()  , "K" , "1");
181                for (Route r : currentNetPlan.getRoutes ()) r.setCarriedTraffic(r.getDemand().getOfferedTraffic() , r.getDemand().getOfferedTraffic());
182                this.D = this.currentNetPlan.getNumberOfDemands();
183                
184                /* Initialize the gains between links, normalizing them so that the maximum gain is one */
185                this.mac_g_nu_ee = WirelessUtils.computeInterferenceMatrixNaturalUnits (currentNetPlan.getLinks () , mac_interferenceAttenuationFactor_nu.getDouble() , mac_pathLossExponent.getDouble());
186                final double maxGainValue = mac_g_nu_ee.getMaxLocation() [0];
187                mac_g_nu_ee.assign (DoubleFunctions.div(maxGainValue));
188                
189                /* Initialize the thermal noise at the receivers, to have a worse case ROT (rise over thermal) */
190                double worseInterferenceReceivedAtMaxPower_nu = WirelessUtils.computeWorseReceiverInterferencePower_nu (mac_maxTransmissionPower_logu.getDouble()  , mac_g_nu_ee);
191
192                /* Adjust the thermal noise in the receivers so that we have a given ROT */
193                this.mac_receptionThermalNoise_nu = worseInterferenceReceivedAtMaxPower_nu / mac_worseRiseOverThermal_nu.getDouble();
194
195                /* Initialize the transmission power variables */
196                this.mac_transmissionPower_logu_e = DoubleFactory1D.dense.make (E , mac_minTransmissionPower_logu.getDouble());
197                this.mac_previousTransmissionPower_logu_e = DoubleFactory1D.dense.make (E , mac_minTransmissionPower_logu.getDouble());
198                        
199                /* Update the netplan object with the resulting capacities */
200                for (Link e : currentNp.getLinks())
201                        e.setCapacity(Math.log(computeSINR_e (e)));
202
203                this.mac_mostUpdatedMe2ValueKnownByLink1_e1e2 = DoubleFactory2D.dense.make (E,E);
204                for (Link e1 : currentNp.getLinks())
205                {
206                        for (Link otherLink : currentNp.getLinks())
207                                if (e1 != otherLink)
208                                        mac_mostUpdatedMe2ValueKnownByLink1_e1e2.set(e1.getIndex () , otherLink.getIndex(), computeMeFactor_e (otherLink));
209                }
210                
211                /* Initially all nodes receive a "wake up to transmit" event, aligned at time zero or y asynchr => randomly chosen */
212                for (Link e : currentNp.getLinks())
213                {
214                        final double signalingTime = (mac_signaling_isSynchronous.getBoolean())? mac_signaling_averageInterMessageTime.getDouble() : Math.max(0 , mac_signaling_averageInterMessageTime.getDouble() + mac_signaling_maxFluctuationInterMessageTime.getDouble() * (rng.nextDouble() - 0.5));
215                        this.scheduleEvent(new SimEvent (signalingTime , SimEvent.DestinationModule.EVENT_PROCESSOR , MAC_SIGNALING_WAKEUPTOSENDMESSAGE , e));
216                        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));
217                        this.scheduleEvent(new SimEvent (updateTime , SimEvent.DestinationModule.EVENT_PROCESSOR , MAC_UPDATE_WAKEUPTOUPDATE , e));
218                }
219                
220                /* INITIALIZATION CONGESTION CONTROL LAYER */
221
222                /* Initialize the information each demand knows of the prices of all the links */
223                this.cc_mostUpdatedLinkPriceKnownByDemand_de = DoubleFactory2D.dense.make (D,E,cc_initialLinkPrices.getDouble());
224                
225                /* Compute the traffic each demand injects, set it in net2plan, and tell routing layer there was an update in this */
226                for (Demand d : currentNetPlan.getDemands())
227                {
228                        final double h_d = this.computeHdFromPrices(d);
229                        d.setOfferedTraffic(h_d);
230                        d.getRoutes ().iterator().next().setCarriedTraffic(h_d , h_d);
231                }
232                
233                /* Initially all nodes receive a "wake up to transmit" event, aligned at time zero or y asynchr => randomly chosen */
234                for (Link e : currentNetPlan.getLinks())
235                {
236                        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));
237                        this.scheduleEvent(new SimEvent (signalingTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_SIGNALING_WAKEUPTOSENDMESSAGE , e));
238                }
239                for (Demand d : currentNetPlan.getDemands())
240                {
241                        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));
242                        this.scheduleEvent(new SimEvent (updateTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_UPDATE_WAKEUPTOUPDATE , d));
243                }
244
245                /* INITIALIZE STATISTIC VARIABLES FOR BOTH LAYERS */
246                this.traceOf_u_e = new TimeTrace  ();
247                this.traceOf_y_e = new TimeTrace  ();
248                this.traceOf_p_e = new TimeTrace ();
249                this.traceOf_pi_e = new TimeTrace ();
250                this.traceOf_h_d = new TimeTrace (); 
251                this.traceOf_objFunction = new TimeTrace ();
252
253                this.traceOf_u_e.add(0.0, this.currentNetPlan.getVectorLinkCapacity());
254                this.traceOf_y_e.add(0.0, this.currentNetPlan.getVectorLinkTotalCarriedTraffic());
255                this.traceOf_p_e.add(0.0, this.mac_transmissionPower_logu_e.copy ());
256                this.traceOf_pi_e.add(0.0, this.cc_price_e.copy ());
257                this.traceOf_h_d.add(0.0, this.currentNetPlan.getVectorDemandCarriedTraffic());
258                this.traceOf_objFunction.add(0.0 , NetworkPerformanceMetrics.alphaUtility(this.currentNetPlan.getVectorDemandOfferedTraffic() , cc_fairnessFactor.getDouble()));
259                
260        }
261
262        @Override
263        public void processEvent(NetPlan currentNetPlan, SimEvent event)
264        {
265                final double t = event.getEventTime();
266                switch (event.getEventType())
267                {
268                case MAC_SIGNALING_RECEIVEDMESSAGE: // A node broadcasts signaling info to its 1 hop neighbors
269                {
270                        final Pair<Link,Pair<Link,Double>> signalInfo = (Pair<Link,Pair<Link,Double>>) event.getEventObject();
271                        final Link eMe = signalInfo.getFirst();
272                        final Link otherLink = signalInfo.getSecond().getFirst();
273                        final double otherLink_me = signalInfo.getSecond().getSecond();
274                        mac_mostUpdatedMe2ValueKnownByLink1_e1e2.set (eMe.getIndex () , otherLink.getIndex () , otherLink_me);
275                        break;
276                }
277                
278                case MAC_SIGNALING_WAKEUPTOSENDMESSAGE: // A node broadcasts signaling info to its 1 hop neighbors
279                {
280                        final Link eMe = (Link) event.getEventObject();
281
282                        Pair<Link,Double> infoToSignal = Pair.of(eMe , this.computeMeFactor_e(eMe));
283                        if (rng.nextDouble() >= this.mac_signaling_signalingLossProbability.getDouble()) // the signaling may be lost => lost to all nodes
284                                for (Link otherLink : currentNetPlan.getLinks ())
285                                        if (otherLink != eMe)
286                                        {
287                                                final double signalingReceptionTime = t + Math.max(0 , mac_signaling_averageDelay.getDouble() + mac_signaling_maxFluctuationInDelay.getDouble() * (rng.nextDouble() - 0.5));
288                                                this.scheduleEvent(new SimEvent (signalingReceptionTime , SimEvent.DestinationModule.EVENT_PROCESSOR , MAC_SIGNALING_RECEIVEDMESSAGE , Pair.of(otherLink , infoToSignal)));
289                                        }
290
291                        /* Re-schedule when to wake up again */
292                        final double signalingTime = mac_signaling_isSynchronous.getBoolean()? t + mac_signaling_averageInterMessageTime.getDouble() : Math.max(t , t + mac_signaling_averageInterMessageTime.getDouble() + mac_signaling_maxFluctuationInterMessageTime.getDouble() * (rng.nextDouble() - 0.5));
293                        this.scheduleEvent(new SimEvent (signalingTime , SimEvent.DestinationModule.EVENT_PROCESSOR , MAC_SIGNALING_WAKEUPTOSENDMESSAGE , eMe));
294                        break;
295                }
296
297                case MAC_UPDATE_WAKEUPTOUPDATE: 
298                {
299                        final Link eMe = (Link) event.getEventObject(); 
300                        
301                        final double currentTransmissionPower_logu = this.mac_transmissionPower_logu_e.get(eMe.getIndex ());
302                        double gradientThisLink = computeGradient (eMe) + 2*mac_gradient_maxGradientAbsoluteNoise.getDouble()*(rng.nextDouble()-0.5);
303                        double nextTransmissionPower_logu = currentTransmissionPower_logu + this.mac_gradient_gammaStep.getDouble() * gradientThisLink;
304
305                        /* heavy ball */
306                        nextTransmissionPower_logu += this.mac_gradient_heavyBallBetaParameter.getDouble() * (currentTransmissionPower_logu - mac_previousTransmissionPower_logu_e.get(eMe.getIndex ()));
307                        /* projection */
308                        nextTransmissionPower_logu = GradientProjectionUtils.euclideanProjection_boxLike(nextTransmissionPower_logu, this.mac_minTransmissionPower_logu.getDouble(), this.mac_maxTransmissionPower_logu.getDouble());
309
310                        if (mac_gradient_maxGradientCoordinateChange.getDouble() > 0)
311                                nextTransmissionPower_logu = GradientProjectionUtils.scaleDown_maxAbsoluteCoordinateChange (currentTransmissionPower_logu , nextTransmissionPower_logu ,  mac_gradient_maxGradientCoordinateChange.getDouble());
312                        
313                        this.mac_previousTransmissionPower_logu_e.set(eMe.getIndex (), this.mac_transmissionPower_logu_e.get(eMe.getIndex ()));
314                        this.mac_transmissionPower_logu_e.set(eMe.getIndex (), nextTransmissionPower_logu);
315
316                        /* Send event next recomputing time */
317                        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));
318                        this.scheduleEvent(new SimEvent (updateTime , SimEvent.DestinationModule.EVENT_PROCESSOR , MAC_UPDATE_WAKEUPTOUPDATE,  eMe));
319
320                        /* Update in currentNetPlan the capacity values */
321                        for (Link e : this.currentNetPlan.getLinks())
322                                e.setCapacity(Math.log(computeSINR_e (e)));
323
324                        this.traceOf_u_e.add(t, this.currentNetPlan.getVectorLinkCapacity());
325                        this.traceOf_y_e.add(t, this.currentNetPlan.getVectorLinkTotalCarriedTraffic());
326                        this.traceOf_p_e.add(t, this.mac_transmissionPower_logu_e.copy ());
327                        this.traceOf_pi_e.add(t, this.cc_price_e.copy ());
328                        this.traceOf_h_d.add(t, this.currentNetPlan.getVectorDemandCarriedTraffic());
329                        this.traceOf_objFunction.add(t , NetworkPerformanceMetrics.alphaUtility(this.currentNetPlan.getVectorDemandOfferedTraffic() , cc_fairnessFactor.getDouble()));
330
331                        if (t > this.mac_simulation_maxNumberOfUpdateIntervals.getDouble() * this.mac_update_averageInterUpdateTime.getDouble()) { this.endSimulation (); }
332                        break;
333                }
334
335                case CC_SIGNALING_RECEIVEDMESSAGE: // A node receives from an out neighbor the q_nt for any destination
336                {
337                        final Pair<Demand,Pair<Link,Double>> signalInfo = (Pair<Demand,Pair<Link,Double>>) event.getEventObject();
338                        final Demand d = signalInfo.getFirst();
339                        final Pair<Link,Double> receivedInfo_price_e = signalInfo.getSecond();
340                        this.cc_mostUpdatedLinkPriceKnownByDemand_de.set(d.getIndex() , receivedInfo_price_e.getFirst().getIndex () ,receivedInfo_price_e.getSecond());  
341                        break;
342                }
343                
344                case CC_SIGNALING_WAKEUPTOSENDMESSAGE: // A node broadcasts signaling info to its 1 hop neighbors
345                {
346                        final Link eMe = (Link) event.getEventObject();
347
348                        /* Update the new price with the gradient approach */
349                        final double u_e = eMe.getCapacity();
350                        final double y_e = eMe.getCarriedTrafficIncludingProtectionSegments();
351                        final double old_pie = this.cc_price_e.get(eMe.getIndex());
352                        final double new_pie = Math.max(0, old_pie - this.cc_gradient_gammaStep.getDouble() * (u_e - y_e) + 2*cc_gradient_maxGradientAbsoluteNoise.getDouble()*(rng.nextDouble()-0.5));
353                        this.cc_price_e.set(eMe.getIndex(), new_pie);
354
355                        /* Create the info I will signal */
356                        Pair<Link,Double> infoToSignal = Pair.of(eMe ,  new_pie);
357
358                        /* Send the events of the signaling information messages to all the nodes */
359                        for (Route route : eMe.getTraversingRoutes())
360                        {
361                                if (rng.nextDouble() < this.cc_signaling_signalingLossProbability.getDouble()) continue; // the signaling may be lost => lost to all demands
362                                final Demand d = route.getDemand();
363                                final double signalingReceptionTime = t + Math.max(0 , cc_signaling_averageDelay.getDouble() + cc_signaling_maxFluctuationInDelay.getDouble() * (rng.nextDouble() - 0.5));
364                                this.scheduleEvent(new SimEvent (signalingReceptionTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_SIGNALING_RECEIVEDMESSAGE , Pair.of(d , infoToSignal)));
365                        }
366
367                        /* Re-schedule when to wake up again */
368                        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));
369                        this.scheduleEvent(new SimEvent (signalingTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_SIGNALING_WAKEUPTOSENDMESSAGE , eMe));
370                        break;
371                }
372
373                case CC_UPDATE_WAKEUPTOUPDATE: // a node updates its p_n, p_e values, using most updated info available
374                {
375                        final Demand dMe = (Demand) event.getEventObject();
376                        
377                        /* compute the new h_d and apply it */
378                        final double new_hd = computeHdFromPrices (dMe);
379                        dMe.setOfferedTraffic(new_hd);
380                        dMe.getRoutes ().iterator().next().setCarriedTraffic(new_hd , new_hd);
381
382                        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));
383                        this.scheduleEvent(new SimEvent (updateTime , SimEvent.DestinationModule.EVENT_PROCESSOR , CC_UPDATE_WAKEUPTOUPDATE,  dMe));
384
385                        this.traceOf_objFunction.add(t, NetworkPerformanceMetrics.alphaUtility(this.currentNetPlan.getVectorDemandOfferedTraffic() , cc_fairnessFactor.getDouble()));
386
387                        if (t > this.cc_simulation_maxNumberOfUpdateIntervals.getDouble() * this.cc_update_averageInterUpdateTime.getDouble()) { this.endSimulation (); }
388                        break;
389                }
390                
391                
392                default: throw new RuntimeException ("Unexpected received event");
393                }
394                
395                
396        }
397
398        public String finish (StringBuilder st , double simTime)
399        {
400                if (simulation_outFileNameRoot.getString().equals("")) return null;
401                this.traceOf_objFunction.printToFile(new File (simulation_outFileNameRoot.getString() + "_objFunc.txt"));
402                this.traceOf_u_e.printToFile(new File (simulation_outFileNameRoot.getString() + "_ue.txt"));
403                this.traceOf_y_e.printToFile(new File (simulation_outFileNameRoot.getString() + "_ye.txt"));
404                this.traceOf_p_e.printToFile(new File (simulation_outFileNameRoot.getString() + "_pe.txt"));
405                this.traceOf_pi_e.printToFile(new File (simulation_outFileNameRoot.getString() + "_pie.txt"));
406                this.traceOf_h_d.printToFile(new File (simulation_outFileNameRoot.getString() + "_hd.txt"));
407
408                Map<String,String> param = new HashMap<String,String> (algorithmParameters);
409                param.put("solverName", "ipopt");
410                param.put("solverLibraryName", "");
411                param.put("maxSolverTimeInSeconds", "-1");
412                new Offline_cba_wirelessCongControlTransmissionPowerAssignment ().executeAlgorithm(copyInitialNetPlan , param , this.net2planParameters);
413                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_objFunc.txt") , NetworkPerformanceMetrics.alphaUtility(copyInitialNetPlan.getVectorDemandOfferedTraffic() , cc_fairnessFactor.getDouble()));
414                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_ue.txt") , copyInitialNetPlan.getVectorLinkCapacity());
415                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_ye.txt") , copyInitialNetPlan.getVectorLinkTotalCarriedTraffic());
416                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_hd.txt") , copyInitialNetPlan.getVectorDemandCarriedTraffic());
417                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_pe.txt") , copyInitialNetPlan.getVectorAttributeValues(copyInitialNetPlan.getLinks () , "p_e"));
418                TimeTrace.printToFile(new File (simulation_outFileNameRoot.getString() + "_jom_pie.txt") , copyInitialNetPlan.getVectorAttributeValues(copyInitialNetPlan.getLinks () , "pi_e"));
419                return null;
420        }
421        
422
423        private double computeGradient (Link thisLink)
424        {
425                final DoubleMatrix1D infoIKnow_me = mac_mostUpdatedMe2ValueKnownByLink1_e1e2.viewRow (thisLink.getIndex ());
426                
427                double gradient = cc_price_e.get(thisLink.getIndex ()); 
428                double accumFactor = 0;
429                for (Link epp : this.currentNetPlan.getLinks())
430                        if (epp != thisLink)
431                                accumFactor += this.mac_g_nu_ee.get(thisLink.getIndex (),epp.getIndex ()) * infoIKnow_me.get(epp.getIndex ());
432                accumFactor *= Math.exp(this.mac_transmissionPower_logu_e.get(thisLink.getIndex ()));
433                gradient -= accumFactor;
434                return gradient;
435        }
436
437        private double computeMeFactor_e (Link e)
438        {
439                final double snr_e = Math.exp(e.getCapacity());
440                final double pi_e = this.cc_price_e.get(e.getIndex ());
441                final double g_ee = this.mac_g_nu_ee.get(e.getIndex (),e.getIndex ());
442                return  pi_e * snr_e / (Math.exp(this.mac_transmissionPower_logu_e.get(e.getIndex ())) * g_ee);
443        }
444        
445        private double computeSINR_e (Link e)
446        {
447                final double receivedPower_nu = Math.exp(this.mac_transmissionPower_logu_e.get(e.getIndex())) * this.mac_g_nu_ee.get(e.getIndex (),e.getIndex ());
448                double interferencePower_nu = this.mac_receptionThermalNoise_nu;
449                for (Link eInt : this.currentNetPlan.getLinks ())
450                        if (eInt != e) interferencePower_nu += Math.exp(this.mac_transmissionPower_logu_e.get(eInt.getIndex())) * this.mac_g_nu_ee.get(eInt.getIndex (),e.getIndex ());
451                return receivedPower_nu / interferencePower_nu;
452        }
453
454        private double computeHdFromPrices (Demand d)
455        {
456                DoubleMatrix1D infoIKnow_price_e = this.cc_mostUpdatedLinkPriceKnownByDemand_de.viewRow(d.getIndex());
457
458                /* compute the demand price as weighted sum in the routes of route prices  */
459                double demandWeightedSumLinkPrices = 0;
460                double demandCarriedTraffic = 0; 
461                for (Route r : d.getRoutes ())
462                {
463                        final double h_r = r.getCarriedTraffic();
464                        demandCarriedTraffic += h_r;
465                        for (Link e : r.getSeqLinksRealPath())
466                                demandWeightedSumLinkPrices += h_r * infoIKnow_price_e.get(e.getIndex());
467                }
468                //if (Math.abs(demandCarriedTraffic - this.currentNetPlan.getDemandCarriedTraffic(dIdMe)) > 1E-3) throw new RuntimeException ("Not all the traffic is carried");
469                demandWeightedSumLinkPrices /= demandCarriedTraffic;
470
471                /* compute the new h_d */
472                final double new_hd = Math.max(this.cc_minHd.getDouble() , Math.min(this.cc_maxHd.getDouble(), Math.pow(demandWeightedSumLinkPrices, -1/this.cc_fairnessFactor.getDouble())));
473
474                return new_hd;
475        }
476
477}