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.io.FileOutputStream; 015import java.io.PrintStream; 016import java.util.HashMap; 017import java.util.HashSet; 018import java.util.List; 019import java.util.Map; 020import java.util.Random; 021import java.util.Set; 022 023import cern.colt.matrix.tdouble.DoubleFactory1D; 024import cern.colt.matrix.tdouble.DoubleFactory2D; 025import cern.colt.matrix.tdouble.DoubleMatrix1D; 026import cern.colt.matrix.tdouble.DoubleMatrix2D; 027import cern.jet.math.tdouble.DoublePlusMultFirst; 028import cern.jet.math.tdouble.DoublePlusMultSecond; 029 030import com.jom.DoubleMatrixND; 031import com.jom.OptimizationProblem; 032import com.net2plan.examples.ocnbook.offline.Offline_fa_xteFormulations; 033import com.net2plan.interfaces.networkDesign.Demand; 034import com.net2plan.interfaces.networkDesign.Link; 035import com.net2plan.interfaces.networkDesign.Net2PlanException; 036import com.net2plan.interfaces.networkDesign.NetPlan; 037import com.net2plan.interfaces.networkDesign.Node; 038import com.net2plan.interfaces.simulation.IEventProcessor; 039import com.net2plan.interfaces.simulation.SimEvent; 040import com.net2plan.libraries.GraphUtils.ClosedCycleRoutingException; 041import com.net2plan.utils.Constants.RoutingType; 042import com.net2plan.utils.Constants.SearchType; 043import com.net2plan.utils.DoubleUtils; 044import com.net2plan.utils.InputParameter; 045import com.net2plan.utils.IntUtils; 046import com.net2plan.utils.TimeTrace; 047import com.net2plan.utils.Triple; 048 049/** 050 * This module implements a distributed primal-decomposition-based gradient algorithm, for a coordinated adjustment of the routing in multiple domains (or cluster, or autonomous systems) in a network, so that domains do not need to exchange sensitive internal information, and minimize the average number of hops in the network. 051 * 052 * Ths event processor is adapted to permit observing the algorithm performances under user-defined conditions, 053 * including asynchronous distributed executions, where signaling can be affected by losses and/or delays, and/or measurement errors. 054 * The time evolution of different metrics can be stored in output files, for later processing. 055 * As an example, see the <a href="../../../../../../graphGeneratorFiles/fig_sec11_6_multidomainRouting_primalDecomp.m">{@code fig_sec11_6_multidomainRouting_primalDecomp.m}</a> MATLAB file used for generating the graph/s of the case study in the 056 * <a href="http://eu.wiley.com/WileyCDA/WileyTitle/productCd-1119013356.html">book</a> using this algorithm. 057 * 058 * To simulate a network with this module, use the {@code Online_evGen_doNothing} generator. 059 * 060 * @net2plan.keywords Multidomain network, Primal decomposition, Distributed algorithm, Flow assignment (FA), Destination-based routing, Destination-link formulation 061 * @net2plan.ocnbooksections Section 11.6 062 * @net2plan.inputParameters 063 * @author Pablo Pavon-Marino 064 */ 065public class Online_evProc_multidomainRoutingPrimalDecomp extends IEventProcessor 066{ 067 private InputParameter routing_isDomainRoutingRecomputationSynchronous = new InputParameter("routing_isDomainRoutingRecomputationSynchronous", false, "true if all domains wake up synchronously to recompute its routing"); 068 private InputParameter routing_averageDomainRoutingRecomputationInterval = new InputParameter("routing_averageDomainRoutingRecomputationInterval", 1.0, "Average time between two routing recomputation of a domain", 0, false, Double.MAX_VALUE, true); 069 private InputParameter routing_maxFluctuationDomainRoutingRecomputationInterval = new InputParameter("routing_maxFluctuationDomainRoutingRecomputationInterval", 0.5, "Max fluctuation in time between two routing recomputation of a domain", 0, true, Double.MAX_VALUE, true); 070 private InputParameter routing_averageMasterUpdateInterval = new InputParameter("routing_averageMasterUpdateInterval", 1.0, "Average time between two master recomputation", 0, false, Double.MAX_VALUE, true); 071 private InputParameter routing_maxFluctuationMasterUpdateInterval = new InputParameter("routing_maxFluctuationMasterUpdateInterval", 0.5, "Max fluctuation in time between two master updates", 0, true, Double.MAX_VALUE, true); 072 private InputParameter gradient_gammaStep = new InputParameter("gradient_gammaStep", 0.1, "Gamma step in the gradient algorithm", 0, false, Double.MAX_VALUE, true); 073 private InputParameter simulation_maxNumberOfMasterUpdateIntervals = new InputParameter("simulation_maxNumberOfMasterUpdateIntervals", 100.0, "Maximum number of update intervals in average per agent", 0, false, Double.MAX_VALUE, true); 074 private InputParameter simulation_randomSeed = new InputParameter("simulation_randomSeed", (long) 1, "Seed of the random number generator"); 075 private InputParameter simulation_outFileNameRoot = new InputParameter("simulation_outFileNameRoot", "multidomainRoutingPrimalDecompostion", "Root of the file name to be used in the output files. If blank, no output"); 076 077 private int N, E, C, Ef, D; // total number of nodes and links in the network. Ef is the number of frontier links 078 private int[] allT; // array of constants 0...N-1 079 private int[] allE; // array of constants 0...E-1 080 private double PENALIZATIONOVERSUBSCRIBED = 1E1; 081 private double PRECISIONFACTOR; 082 083 Random rng; 084 085 private Map<String, int[]> nodeGlobalIndexes_c; 086 private Map<String, int[]> frontierGlobalLinkIndexes_c; 087 private Map<String, int[]> internalAndFrontierGlobalLinkIndexes_c; 088 089 private Set<String> clusterIds; 090 private Map<Node, String> clusterIdOfNode; 091 // private Map<Long,Integer> globalLinkIndex; 092 // private long [] globalLinkId; 093 // private Map<Long,Integer> globalNodeIndex; 094 // private long [] globalNodeId; 095 private Map<String, Integer> globalClusterIndex; 096 private String[] globalClusterId; 097 098 private int[] frontierLinkGlobalIndex_ef; // one element per frontier link (in m_te and A_ce), returns its global index 099 100 private DoubleMatrix2D global_x_te; // indexes: global node, global link 101 private DoubleMatrix1D global_u_c; // indexes: global cluster 102 private DoubleMatrix2D global_m_te;// indexes: global node, global link (although only frontier link columns are used) 103 104 private Map<String, DoubleMatrix2D> mostUpdated_vte_fromC; 105 106 private DoubleMatrix2D global_h_nt; 107 private DoubleMatrix2D global_A_ne; 108 private DoubleMatrix1D global_u_e; 109 private DoubleMatrix2D global_h_ct; 110 private DoubleMatrix2D global_A_cef; 111 112 private Map<String, String> net2planParameters; 113 114 static final int ROUTING_MASTERUPDATE = 300; 115 static final int ROUTING_DOMAINRECOMPUTE = 301; 116 117 private TimeTrace stat_traceOf_objFun; 118 private TimeTrace stat_traceOf_x_te; 119 private TimeTrace stat_traceOf_y_e; 120 121 private NetPlan currentNetPlan; 122 123 @Override 124 public String getDescription() 125 { 126 return "This module implements a distributed primal-decomposition-based gradient algorithm, for a coordinated adjustment of the routing in multiple domains (or cluster, or autonomous systems) in a network, so that domains do not need to exchange sensitive internal information, and minimize the average number of hops in the network."; 127 } 128 129 @Override 130 public List<Triple<String, String, String>> getParameters() 131 { 132 /* Returns the parameter information for all the InputParameter objects defined in this object (uses Java reflection) */ 133 return InputParameter.getInformationAllInputParameterFieldsOfObject(this); 134 } 135 136 @Override 137 public void initialize(NetPlan currentNetPlan, Map<String, String> algorithmParameters, Map<String, String> simulationParameters, Map<String, String> net2planParameters) 138 { 139 try 140 { 141 /* Initialize all InputParameter objects defined in this object (this uses Java reflection) */ 142 InputParameter.initializeAllInputParameterFieldsOfObject(this, algorithmParameters); 143 144 if (currentNetPlan.getNumberOfLayers() != 1) throw new Net2PlanException("This algorithm works in single layer networks"); 145 146 this.rng = new Random(simulation_randomSeed.getLong()); 147 this.currentNetPlan = currentNetPlan; 148 this.net2planParameters = net2planParameters; 149 150 this.currentNetPlan.removeAllUnicastRoutingInformation(); 151 this.currentNetPlan.setRoutingType(RoutingType.HOP_BY_HOP_ROUTING); 152 153 this.N = this.currentNetPlan.getNumberOfNodes(); // total number of nodes in the network 154 this.E = this.currentNetPlan.getNumberOfLinks(); // total number of nodes in the network 155 this.D = this.currentNetPlan.getNumberOfDemands(); // total number of nodes in the network 156 if ((E == 0) || (D == 0)) throw new Net2PlanException ("The input design should have links and demands"); 157 this.allT = new int[N]; 158 for (int n = 0; n < N; n++) 159 allT[n] = n; 160 this.allE = new int[E]; 161 for (int e = 0; e < E; e++) 162 allE[e] = e; 163 this.PRECISIONFACTOR = 1e-2;//Double.parseDouble (net2planParameters.get("precisionFactor")); 164 165 this.clusterIds = new HashSet<String>(); 166 this.clusterIdOfNode = new HashMap<Node, String>(); 167 this.nodeGlobalIndexes_c = new HashMap<String, int[]>(); 168 this.frontierGlobalLinkIndexes_c = new HashMap<String, int[]>(); 169 this.internalAndFrontierGlobalLinkIndexes_c = new HashMap<String, int[]>(); 170 for (Node n : this.currentNetPlan.getNodes()) 171 { 172 final int indexN = n.getIndex(); 173 final String stClusterId = n.getAttribute("clusterId"); 174 if (stClusterId == null) throw new Net2PlanException("Node " + n + " does not belong to any cluster"); 175 clusterIds.add(stClusterId); 176 clusterIdOfNode.put(n, stClusterId); 177 if (!nodeGlobalIndexes_c.containsKey(stClusterId)) 178 { 179 nodeGlobalIndexes_c.put(stClusterId, new int[] {}); 180 frontierGlobalLinkIndexes_c.put(stClusterId, new int[] {}); 181 internalAndFrontierGlobalLinkIndexes_c.put(stClusterId, new int[] {}); 182 } 183 nodeGlobalIndexes_c.put(stClusterId, IntUtils.concatenate(nodeGlobalIndexes_c.get(stClusterId), new int[] { indexN })); 184 } 185 186 this.C = nodeGlobalIndexes_c.size(); 187 this.globalClusterIndex = new HashMap<String, Integer>(); 188 this.globalClusterId = new String[C]; 189 int counter = 0; 190 for (String cId : clusterIds) 191 { 192 globalClusterId[counter] = cId; 193 globalClusterIndex.put(cId, counter); 194 counter++; 195 } 196 197 /* For each cluster, the first links in xte are the frontier links, in the same order as they appear in cluster mte */ 198 this.frontierLinkGlobalIndex_ef = new int[0]; 199 for (Link e : this.currentNetPlan.getLinks()) 200 { 201 final int indexE = e.getIndex(); 202 final Node a_e = e.getOriginNode(); 203 final Node b_e = e.getDestinationNode(); 204 final String clusterId_ae = clusterIdOfNode.get(a_e); 205 final String clusterId_be = clusterIdOfNode.get(b_e); 206 if (!clusterId_ae.equals(clusterId_be)) 207 { 208 frontierGlobalLinkIndexes_c.put(clusterId_ae, IntUtils.concatenate(frontierGlobalLinkIndexes_c.get(clusterId_ae), new int[] { indexE })); 209 frontierGlobalLinkIndexes_c.put(clusterId_be, IntUtils.concatenate(frontierGlobalLinkIndexes_c.get(clusterId_be), new int[] { indexE })); 210 internalAndFrontierGlobalLinkIndexes_c.put(clusterId_ae, IntUtils.concatenate(internalAndFrontierGlobalLinkIndexes_c.get(clusterId_ae), new int[] { indexE })); 211 internalAndFrontierGlobalLinkIndexes_c.put(clusterId_be, IntUtils.concatenate(internalAndFrontierGlobalLinkIndexes_c.get(clusterId_be), new int[] { indexE })); 212 frontierLinkGlobalIndex_ef = IntUtils.concatenate(frontierLinkGlobalIndex_ef, new int[] { indexE }); 213 } 214 } 215 for (Link e : this.currentNetPlan.getLinks()) 216 { 217 final int indexE = e.getIndex(); 218 final Node a_e = e.getOriginNode(); 219 final Node b_e = e.getDestinationNode(); 220 final String clusterId_ae = clusterIdOfNode.get(a_e); 221 final String clusterId_be = clusterIdOfNode.get(b_e); 222 if (clusterId_ae.equals(clusterId_be)) internalAndFrontierGlobalLinkIndexes_c.put(clusterId_ae, IntUtils.concatenate(internalAndFrontierGlobalLinkIndexes_c.get(clusterId_ae), new int[] { indexE })); 223 } 224 225 this.Ef = frontierLinkGlobalIndex_ef.length; 226 227 this.global_x_te = DoubleFactory2D.sparse.make(N, E); 228 this.global_u_c = DoubleFactory1D.dense.make(C); 229 this.global_m_te = DoubleFactory2D.dense.make(N, E); 230 this.mostUpdated_vte_fromC = new HashMap<String, DoubleMatrix2D>(); 231 for (String cId : clusterIds) 232 mostUpdated_vte_fromC.put(cId, null); 233 this.global_A_ne = this.currentNetPlan.getMatrixNodeLinkIncidence(); 234 this.global_u_e = this.currentNetPlan.getVectorLinkCapacity(); 235 this.global_h_nt = this.currentNetPlan.getMatrixNode2NodeOfferedTraffic(); 236 for (int t = 0; t < N; t++) 237 { 238 double h_t = 0; 239 for (int n = 0; n < N; n++) 240 h_t += global_h_nt.get(n, t); 241 global_h_nt.set(t, t, -h_t); 242 } 243 244 this.global_h_ct = DoubleFactory2D.dense.make(C, N); 245 for (String clusterId : clusterIds) 246 { 247 final int index_c = globalClusterIndex.get(clusterId); 248 for (Node t : this.currentNetPlan.getNodes()) 249 { 250 final int index_t = t.getIndex(); 251 final int index_cluster_t = globalClusterIndex.get(clusterIdOfNode.get(t)); 252 final boolean targetInsideTheCluster = (index_c == index_cluster_t); 253 if (targetInsideTheCluster) 254 { 255 /* h_ct is what other nodes outside cluster generate towards t (in negative) */ 256 for (Node n : this.currentNetPlan.getNodes()) 257 { 258 final int index_n = n.getIndex(); 259 final int index_cluster_n = globalClusterIndex.get(clusterIdOfNode.get(n)); 260 if (index_cluster_n == index_c) continue; // only outside cluster nodes matter 261 global_h_ct.set(index_c, index_t, global_h_ct.get(index_c, index_t) - global_h_nt.get(index_n, index_t)); 262 } 263 } else 264 { 265 /* h_ct is what this cluster generate towards t (in negative) */ 266 for (int index_n : nodeGlobalIndexes_c.get(clusterId)) 267 { 268 final int index_cluster_n = globalClusterIndex.get(clusterIdOfNode.get(currentNetPlan.getNode(index_n))); 269 if (index_cluster_n != index_c) throw new RuntimeException("Bad"); 270 global_h_ct.set(index_c, index_t, global_h_ct.get(index_c, index_t) + global_h_nt.get(index_n, index_t)); 271 } 272 } 273 } 274 } 275 276 this.global_A_cef = DoubleFactory2D.dense.make(C, Ef); 277 for (int index_ef = 0; index_ef < frontierLinkGlobalIndex_ef.length; index_ef++) 278 { 279 final int index_e = frontierLinkGlobalIndex_ef[index_ef]; 280 final Link e = currentNetPlan.getLink(index_e); 281 final Node a_e = e.getOriginNode(); 282 final Node b_e = e.getDestinationNode(); 283 final int indexC_ae = this.globalClusterIndex.get(clusterIdOfNode.get(a_e)); 284 final int indexC_be = this.globalClusterIndex.get(clusterIdOfNode.get(b_e)); 285 if (indexC_ae != indexC_be) 286 { 287 global_A_cef.set(indexC_ae, index_ef, 1.0); 288 global_A_cef.set(indexC_be, index_ef, -1.0); 289 } 290 } 291 292 stat_traceOf_objFun = new TimeTrace(); 293 stat_traceOf_x_te = new TimeTrace(); 294 stat_traceOf_y_e = new TimeTrace(); 295 296 /* Initialize the master variables. The feasible, closest to zero */ 297 this.global_m_te = projectMaster(global_m_te, null); 298 299 /* Create the cluster topologies */ 300 /* Initially all nodes receive a "wake up to transmit" event, aligned at time zero or y asynchr => randomly chosen */ 301 for (String cId : clusterIds) 302 scheduleEvent(new SimEvent(0, SimEvent.DestinationModule.EVENT_PROCESSOR, ROUTING_DOMAINRECOMPUTE, cId)); 303 scheduleEvent(new SimEvent(0.0001, SimEvent.DestinationModule.EVENT_PROCESSOR, ROUTING_MASTERUPDATE, -1)); 304 305 checksInitialize(); 306 } catch (Exception e) 307 { 308 e.printStackTrace(); 309 throw e; 310 } 311 } 312 313 @Override 314 public void processEvent(NetPlan currentNetPlan, SimEvent event) 315 { 316 final double time = event.getEventTime(); 317 switch (event.getEventType()) 318 { 319 case ROUTING_DOMAINRECOMPUTE: 320 { 321 final String clusterId = (String) event.getEventObject(); 322 final int index_c = globalClusterIndex.get(clusterId); 323 324 /* Takes the most updated master information it has */ 325 Triple<DoubleMatrix2D, DoubleMatrix2D, Double> p = computeSubproblem(clusterId, this.currentNetPlan); 326 DoubleMatrix2D cluster_xte = p.getFirst(); 327 DoubleMatrix2D cluster_vte = p.getSecond(); 328 double u_c = p.getThird(); 329 330 final int[] cluster_Eif = internalAndFrontierGlobalLinkIndexes_c.get(clusterId); 331 332 /* global x_te can be not feasible, if different domains operate asynchornously */ 333 global_x_te.viewSelection(allT, cluster_Eif).assign(cluster_xte.copy().viewSelection(allT, cluster_Eif)); 334 global_u_c.set(index_c, u_c); 335 mostUpdated_vte_fromC.put(clusterId, cluster_vte.copy()); 336 337 /* They will be the new forwarding rules */ 338 DoubleMatrix2D new_fde = this.currentNetPlan.getMatrixDemandBasedForwardingRules(); 339 340 /* Remove all the forwarding rules of the cluster nodes */ 341 for (int index_n : nodeGlobalIndexes_c.get(clusterId)) 342 for (Link e : currentNetPlan.getNode(index_n).getOutgoingLinks()) 343 new_fde.viewColumn(e.getIndex()).assign(0); 344 345 /* Create new forwarding rules according to the routing computed */ 346 for (int index_n : nodeGlobalIndexes_c.get(clusterId)) 347 for (Node t : this.currentNetPlan.getNodes()) 348 { 349 final Node n = this.currentNetPlan.getNode(index_n); 350 if (n == t) continue; 351 final int index_t = t.getIndex(); 352 final Set<Demand> demandIdsEndingIn_t = t.getIncomingDemands(); 353 double outTraffic_nt = 0; 354 for (Link e : n.getOutgoingLinks()) 355 outTraffic_nt += global_x_te.get(index_t, e.getIndex()); 356 357 Link outLinkMaxFte = null; 358 double maxFte = -1; 359 double accumFteToForward = 0; 360 for (Link e : n.getOutgoingLinks()) 361 { 362 final double f_te = (outTraffic_nt == 0) ? 0 : global_x_te.get(index_t, e.getIndex()) / outTraffic_nt; 363 if (f_te > maxFte) 364 { 365 maxFte = f_te; 366 outLinkMaxFte = e; 367 } 368 if (f_te > PRECISIONFACTOR) 369 for (Demand d : demandIdsEndingIn_t) 370 new_fde.set(d.getIndex(), e.getIndex(), f_te); //changingFR.put (Pair.of (d,e),f_te); 371 else 372 accumFteToForward += f_te; 373 } 374 for (Demand d : demandIdsEndingIn_t) 375 new_fde.set(d.getIndex(), outLinkMaxFte.getIndex(), new_fde.get(d.getIndex(), outLinkMaxFte.getIndex()) + accumFteToForward); 376 } 377 378 DoubleMatrix2D outFR_dn = new_fde.zMult(currentNetPlan.getMatrixNodeLinkOutgoingIncidence().viewDice(), null); 379 for (int n = 0; n < N; n++) 380 for (int d = 0; d < D; d++) 381 if ((Math.abs(outFR_dn.get(d, n) - 1) > 1e-3) && (Math.abs(outFR_dn.get(d, n)) > 1e-3)) throw new RuntimeException("Bad. outFR_dn: " + outFR_dn); 382 383 /* In some rare occasions, the routing may create loops. These routings are not applied to the network */ 384 try { this.currentNetPlan.setForwardingRules(new_fde); } catch (ClosedCycleRoutingException e) { System.out.println ("The routing after this domain update would produce a closed loop. It is not applied."); } 385 386 if (time > 0) 387 { 388 DoubleMatrix1D blockedTraffic = currentNetPlan.getVectorDemandBlockedTraffic(); 389 DoubleMatrix1D linkOversubscription = currentNetPlan.getVectorLinkOversubscribedTraffic(); 390 if (blockedTraffic.getMaxLocation()[0] > 1e-3) System.out.println ("Some traffic is lost. Sum blocked traffic of the demands: " + blockedTraffic.zSum()); 391 if (linkOversubscription.getMaxLocation()[0] > 1e-3) System.out.println ("Some links are ovsersubscribed. Sum oversubscription traffic in the links: " + linkOversubscription.zSum()); 392 stat_traceOf_objFun.add(time, computeObjectiveFunctionFromNetPlan(this.currentNetPlan)); 393 stat_traceOf_x_te.add(time, this.currentNetPlan.getMatrixDestination2LinkTrafficCarried()); 394 stat_traceOf_y_e.add(time, this.currentNetPlan.getVectorLinkTotalCarriedTraffic()); 395 } 396 397 if (time > this.simulation_maxNumberOfMasterUpdateIntervals.getDouble() * routing_averageMasterUpdateInterval.getDouble()) 398 { 399 this.endSimulation(); 400 } 401 402 final double nextWakeUpTime = (routing_isDomainRoutingRecomputationSynchronous.getBoolean()) ? time + routing_averageDomainRoutingRecomputationInterval.getDouble() : Math.max(time, time + routing_averageDomainRoutingRecomputationInterval.getDouble() + routing_maxFluctuationDomainRoutingRecomputationInterval.getDouble() * (rng.nextDouble() - 0.5)); 403 this.scheduleEvent(new SimEvent(nextWakeUpTime, SimEvent.DestinationModule.EVENT_PROCESSOR, ROUTING_DOMAINRECOMPUTE, clusterId)); 404 405 break; 406 } 407 case ROUTING_MASTERUPDATE: // A node receives from an out neighbor the q_nt for any destination 408 { 409 DoubleMatrix2D gradient_te = DoubleFactory2D.dense.make(N, Ef, 1.0); // initially all ones 410 411 for (String clusterId : clusterIds) 412 { 413 DoubleMatrix2D v_te_fromC = mostUpdated_vte_fromC.get(clusterId); 414 gradient_te.assign(v_te_fromC.viewSelection(allT, frontierLinkGlobalIndex_ef), DoublePlusMultFirst.plusMult(1)); // -1 no 415 } 416 417 //final double gamma = gradient_gammaStep.getDouble() / (1+time*time); 418 final double gamma = gradient_gammaStep.getDouble(); 419 // final double gamma = gradient_gammaStep.getDouble() / (1+time); 420 // final double gamma = time < 20? gradient_gammaStep.getDouble() : gradient_gammaStep.getDouble() / (Math.sqrt(1+time)); 421 DoubleMatrix2D candidate_mte = global_m_te.copy(); 422 candidate_mte.viewSelection(allT, frontierLinkGlobalIndex_ef).assign(gradient_te, DoublePlusMultSecond.plusMult(-gamma)); 423 424 /* Now we have to project the candidate master variables into the set of feasible master variables */ 425 this.global_m_te = projectMaster(candidate_mte, global_m_te); 426 427 if (time > this.simulation_maxNumberOfMasterUpdateIntervals.getDouble() * routing_averageDomainRoutingRecomputationInterval.getDouble()) 428 { 429 endSimulation(); 430 } 431 432 final double nextMasterTime = Math.max(time, time + routing_averageMasterUpdateInterval.getDouble() + routing_maxFluctuationMasterUpdateInterval.getDouble() * (rng.nextDouble() - 0.5)); 433 scheduleEvent(new SimEvent(nextMasterTime, SimEvent.DestinationModule.EVENT_PROCESSOR, ROUTING_MASTERUPDATE, -1)); 434 break; 435 } 436 437 default: 438 throw new RuntimeException("Unexpected received event"); 439 } 440 } 441 442 public String finish(StringBuilder st, double simTime) 443 { 444 if (simulation_outFileNameRoot.getString().equals("")) return null; 445 stat_traceOf_objFun.printToFile(new File(simulation_outFileNameRoot.getString() + "_objFunc.txt")); 446 stat_traceOf_x_te.printToFile(new File(simulation_outFileNameRoot.getString() + "_xte.txt")); 447 stat_traceOf_y_e.printToFile(new File(simulation_outFileNameRoot.getString() + "_ye.txt")); 448 final NetPlan optNetPlan = computeOptimumSolution(); 449 TimeTrace.printToFile(new File(simulation_outFileNameRoot.getString() + "_jom_objFunc.txt"), new double[] { computeObjectiveFunctionFromNetPlan(optNetPlan) }); 450 TimeTrace.printToFile(new File(simulation_outFileNameRoot.getString() + "_jom_xte.txt"), optNetPlan.getMatrixDestination2LinkTrafficCarried()); 451 TimeTrace.printToFile(new File(simulation_outFileNameRoot.getString() + "_jom_ye.txt"), optNetPlan.getVectorLinkTotalCarriedTraffic()); 452 return null; 453 } 454 455 /* They use the most updated m_te information */ 456 private Triple<DoubleMatrix2D, DoubleMatrix2D, Double> computeSubproblem(String clusterId, NetPlan np) 457 { 458 final int[] cluster_Ef = frontierGlobalLinkIndexes_c.get(clusterId); 459 final int[] cluster_N = nodeGlobalIndexes_c.get(clusterId); 460 461 DoubleMatrix2D h_nt = global_h_nt.copy().viewSelection(cluster_N, allT); 462 DoubleMatrix2D A_ne = global_A_ne.copy().viewSelection(cluster_N, allE); 463 DoubleMatrix1D u_e = global_u_e.copy(); 464 465 OptimizationProblem op = new OptimizationProblem(); 466 467 op.setInputParameter("h_nt", h_nt.toArray()); 468 op.setInputParameter("A_ne", A_ne.toArray()); 469 op.setInputParameter("u_e", u_e.toArray(), "row"); 470 op.setInputParameter("m_te", global_m_te.toArray()); 471 op.setInputParameter("allT", allT, "row"); 472 op.setInputParameter("frontC", cluster_Ef, "row"); 473 op.setInputParameter("M", PENALIZATIONOVERSUBSCRIBED); 474 475 op.addDecisionVariable("x_te", false, new int[] { N, E }, 0, Double.MAX_VALUE); 476 op.addDecisionVariable("u", false, new int[] { 1, 1 }, 0, Double.MAX_VALUE); 477 478 /* Set some input parameters */ 479 op.setObjectiveFunction("minimize", "sum(x_te) + M*u"); 480 481 op.addConstraint("A_ne * (x_te') == h_nt"); // the flow-conservation constraints (NxN constraints) 482 op.addConstraint("sum(x_te,1) <= u_e + u"); // the capacity constraints (E constraints) 483 op.addConstraint("x_te(allT,frontC) == m_te(allT,frontC)", "v_te"); // the capacity constraints (E constraints) 484 485 /* Call the solver to solve the problem */ 486 op.solve("cplex"); 487 488 /* If an optimal solution was not found, quit */ 489 if (!op.solutionIsOptimal()) throw new RuntimeException("An optimal solution was not found"); 490 491 /* Retrieve the optimum solutions. Convert the bps into Erlangs */ 492 DoubleMatrix2D x_te = op.getPrimalSolution("x_te").view2D(); 493 double u = op.getPrimalSolution("u").toValue(); 494 495 /* v_te: zeros in no frontier links, the appropriate values in this cluster frontier links */ 496 DoubleMatrix2D v_te = DoubleFactory2D.dense.make(N, E); 497 v_te.viewSelection(allT, cluster_Ef).assign(op.getMultipliersOfConstraint("v_te").view2D()); 498 499 return Triple.of(x_te, v_te, u); 500 } 501 502 private DoubleMatrix2D projectMaster(DoubleMatrix2D can_mte, DoubleMatrix2D old_mte) 503 { 504 OptimizationProblem op = new OptimizationProblem(); 505 506 op.addDecisionVariable("m_te", false, new int[] { N, Ef }, 0, Double.MAX_VALUE); 507 op.setInputParameter("h_ct", global_h_ct.toArray()); 508 op.setInputParameter("A_ce", global_A_cef.toArray()); 509 op.setInputParameter("can_mte", can_mte.viewSelection(allT, frontierLinkGlobalIndex_ef).toArray()); 510 op.setInputParameter("u_e", global_u_e.viewSelection(frontierLinkGlobalIndex_ef).toArray(), "row"); 511 512 /* Set some input parameters */ 513 op.setObjectiveFunction("minimize", "sum((m_te-can_mte)^2)"); 514 op.addConstraint("A_ce * (m_te') == h_ct"); // the flow-conservation constraints (NxN constraints) 515 op.addConstraint("sum(m_te,1) <= u_e"); // the capacity constraints (E constraints) 516 517 for (String cId : clusterIds) 518 { 519 final int indexC = this.globalClusterIndex.get(cId); 520 final int[] nIdsCluster = nodeGlobalIndexes_c.get(cId); 521 final int[] frontierOutLinkIndexesInMte = DoubleUtils.find(global_A_cef.viewRow(indexC).toArray(), 1.0, SearchType.ALL); 522 op.setInputParameter("destCluster", nIdsCluster, "row"); 523 op.setInputParameter("clusterOutLinks", frontierOutLinkIndexesInMte, "row"); 524 op.addConstraint("m_te(destCluster,clusterOutLinks) == 0"); // traffic to a cluster does not leave the cluster 525 } 526 if (old_mte != null) op.setInitialSolution("m_te", new DoubleMatrixND(old_mte.viewSelection(allT, frontierLinkGlobalIndex_ef).copy())); 527 528 /* Call the solver to solve the problem */ 529 try 530 { 531 op.solve("ipopt"); 532 } catch (Exception e) 533 { 534 e.printStackTrace(); 535 } // sometimes IPOPT is unstable. do not stop the algorithm then 536 537 /* If an optimal solution was not found, quit */ 538 if (!op.solutionIsOptimal()) return old_mte.copy(); // sometimes IPOPT is unstable. do not stop the algorithm then 539 540 /* Retrieve the optimum solutions. Convert the bps into Erlangs */ 541 double[][] m_te_array = (double[][]) op.getPrimalSolution("m_te").toArray(); 542 543 for (String cId : clusterIds) 544 { 545 final int indexC = this.globalClusterIndex.get(cId); 546 final int[] nIdsCluster = nodeGlobalIndexes_c.get(cId); 547 final int[] frontierOutLinkIndexesInMte = DoubleUtils.find(global_A_cef.viewRow(indexC).toArray(), 1.0, SearchType.ALL); 548 for (int t : nIdsCluster) 549 for (int e : frontierOutLinkIndexesInMte) 550 if (m_te_array[t][e] > PRECISIONFACTOR) 551 throw new RuntimeException("Bad: m_te_array[t][e]: " + m_te_array[t][e]); 552 else 553 m_te_array[t][e] = 0; 554 } 555 556 DoubleMatrix2D m_te = DoubleFactory2D.dense.make(N, E); 557 m_te.viewSelection(allT, frontierLinkGlobalIndex_ef).assign(m_te_array); 558 return m_te; 559 } 560 561 private double computeObjectiveFunctionFromNetPlan(NetPlan np) 562 { 563 double optCost = 0; 564 double maxOverssubs = 0; 565 for (Link e : np.getLinks()) 566 { 567 final double y_e = e.getCarriedTrafficIncludingProtectionSegments(); 568 final double u_e = e.getCapacity(); 569 maxOverssubs = Math.max(maxOverssubs, y_e - u_e); 570 optCost += Math.min(y_e, u_e); 571 } 572 optCost += PENALIZATIONOVERSUBSCRIBED * maxOverssubs; 573 return optCost; 574 } 575 576 private NetPlan computeOptimumSolution() 577 { 578 /* Compute the optimum solution */ 579 Map<String, String> m = new HashMap<String, String>(); 580 m.put("optimizationTarget", "min-av-num-hops"); 581 m.put("solverName", "cplex"); 582 m.put("solverLibraryName", ""); 583 m.put("maxSolverTimeInSeconds", "-1"); 584 m.put("binaryRatePerTrafficUnit_bps", "1"); 585 m.put("averagePacketLengthInBytes", "1"); 586 m.put("nonBifurcatedRouting", "false"); 587 NetPlan npCopy = this.currentNetPlan.copy(); 588 new Offline_fa_xteFormulations().executeAlgorithm(npCopy, m, this.net2planParameters); 589 return npCopy; 590 } 591 592 private void checksInitialize() 593 { 594 /* Check in h_ct, the columns sum zero */ 595 for (int t = 0; t < this.currentNetPlan.getNumberOfNodes(); t++) 596 if (Math.abs(global_h_ct.viewColumn(t).zSum()) > PRECISIONFACTOR) throw new RuntimeException("Bad"); 597 /* Check in h_nt, the columns sum zero */ 598 for (int t = 0; t < this.currentNetPlan.getNumberOfNodes(); t++) 599 if (Math.abs(global_h_nt.viewColumn(t).zSum()) > PRECISIONFACTOR) throw new RuntimeException("Bad"); 600 } 601 602}