1# graphops.tcl -- 2# 3# Operations on and algorithms for graph data structures. 4# 5# Copyright (c) 2008 Alejandro Paz <vidriloco@gmail.com>, algorithm implementation 6# Copyright (c) 2008 Andreas Kupries, integration with Tcllib's struct::graph 7# 8# See the file "license.terms" for information on usage and redistribution 9# of this file, and for a DISCLAIMER OF ALL WARRANTIES. 10# 11# RCS: @(#) $Id: graphops.tcl,v 1.19 2009/09/24 19:30:10 andreas_kupries Exp $ 12 13# ### ### ### ######### ######### ######### 14## Requisites 15 16package require Tcl 8.5 17 18package require struct::disjointset ; # Used by kruskal 19package require struct::prioqueue ; # Used by kruskal, prim 20package require struct::queue ; # Used by isBipartite?, connectedComponent(Of) 21package require struct::stack ; # Used by tarjan 22package require struct::graph ; # isBridge, isCutVertex 23package require struct::tree ; # Used by BFS 24 25# ### ### ### ######### ######### ######### 26## 27 28namespace eval ::struct::graph::op {} 29 30# ### ### ### ######### ######### ######### 31## 32 33# This command constructs an adjacency matrix representation of the 34# graph argument. 35 36# Reference: http://en.wikipedia.org/wiki/Adjacency_matrix 37# 38# Note: The reference defines the matrix in such a way that some of 39# the limitations of the code here are not present. I.e. the 40# definition at wikipedia deals properly with arc directionality 41# and parallelism. 42# 43# TODO: Rework the code so that the result is in line with the reference. 44# Add features to handle weights as well. 45 46proc ::struct::graph::op::toAdjacencyMatrix {g} { 47 set nodeList [lsort -dict [$g nodes]] 48 # Note the lsort. This is used to impose some order on the matrix, 49 # for comparability of results. Otherwise different versions of 50 # Tcl and struct::graph (critcl) may generate different, yet 51 # equivalent matrices, dependent on things like the order a hash 52 # search is done, or nodes have been added to the graph, or ... 53 54 # Fill an array for index tracking later. Note how we start from 55 # index 1. This allows us avoid multiple expr+1 later on when 56 # iterating over the nodes and converting the names to matrix 57 # indices. See (*). 58 59 set i 1 60 foreach n $nodeList { 61 set nodeDict($n) $i 62 incr i 63 } 64 65 set matrix {} 66 lappend matrix [linsert $nodeList 0 {}] 67 68 # Setting up a template row with all of it's elements set to zero. 69 70 set baseRow 0 71 foreach n $nodeList { 72 lappend baseRow 0 73 } 74 75 foreach node $nodeList { 76 77 # The first element in every row is the name of its 78 # corresponding node. Using lreplace to overwrite the initial 79 # data in the template we get a copy apart from the template, 80 # which we can then modify further. 81 82 set currentRow [lreplace $baseRow 0 0 $node] 83 84 # Iterate over the neighbours, also known as 'adjacent' 85 # rows. The exact set of neighbours depends on the mode. 86 87 foreach neighbour [$g nodes -adj $node] { 88 # Set value for neighbour on this node list 89 set at $nodeDict($neighbour) 90 91 # (*) Here we avoid +1 due to starting from index 1 in the 92 # initialization of nodeDict. 93 set currentRow [lreplace $currentRow $at $at 1] 94 } 95 lappend matrix $currentRow 96 } 97 98 # The resulting matrix is a list of lists, size (n+1)^2 where n = 99 # number of nodes. First row and column (index 0) are node 100 # names. The other entries are boolean flags. True when an arc is 101 # present, False otherwise. The matrix represents an 102 # un-directional form of the graph with parallel arcs collapsed. 103 104 return $matrix 105} 106 107#Adjacency List 108#------------------------------------------------------------------------------------- 109#Procedure creates for graph G, it's representation as Adjacency List. 110# 111#In comparison to Adjacency Matrix it doesn't force using array with quite big 112#size - V^2, where V is a number of vertices ( instead, memory we need is about O(E) ). 113#It's especially important when concerning rare graphs ( graphs with amount of vertices 114#far bigger than amount of edges ). In practise, it turns out that generally, 115#Adjacency List is more effective. Moreover, going through the set of edges take 116#less time ( O(E) instead of O(E^2) ) and adding new edges is rapid. 117#On the other hand, checking if particular edge exists in graph G takes longer 118#( checking if edge {v1,v2} belongs to E(G) in proportion to min{deg(v1,v2)} ). 119#Deleting an edge is also longer - in proportion to max{ deg(v1), deg(v2) }. 120# 121#Input: 122# graph G ( directed or undirected ). Default is undirected. 123# 124#Output: 125# Adjacency List for graph G, represented by dictionary containing lists of adjacent nodes 126#for each node in G (key). 127# 128#Options: 129# -weights - adds to returning dictionary arc weights for each connection between nodes, so 130#each node returned by list as adjacent has additional parameter - weight of arc between him and 131#current node. 132# -directed - sets graph G to be interpreted as directed graph. 133# 134#Reference: 135#http://en.wikipedia.org/wiki/Adjacency_list 136# 137 138proc ::struct::graph::op::toAdjacencyList {G args} { 139 140 set arcTraversal "undirected" 141 set weightsOn 0 142 143 #options for procedure 144 foreach option $args { 145 switch -exact -- $option { 146 -directed { 147 set arcTraversal "directed" 148 } 149 -weights { 150 #checking if all edges have their weights set 151 VerifyWeightsAreOk $G 152 set weightsOn 1 153 } 154 default { 155 return -code error "Bad option \"$option\". Expected -directed or -weights" 156 } 157 } 158 } 159 160 set V [lsort -dict [$G nodes]] 161 162 #mainloop 163 switch -exact -- $arcTraversal { 164 undirected { 165 #setting up the Adjacency List with nodes 166 foreach v [lsort -dict [$G nodes]] { 167 dict set AdjacencyList $v {} 168 } 169 #appending the edges adjacent to nodes 170 foreach e [$G arcs] { 171 172 set v [$G arc source $e] 173 set u [$G arc target $e] 174 175 if { !$weightsOn } { 176 dict lappend AdjacencyList $v $u 177 dict lappend AdjacencyList $u $v 178 } else { 179 dict lappend AdjacencyList $v [list $u [$G arc getweight $e]] 180 dict lappend AdjacencyList $u [list $v [$G arc getweight $e]] 181 } 182 } 183 #deleting duplicated edges 184 foreach x [dict keys $AdjacencyList] { 185 dict set AdjacencyList $x [lsort -unique [dict get $AdjacencyList $x]] 186 } 187 } 188 directed { 189 foreach v $V { 190 set E [$G arcs -out $v] 191 set adjNodes {} 192 foreach e $E { 193 if { !$weightsOn } { 194 lappend adjNodes [$G arc target $e] 195 } else { 196 lappend adjNodes [list [$G arc target $e] [$G arc getweight $e]] 197 } 198 } 199 dict set AdjacencyList $v $adjNodes 200 } 201 } 202 default { 203 return -code error "Error while executing procedure" 204 } 205 } 206 207 return $AdjacencyList 208} 209 210#Bellman's Ford Algorithm 211#------------------------------------------------------------------------------------- 212#Searching for shortest paths between chosen node and 213#all other nodes in graph G. Based on relaxation method. In comparison to Dijkstra 214#it doesn't assume that all weights on edges are positive. However, this generality 215#costs us time complexity - O(V*E), where V is number of vertices and E is number 216#of edges. 217# 218#Input: 219#Directed graph G, weighted on edges and not containing 220#any cycles with negative sum of weights ( the presence of such cycles means 221#there is no shortest path, since the total weight becomes lower each time the 222#cycle is traversed ). Possible negative weights on edges. 223# 224#Output: 225#dictionary d[u] - distances from start node to each other node in graph G. 226# 227#Reference: http://en.wikipedia.org/wiki/Bellman-Ford_algorithm 228# 229 230proc ::struct::graph::op::BellmanFord { G startnode } { 231 232 #checking if all edges have their weights set 233 VerifyWeightsAreOk $G 234 235 #checking if the startnode exists in given graph G 236 if {![$G node exists $startnode]} { 237 return -code error "node \"$startnode\" does not exist in graph \"$G\"" 238 } 239 240 #sets of nodes and edges for graph G 241 set V [$G nodes] 242 set E [$G arcs] 243 244 #initialization 245 foreach i $V { 246 dict set distances $i Inf 247 } 248 249 dict set distances $startnode 0 250 251 #main loop (relaxation) 252 for { set i 1 } { $i <= ([dict size $distances]-1) } { incr i } { 253 254 foreach j $E { 255 set u [$G arc source $j] ;# start node of edge j 256 set v [$G arc target $j] ;# end node of edge j 257 258 if { [ dict get $distances $v ] > [ dict get $distances $u ] + [ $G arc getweight $j ]} { 259 dict set distances $v [ expr {[dict get $distances $u] + [$G arc getweight $j]} ] 260 } 261 } 262 } 263 264 #checking if there exists cycle with negative sum of weights 265 foreach i $E { 266 set u [$G arc source $i] ;# start node of edge i 267 set v [$G arc target $i] ;# end node of edge i 268 269 if { [dict get $distances $v] > [ dict get $distances $u ] + [$G arc getweight $i] } { 270 return -code error "Error. Given graph \"$G\" contains cycle with negative sum of weights." 271 } 272 } 273 274 return $distances 275 276} 277 278 279#Johnson's Algorithm 280#------------------------------------------------------------------------------------- 281#Searching paths between all pairs of vertices in graph. For rare graphs 282#asymptotically quicker than Floyd-Warshall's algorithm. Johnson's algorithm 283#uses Bellman-Ford's and Dijkstra procedures. 284# 285#Input: 286#Directed graph G, weighted on edges and not containing 287#any cycles with negative sum of weights ( the presence of such cycles means 288#there is no shortest path, since the total weight becomes lower each time the 289#cycle is traversed ). Possible negative weights on edges. 290#Possible options: 291# -filter ( returns only existing distances, cuts all Inf values for 292# non-existing connections between pairs of nodes ) 293# 294#Output: 295# Dictionary containing distances between all pairs of vertices 296# 297#Reference: http://en.wikipedia.org/wiki/Johnson_algorithm 298# 299 300proc ::struct::graph::op::Johnsons { G args } { 301 302 #options for procedure 303 set displaymode 0 304 foreach option $args { 305 switch -exact -- $option { 306 -filter { 307 set displaymode 1 308 } 309 default { 310 return -code error "Bad option \"$option\". Expected -filter" 311 } 312 } 313 } 314 315 #checking if all edges have their weights set 316 VerifyWeightsAreOk $G 317 318 #Transformation of graph G - adding one more node connected with 319 #each existing node with an edge, which weight is 0 320 set V [$G nodes] 321 set s [$G node insert] 322 323 foreach i $V { 324 if { $i ne $s } { 325 $G arc insert $s $i 326 } 327 } 328 329 $G arc setunweighted 330 331 #set potential values with Bellman-Ford's 332 set h [BellmanFord $G $s] 333 334 #transformed graph no needed longer - deleting added node and edges 335 $G node delete $s 336 337 #setting new weights for edges in graph G 338 foreach i [$G arcs] { 339 set u [$G arc source $i] 340 set v [$G arc target $i] 341 342 lappend weights [$G arc getweight $i] 343 $G arc setweight $i [ expr { [$G arc getweight $i] + [dict get $h $u] - [dict get $h $v] } ] 344 } 345 346 #finding distances between all pair of nodes with Dijkstra started from each node 347 foreach i [$G nodes] { 348 set dijkstra [dijkstra $G $i -arcmode directed -outputformat distances] 349 350 foreach j [$G nodes] { 351 if { $i ne $j } { 352 if { $displaymode eq 1 } { 353 if { [dict get $dijkstra $j] ne "Inf" } { 354 dict set values [list $i $j] [ expr {[ dict get $dijkstra $j] - [dict get $h $i] + [dict get $h $j]} ] 355 } 356 } else { 357 dict set values [list $i $j] [ expr {[ dict get $dijkstra $j] - [dict get $h $i] + [dict get $h $j]} ] 358 } 359 } 360 } 361 } 362 363 #setting back edge weights for graph G 364 set k 0 365 foreach i [$G arcs] { 366 $G arc setweight $i [ lindex $weights $k ] 367 incr k 368 } 369 370 return $values 371} 372 373 374#Floyd-Warshall's Algorithm 375#------------------------------------------------------------------------------------- 376#Searching shortest paths between all pairs of edges in weighted graphs. 377#Time complexity: O(V^3) - where V is number of vertices. 378#Memory complexity: O(V^2) 379#Input: directed weighted graph G 380#Output: dictionary containing shortest distances to each node from each node 381# 382#Algorithm finds solutions dynamically. It compares all possible paths through the graph 383#between each pair of vertices. Graph shouldn't possess any cycle with negative 384#sum of weights ( the presence of such cycles means there is no shortest path, 385#since the total weight becomes lower each time the cycle is traversed ). 386#On the other hand algorithm can be used to find those cycles - if any shortest distance 387#found by algorithm for any nodes v and u (when v is the same node as u) is negative, 388#that node surely belong to at least one negative cycle. 389# 390#Reference: http://en.wikipedia.org/wiki/Floyd-Warshall_algorithm 391# 392 393proc ::struct::graph::op::FloydWarshall { G } { 394 395 VerifyWeightsAreOk $G 396 397 foreach v1 [$G nodes] { 398 foreach v2 [$G nodes] { 399 dict set values [list $v1 $v2] Inf 400 } 401 dict set values [list $v1 $v1] 0 402 } 403 404 foreach e [$G arcs] { 405 set v1 [$G arc source $e] 406 set v2 [$G arc target $e] 407 dict set values [list $v1 $v2] [$G arc getweight $e] 408 } 409 410 foreach u [$G nodes] { 411 foreach v1 [$G nodes] { 412 foreach v2 [$G nodes] { 413 414 set x [dict get $values [list $v1 $u]] 415 set y [dict get $values [list $u $v2]] 416 set d [ expr {$x + $y}] 417 418 if { [dict get $values [list $v1 $v2]] > $d } { 419 dict set values [list $v1 $v2] $d 420 } 421 } 422 } 423 } 424 #finding negative cycles 425 foreach v [$G nodes] { 426 if { [dict get $values [list $v $v]] < 0 } { 427 return -code error "Error. Given graph \"$G\" contains cycle with negative sum of weights." 428 } 429 } 430 431 return $values 432} 433 434#Metric Travelling Salesman Problem (TSP) - 2 approximation algorithm 435#------------------------------------------------------------------------------------- 436#Travelling salesman problem is a very popular problem in graph theory, where 437#we are trying to find minimal Hamilton cycle in weighted complete graph. In other words: 438#given a list of cities (nodes) and their pairwise distances (edges), the task is to find 439#a shortest possible tour that visits each city exactly once. 440#TSP problem is NP-Complete, so there is no efficient algorithm to solve it. Greedy methods 441#are getting extremely slow, with the increase in the set of nodes. 442# 443#For this algorithm we consider a case when for given graph G, the triangle inequality is 444#satisfied. So for example, for any three nodes A, B and C the distance between A and C must 445#be at most the distance from A to B plus the distance from B to C. What's important 446#most of the considered cases in TSP problem will satisfy this condition. 447# 448#Input: undirected, weighted graph G 449#Output: approximated solution of minimum Hamilton Cycle - closed path visiting all nodes, 450#each exactly one time. 451# 452#Reference: http://en.wikipedia.org/wiki/Travelling_salesman_problem 453# 454 455proc ::struct::graph::op::MetricTravellingSalesman { G } { 456 457 #checking if graph is connected 458 if { ![isConnected? $G] } { 459 return -code error "Error. Given graph \"$G\" is not a connected graph." 460 } 461 #checking if all weights are set 462 VerifyWeightsAreOk $G 463 464 # Extend graph to make it complete. 465 # NOTE: The graph is modified in place. 466 createCompleteGraph $G originalEdges 467 468 #create minimum spanning tree for graph G 469 set T [prim $G] 470 471 #TGraph - spanning tree of graph G 472 #filling TGraph with edges and nodes 473 set TGraph [createTGraph $G $T 0] 474 475 #finding Hamilton cycle 476 set result [findHamiltonCycle $TGraph $originalEdges $G] 477 478 $TGraph destroy 479 480 # Note: Fleury, which is the algorithm used to find our the cycle 481 # (inside of isEulerian?) is inherently directionless, i.e. it 482 # doesn't care about arc direction. This does not matter if our 483 # input is a symmetric graph, i.e. u->v and v->u have the same 484 # weight for all nodes u, v in G, u != v. But for an asymmetric 485 # graph as our input we really have to check the two possible 486 # directions of the returned tour for the one with the smaller 487 # weight. See test case MetricTravellingSalesman-1.1 for an 488 # exmaple. 489 490 set w {} 491 foreach a [$G arcs] { 492 set u [$G arc source $a] 493 set v [$G arc target $a] 494 set uv [list $u $v] 495 # uv = <$G arc nodes $arc> 496 dict set w $uv [$G arc getweight $a] 497 } 498 foreach k [dict keys $w] { 499 lassign $k u v 500 set vu [list $v $u] 501 if {[dict exists $w $vu]} continue 502 dict set w $vu [dict get $w $k] 503 } 504 505 set reversed [lreverse $result] 506 507 if {[TourWeight $w $result] > [TourWeight $w $reversed]} { 508 return $reversed 509 } 510 return $result 511} 512 513proc ::struct::graph::op::TourWeight {w tour} { 514 set total 0 515 foreach \ 516 u [lrange $tour 0 end-1] \ 517 v [lrange $tour 1 end] { 518 set uv [list $u $v] 519 set total [expr { 520 $total + 521 [dict get $w $uv] 522 }] 523 } 524 return $total 525} 526 527#Christofides Algorithm - for Metric Travelling Salesman Problem (TSP) 528#------------------------------------------------------------------------------------- 529#Travelling salesman problem is a very popular problem in graph theory, where 530#we are trying to find minimal Hamilton cycle in weighted complete graph. In other words: 531#given a list of cities (nodes) and their pairwise distances (edges), the task is to find 532#a shortest possible tour that visits each city exactly once. 533#TSP problem is NP-Complete, so there is no efficient algorithm to solve it. Greedy methods 534#are getting extremely slow, with the increase in the set of nodes. 535# 536#For this algorithm we consider a case when for given graph G, the triangle inequality is 537#satisfied. So for example, for any three nodes A, B and C the distance between A and C must 538#be at most the distance from A to B plus the distance from B to C. What's important 539#most of the considered cases in TSP problem will satisfy this condition. 540# 541#Christofides is a 3/2 approximation algorithm. For a graph given at input, it returns 542#found Hamilton cycle (list of nodes). 543# 544#Reference: http://en.wikipedia.org/wiki/Christofides_algorithm 545# 546 547proc ::struct::graph::op::Christofides { G } { 548 549 #checking if graph is connected 550 if { ![isConnected? $G] } { 551 return -code error "Error. Given graph \"$G\" is not a connected graph." 552 } 553 #checking if all weights are set 554 VerifyWeightsAreOk $G 555 556 createCompleteGraph $G originalEdges 557 558 #create minimum spanning tree for graph G 559 set T [prim $G] 560 561 #setting graph algorithm is working on - spanning tree of graph G 562 set TGraph [createTGraph $G $T 1] 563 564 set oddTGraph [struct::graph] 565 566 foreach v [$TGraph nodes] { 567 if { [$TGraph node degree $v] % 2 == 1 } { 568 $oddTGraph node insert $v 569 } 570 } 571 572 #create complete graph 573 foreach v [$oddTGraph nodes] { 574 foreach u [$oddTGraph nodes] { 575 if { ($u ne $v) && ![$oddTGraph arc exists [list $u $v]] } { 576 $oddTGraph arc insert $v $u [list $v $u] 577 $oddTGraph arc setweight [list $v $u] [distance $G $v $u] 578 } 579 580 } 581 } 582 583 #### 584 # MAX MATCHING HERE!!! 585 #### 586 set M [GreedyMaxMatching $oddTGraph] 587 588 foreach e [$oddTGraph arcs] { 589 if { ![struct::set contains $M $e] } { 590 $oddTGraph arc delete $e 591 } 592 } 593 594 #operation: M + T 595 foreach e [$oddTGraph arcs] { 596 set u [$oddTGraph arc source $e] 597 set v [$oddTGraph arc target $e] 598 set uv [list $u $v] 599 600 # Check if the arc in max-matching is parallel or not, to make 601 # sure that we always insert an anti-parallel arc. 602 603 if {[$TGraph arc exists $uv]} { 604 set vu [list $v $u] 605 $TGraph arc insert $v $u $vu 606 $TGraph arc setweight $vu [$oddTGraph arc getweight $e] 607 } else { 608 $TGraph arc insert $u $v $uv 609 $TGraph arc setweight $uv [$oddTGraph arc getweight $e] 610 } 611 } 612 613 #finding Hamilton Cycle 614 set result [findHamiltonCycle $TGraph $originalEdges $G] 615 $oddTGraph destroy 616 $TGraph destroy 617 return $result 618} 619 620#Greedy Max Matching procedure, which finds maximal ( not maximum ) matching 621#for given graph G. It adds edges to solution, beginning from edges with the 622#lowest cost. 623 624proc ::struct::graph::op::GreedyMaxMatching {G} { 625 626 set maxMatch {} 627 628 foreach e [sortEdges $G] { 629 set v [$G arc source $e] 630 set u [$G arc target $e] 631 set neighbours [$G arcs -adj $v $u] 632 set noAdjacentArcs 1 633 634 lremove neighbours $e 635 636 foreach a $neighbours { 637 if { $a in $maxMatch } { 638 set noAdjacentArcs 0 639 break 640 } 641 } 642 if { $noAdjacentArcs } { 643 lappend maxMatch $e 644 } 645 } 646 647 return $maxMatch 648} 649 650#Subprocedure which for given graph G, returns the set of edges 651#sorted with their costs. 652proc ::struct::graph::op::sortEdges {G} { 653 set weights [$G arc weights] 654 655 # NOTE: Look at possible rewrite, simplification. 656 657 set sortedEdges {} 658 659 foreach val [lsort [dict values $weights]] { 660 foreach x [dict keys $weights] { 661 if { [dict get $weights $x] == $val } { 662 set weights [dict remove $weights $x] 663 lappend sortedEdges $x ;#[list $val $x] 664 } 665 } 666 } 667 668 return $sortedEdges 669} 670 671#Subprocedure, which for given graph G, returns the dictionary 672#containing edges sorted by weights (sortMode -> weights) or 673#nodes sorted by degree (sortMode -> degrees). 674 675proc ::struct::graph::op::sortGraph {G sortMode} { 676 677 switch -exact -- $sortMode { 678 weights { 679 set weights [$G arc weights] 680 foreach val [lsort [dict values $weights]] { 681 foreach x [dict keys $weights] { 682 if { [dict get $weights $x] == $val } { 683 set weights [dict remove $weights $x] 684 dict set sortedVals $x $val 685 } 686 } 687 } 688 } 689 degrees { 690 foreach v [$G nodes] { 691 dict set degrees $v [$G node degree $v] 692 } 693 foreach x [lsort -integer -decreasing [dict values $degrees]] { 694 foreach y [dict keys $degrees] { 695 if { [dict get $degrees $y] == $x } { 696 set degrees [dict remove $degrees $y] 697 dict set sortedVals $y $x 698 } 699 } 700 } 701 } 702 default { 703 return -code error "Unknown sort mode \"$sortMode\", expected weights, or degrees" 704 } 705 } 706 707 return $sortedVals 708} 709 710#Finds Hamilton cycle in given graph G 711#Procedure used by Metric TSP Algorithms: 712#Christofides and Metric TSP 2-approximation algorithm 713 714proc ::struct::graph::op::findHamiltonCycle {G originalEdges originalGraph} { 715 716 isEulerian? $G tourvar tourstart 717 718 # Note: The start node is not necessarily the source node of the 719 # first arc in the tour. The Fleury in isEulerian? may have walked 720 # the arcs against! their direction. See also the note in our 721 # caller (MetricTravellingSalesman). 722 723 # Instead of reconstructing the start node by intersecting the 724 # node-set for first and last arc, we are taking the easy and get 725 # it directly from isEulerian?, as that command knows which node 726 # it had chosen for this. 727 728 lappend result $tourstart 729 lappend tourvar [lindex $tourvar 0] 730 731 set v $tourstart 732 foreach i $tourvar { 733 set u [$G node opposite $v $i] 734 735 if { $u ni $result } { 736 set va [lindex $result end] 737 set vb $u 738 739 if { ([list $va $vb] in $originalEdges) || ([list $vb $va] in $originalEdges) } { 740 lappend result $u 741 } else { 742 743 set path [dict get [dijkstra $G $va] $vb] 744 745 #reversing the path 746 set path [lreverse $path] 747 #cutting the start element 748 set path [lrange $path 1 end] 749 750 #adding the path and the target element 751 lappend result {*}$path 752 lappend result $vb 753 } 754 } 755 set v $u 756 } 757 758 set path [dict get [dijkstra $originalGraph [lindex $result 0]] [lindex $result end]] 759 set path [lreverse $path] 760 761 set path [lrange $path 1 end] 762 763 if { [llength $path] } { 764 lappend result {*}$path 765 } 766 767 lappend result $tourstart 768 return $result 769} 770 771#Subprocedure for TSP problems. 772# 773#Creating graph from sets of given nodes and edges. 774#In option doubledArcs we decide, if we want edges to be 775#duplicated or not: 776#0 - duplicated (Metric TSP 2-approximation algorithm) 777#1 - single (Christofides Algorithm) 778# 779#Note that it assumes that graph's edges are properly weighted. That 780#condition is checked before in procedures that use createTGraph, but for 781#other uses it should be taken into consideration. 782# 783 784proc ::struct::graph::op::createTGraph {G Edges doubledArcs} { 785 #checking if given set of edges is proper (all edges are in graph G) 786 foreach e $Edges { 787 if { ![$G arc exists $e] } { 788 return -code error "Edge \"$e\" doesn't exist in graph \"$G\". Set the proper set of edges." 789 } 790 } 791 792 set TGraph [struct::graph] 793 794 #fill TGraph with nodes 795 foreach v [$G nodes] { 796 $TGraph node insert 797 } 798 799 #fill TGraph with arcs 800 foreach e $Edges { 801 set v [$G arc source $e] 802 set u [$G arc target $e] 803 if { ![$TGraph arc exists [list $u $v]] } { 804 $TGraph arc insert $u $v [list $u $v] 805 $TGraph arc setweight [list $u $v] [$G arc getweight $e] 806 } 807 if { !$doubledArcs } { 808 if { ![$TGraph arc exists [list $v $u]] } { 809 $TGraph arc insert $v $u [list $v $u] 810 $TGraph arc setweight [list $v $u] [$G arc getweight $e] 811 } 812 } 813 } 814 815 return $TGraph 816} 817 818#Subprocedure for some algorithms, e.g. TSP algorithms. 819# 820#It returns graph filled with arcs missing to say that graph is complete. 821#Also it sets variable originalEdges with edges, which existed in given 822#graph G at beginning, before extending the set of edges. 823# 824 825proc ::struct::graph::op::createCompleteGraph {G originalEdges} { 826 827 upvar $originalEdges st 828 set st {} 829 foreach e [$G arcs] { 830 set v [$G arc source $e] 831 set u [$G arc target $e] 832 833 lappend st [list $v $u] 834 } 835 836 foreach v [$G nodes] { 837 foreach u [$G nodes] { 838 if { ($u != $v) && ([list $v $u] ni $st) && ([list $u $v] ni $st) && ![$G arc exists [list $u $v]] } { 839 $G arc insert $v $u [list $v $u] 840 $G arc setweight [list $v $u] Inf 841 } 842 } 843 } 844 return $G 845} 846 847 848#Maximum Cut - 2 approximation algorithm 849#------------------------------------------------------------------------------------- 850#Maximum cut problem is a problem finding a cut not smaller than any other cut. In 851#other words, we divide set of nodes for graph G into such 2 sets of nodes U and V, 852#that the amount of edges connecting U and V is as high as possible. 853# 854#Algorithm is a 2-approximation, so for ALG ( solution returned by Algorithm) and 855#OPT ( optimal solution), such inequality is true: OPT <= 2 * ALG. 856# 857#Input: 858#Graph G 859#U - variable storing first set of nodes (cut) given by solution 860#V - variable storing second set of nodes (cut) given by solution 861# 862#Output: 863#Algorithm returns number of edges between found two sets of nodes. 864# 865#Reference: http://en.wikipedia.org/wiki/Maxcut 866# 867 868proc ::struct::graph::op::MaxCut {G U V} { 869 870 upvar $U _U 871 upvar $V _V 872 873 set _U {} 874 set _V {} 875 set counter 0 876 877 foreach {u v} [lsort -dict [$G nodes]] { 878 lappend _U $u 879 if {$v eq ""} continue 880 lappend _V $v 881 } 882 883 set val 1 884 set ALG [countEdges $G $_U $_V] 885 while {$val>0} { 886 set val [cut $G _U _V $ALG] 887 if { $val > $ALG } { 888 set ALG $val 889 } 890 } 891 return $ALG 892} 893 894#procedure replaces nodes between sets and checks if that change is profitable 895proc ::struct::graph::op::cut {G Uvar Vvar param} { 896 897 upvar $Uvar U 898 upvar $Vvar V 899 set _V {} 900 set _U {} 901 set value 0 902 903 set maxValue $param 904 set _U $U 905 set _V $V 906 907 foreach v [$G nodes] { 908 909 if { $v ni $_U } { 910 lappend _U $v 911 lremove _V $v 912 set value [countEdges $G $_U $_V] 913 } else { 914 lappend _V $v 915 lremove _U $v 916 set value [countEdges $G $_U $_V] 917 } 918 919 if { $value > $maxValue } { 920 set U $_U 921 set V $_V 922 set maxValue $value 923 } else { 924 set _V $V 925 set _U $U 926 } 927 } 928 929 set value $maxValue 930 931 if { $value > $param } { 932 return $value 933 } else { 934 return 0 935 } 936} 937 938#Removing element from the list - auxiliary procedure 939proc ::struct::graph::op::lremove {listVariable value} { 940 upvar 1 $listVariable var 941 set idx [lsearch -exact $var $value] 942 set var [lreplace $var $idx $idx] 943} 944 945#procedure counts edges that link two sets of nodes 946proc ::struct::graph::op::countEdges {G U V} { 947 948 set value 0 949 950 foreach u $U { 951 foreach e [$G arcs -out $u] { 952 set v [$G arc target $e] 953 if {$v ni $V} continue 954 incr value 955 } 956 } 957 foreach v $V { 958 foreach e [$G arcs -out $v] { 959 set u [$G arc target $e] 960 if {$u ni $U} continue 961 incr value 962 } 963 } 964 965 return $value 966} 967 968#K-Center Problem - 2 approximation algorithm 969#------------------------------------------------------------------------------------- 970#Input: 971#Undirected complete graph G, which satisfies triangle inequality. 972#k - positive integer 973# 974#Definition: 975#For any set S ( which is subset of V ) and node v, let the connect(v,S) be the 976#cost of cheapest edge connecting v with any node in S. The goal is to find 977#such S, that |S| = k and max_v{connect(v,S)} is possibly small. 978# 979#In other words, we can use it i.e. for finding best locations in the city ( nodes 980#of input graph ) for placing k buildings, such that those buildings will be as close 981#as possible to all other locations in town. 982# 983#Output: 984#set of nodes - k center for graph G 985# 986 987proc ::struct::graph::op::UnweightedKCenter {G k} { 988 989 #checking if all weights for edges in graph G are set well 990 VerifyWeightsAreOk $G 991 992 #checking if proper value of k is given at input 993 if { $k <= 0 } { 994 return -code error "The \"k\" value must be an positive integer." 995 } 996 997 set j [ expr {$k+1} ] 998 999 #variable for holding the graph G(i) in each iteration 1000 set Gi [struct::graph] 1001 #two squared graph G 1002 set GiSQ [struct::graph] 1003 #sorted set of edges for graph G 1004 set arcs [sortEdges $G] 1005 1006 #initializing both graph variables 1007 foreach v [$G nodes] { 1008 $Gi node insert $v 1009 $GiSQ node insert $v 1010 } 1011 1012 #index i for each iteration 1013 1014 #we seek for final solution, as long as the max independent 1015 #set Mi (found in particular iterations), such that |Mi| <= k, is found. 1016 for {set index 0} {$j > $k} {incr index} { 1017 #source node of an edge we add in current iteration 1018 set u [$G arc source [lindex $arcs $index]] 1019 #target node of an edge we add in current iteration 1020 set v [$G arc target [lindex $arcs $index]] 1021 1022 #adding edge Ei to graph G(i) 1023 $Gi arc insert $u $v [list $u $v] 1024 #extending G(i-1)**2 to G(i)**2 using G(i) 1025 set GiSQ [extendTwoSquaredGraph $GiSQ $Gi $u $v] 1026 1027 #finding maximal independent set for G(i)**2 1028 set Mi [GreedyMaxIndependentSet $GiSQ] 1029 1030 #number of nodes in maximal independent set that was found 1031 set j [llength $Mi] 1032 } 1033 1034 $Gi destroy 1035 $GiSQ destroy 1036 return $Mi 1037} 1038 1039#Weighted K-Center - 3 approximation algorithm 1040#------------------------------------------------------------------------------------- 1041# 1042#The variation of unweighted k-center problem. Besides the fact graph is edge-weighted, 1043#there are also weights on vertices of input graph G. We've got also restriction 1044#W. The goal is to choose such set of nodes S ( which is a subset of V ), that it's 1045#total weight is not greater than W and also function: max_v { min_u { cost(u,v) }} 1046#has the smallest possible worth ( v is a node in V and u is a node in S ). 1047# 1048#Note: 1049#For more information about K-Center problem check Unweighted K-Center algorithm 1050#description. 1051 1052proc ::struct::graph::op::WeightedKCenter {G nodeWeights W} { 1053 1054 #checking if all weights for edges in graph G are set well 1055 VerifyWeightsAreOk $G 1056 1057 #checking if proper value of k is given at input 1058 if { $W <= 0 } { 1059 return -code error "The \"W\" value must be an positive integer." 1060 } 1061 #initilization 1062 set j [ expr {$W+1} ] 1063 1064 #graphs G(i) and G(i)**2 1065 set Gi [struct::graph] 1066 set GiSQ [struct::graph] 1067 #the set of arcs for graph G sorted with their weights (increasing) 1068 set arcs [sortEdges $G] 1069 1070 #initialization of graphs G(i) and G(i)**2 1071 foreach v [$G nodes] { 1072 $Gi node insert $v 1073 $GiSQ node insert $v 1074 } 1075 1076 #the main loop - iteration over all G(i)'s and G(i)**2's, 1077 #extended with each iteration till the solution is found 1078 1079 foreach arc $arcs { 1080 #initilization of the set of nodes, which are cheapest neighbours 1081 #for particular nodes in maximal independent set 1082 set Si {} 1083 1084 set u [$G arc source $arc] 1085 set v [$G arc target $arc] 1086 1087 #extending graph G(i) 1088 $Gi arc insert $u $v [list $u $v] 1089 1090 #extending graph G(i)**2 from G(i-1)**2 using G(i) 1091 set GiSQ [extendTwoSquaredGraph $GiSQ $Gi $u $v] 1092 1093 #finding maximal independent set (Mi) for graph G(i)**2 found in the 1094 #previous step. Mi is found using greedy algorithm that also considers 1095 #weights on vertices. 1096 set Mi [GreedyWeightedMaxIndependentSet $GiSQ $nodeWeights] 1097 1098 #for each node u in Maximal Independent set found in previous step, 1099 #we search for its cheapest ( considering costs at vertices ) neighbour. 1100 #Note that node u is considered as it is a neighbour for itself. 1101 foreach u $Mi { 1102 1103 set minWeightOfSi Inf 1104 1105 #the neighbours of u 1106 set neighbours [$Gi nodes -adj $u] 1107 set smallestNeighbour 0 1108 #u is a neighbour for itself 1109 lappend neighbours $u 1110 1111 #finding neighbour with minimal cost 1112 foreach w [lsort -index 1 $nodeWeights] { 1113 lassign $w node weight 1114 if {[struct::set contains $neighbours $node]} { 1115 set minWeightOfSi $weight 1116 set smallestNeighbour $node 1117 break 1118 } 1119 } 1120 1121 lappend Si [list $smallestNeighbour $minWeightOfSi] 1122 } 1123 1124 set totalSiWeight 0 1125 set possibleSolution {} 1126 1127 foreach s $Si { 1128 #counting the total weight of the set of nodes - Si 1129 set totalSiWeight [ expr { $totalSiWeight + [lindex $s 1] } ] 1130 1131 #it's final solution, if weight found in previous step is 1132 #not greater than W 1133 lappend possibleSolution [lindex $s 0] 1134 } 1135 1136 #checking if final solution is found 1137 if { $totalSiWeight <= $W } { 1138 $Gi destroy 1139 $GiSQ destroy 1140 return $possibleSolution 1141 } 1142 } 1143 1144 $Gi destroy 1145 $GiSQ destroy 1146 1147 #no solution found - error returned 1148 return -code error "No k-center found for restriction W = $W" 1149 1150} 1151 1152#Maximal Independent Set - 2 approximation greedy algorithm 1153#------------------------------------------------------------------------------------- 1154# 1155#A maximal independent set is an independent set such that adding any other node 1156#to the set forces the set to contain an edge. 1157# 1158#Note: 1159#Don't confuse it with maximum independent set, which is a largest independent set 1160#for a given graph G. 1161# 1162#Reference: http://en.wikipedia.org/wiki/Maximal_independent_set 1163 1164proc ::struct::graph::op::GreedyMaxIndependentSet {G} { 1165 1166 set result {} 1167 set nodes [$G nodes] 1168 1169 foreach v $nodes { 1170 if { [struct::set contains $nodes $v] } { 1171 lappend result $v 1172 1173 foreach neighbour [$G nodes -adj $v] { 1174 struct::set exclude nodes $neighbour 1175 } 1176 } 1177 } 1178 1179 return $result 1180} 1181 1182#Weighted Maximal Independent Set - 2 approximation greedy algorithm 1183#------------------------------------------------------------------------------------- 1184# 1185#Weighted variation of Maximal Independent Set. It takes as an input argument 1186#not only graph G but also set of weights for all vertices in graph G. 1187# 1188#Note: 1189#Read also Maximal Independent Set description for more info. 1190# 1191#Reference: http://en.wikipedia.org/wiki/Maximal_independent_set 1192 1193proc ::struct::graph::op::GreedyWeightedMaxIndependentSet {G nodeWeights} { 1194 1195 set result {} 1196 set nodes {} 1197 foreach v [lsort -index 1 $nodeWeights] { 1198 lappend nodes [lindex $v 0] 1199 } 1200 1201 foreach v $nodes { 1202 if { [struct::set contains $nodes $v] } { 1203 lappend result $v 1204 1205 set neighbours [$G nodes -adj $v] 1206 1207 foreach neighbour [$G nodes -adj $v] { 1208 struct::set exclude nodes $neighbour 1209 } 1210 } 1211 } 1212 1213 return $result 1214} 1215 1216#subprocedure creating from graph G two squared graph 1217#G^2 - graph in which edge between nodes u and v exists, 1218#if and only if, when distance (in edges, not weights) 1219#between those nodes is not greater than 2 and u != v. 1220 1221proc ::struct::graph::op::createSquaredGraph {G} { 1222 1223 set H [struct::graph] 1224 foreach v [$G nodes] { 1225 $H node insert $v 1226 } 1227 1228 foreach v [$G nodes] { 1229 foreach u [$G nodes -adj $v] { 1230 if { ($v != $u) && ![$H arc exists [list $v $u]] && ![$H arc exists [list $u $v]] } { 1231 $H arc insert $u $v [list $u $v] 1232 } 1233 foreach z [$G nodes -adj $u] { 1234 if { ($v != $z) && ![$H arc exists [list $v $z]] && ![$H arc exists [list $z $v]] } { 1235 $H arc insert $v $z [list $v $z] 1236 } 1237 } 1238 } 1239 } 1240 1241 return $H 1242} 1243 1244#subprocedure for Metric K-Center problem 1245# 1246#Input: 1247#previousGsq - graph G(i-1)**2 1248#currentGi - graph G(i) 1249#u and v - source and target of an edge added in this iteration 1250# 1251#Output: 1252#Graph G(i)**2 used by next steps of K-Center algorithm 1253 1254proc ::struct::graph::op::extendTwoSquaredGraph {previousGsq currentGi u v} { 1255 1256 #adding new edge 1257 if { ![$previousGsq arc exists [list $v $u]] && ![$previousGsq arc exists [list $u $v]]} { 1258 $previousGsq arc insert $u $v [list $u $v] 1259 } 1260 1261 #adding new edges to solution graph: 1262 #here edges, where source is a $u node and targets are neighbours of node $u except for $v 1263 foreach x [$currentGi nodes -adj $u] { 1264 if { ( $x != $v) && ![$previousGsq arc exists [list $v $x]] && ![$previousGsq arc exists [list $x $v]] } { 1265 $previousGsq arc insert $v $x [list $v $x] 1266 } 1267 } 1268 #here edges, where source is a $v node and targets are neighbours of node $v except for $u 1269 foreach x [$currentGi nodes -adj $v] { 1270 if { ( $x != $u ) && ![$previousGsq arc exists [list $u $x]] && ![$previousGsq arc exists [list $x $u]] } { 1271 $previousGsq arc insert $u $x [list $u $x] 1272 } 1273 } 1274 1275 return $previousGsq 1276} 1277 1278#Vertices Cover - 2 approximation algorithm 1279#------------------------------------------------------------------------------------- 1280#Vertices cover is a set o vertices such that each edge of the graph is incident to 1281#at least one vertex of the set. This 2-approximation algorithm searches for minimum 1282#vertices cover, which is a classical optimization problem in computer science and 1283#is a typical example of an NP-hard optimization problem that has an approximation 1284#algorithm. 1285# 1286#Reference: http://en.wikipedia.org/wiki/Vertex_cover_problem 1287# 1288 1289proc ::struct::graph::op::VerticesCover {G} { 1290 #variable containing final solution 1291 set vc {} 1292 #variable containing sorted (with degree) set of arcs for graph G 1293 set arcs {} 1294 1295 #setting the dictionary with degrees for each node 1296 foreach v [$G nodes] { 1297 dict set degrees $v [$G node degree $v] 1298 } 1299 1300 #creating a list containing the sum of degrees for source and 1301 #target nodes for each edge in graph G 1302 foreach e [$G arcs] { 1303 set v [$G arc source $e] 1304 set u [$G arc target $e] 1305 1306 lappend values [list [expr {[dict get $degrees $v]+[dict get $degrees $u]}] $e] 1307 } 1308 #sorting the list of source and target degrees 1309 set values [lsort -integer -decreasing -index 0 $values] 1310 1311 #setting the set of edges in a right sequence 1312 foreach e $values { 1313 lappend arcs [lindex $e 1] 1314 } 1315 1316 #for each node in graph G, we add it to the final solution and 1317 #erase all arcs adjacent to it, so they cannot be 1318 #added to solution in next iterations 1319 foreach e $arcs { 1320 1321 if { [struct::set contains $arcs $e] } { 1322 set v [$G arc source $e] 1323 set u [$G arc target $e] 1324 lappend vc $v $u 1325 1326 foreach n [$G arcs -adj $v $u] { 1327 struct::set exclude arcs $n 1328 } 1329 } 1330 } 1331 1332 return $vc 1333} 1334 1335 1336#Ford's Fulkerson algorithm - computing maximum flow in a flow network 1337#------------------------------------------------------------------------------------- 1338# 1339#The general idea of algorithm is finding augumenting paths in graph G, as long 1340#as they exist, and for each path updating the edge's weights along that path, 1341#with maximum possible throughput. The final (maximum) flow is found 1342#when there is no other augumenting path from source to sink. 1343# 1344#Input: 1345#graph G - weighted and directed graph. Weights at edges are considered as 1346#maximum throughputs that can be carried by that link (edge). 1347#s - the node that is a source for graph G 1348#t - the node that is a sink for graph G 1349# 1350#Output: 1351#Procedure returns the dictionary contaning throughputs for all edges. For 1352#each key ( the edge between nodes u and v in the for of list u v ) there is 1353#a value that is a throughput for that key. Edges where throughput values 1354#are equal to 0 are not returned ( it is like there was no link in the flow network 1355#between nodes connected by such edge). 1356# 1357#Reference: http://en.wikipedia.org/wiki/Ford-Fulkerson_algorithm 1358 1359proc ::struct::graph::op::FordFulkerson {G s t} { 1360 1361 #checking if nodes s and t are in graph G 1362 if { !([$G node exists $s] && [$G node exists $t]) } { 1363 return -code error "Nodes \"$s\" and \"$t\" should be contained in graph's G set of nodes" 1364 } 1365 1366 #checking if all attributes for input network are set well ( costs and throughputs ) 1367 foreach e [$G arcs] { 1368 if { ![$G arc keyexists $e throughput] } { 1369 return -code error "The input network doesn't have all attributes set correctly... Please, check again attributes: \"throughput\" for input graph." 1370 } 1371 } 1372 1373 #initilization 1374 foreach e [$G arcs] { 1375 set u [$G arc source $e] 1376 set v [$G arc target $e] 1377 dict set f [list $u $v] 0 1378 dict set f [list $v $u] 0 1379 } 1380 1381 #setting the residual graph for the first iteration 1382 set residualG [createResidualGraph $G $f] 1383 1384 #deleting the arcs that are 0-weighted 1385 foreach e [$residualG arcs] { 1386 if { [$residualG arc set $e throughput] == 0 } { 1387 $residualG arc delete $e 1388 } 1389 } 1390 1391 #the main loop - works till the path between source and the sink can be found 1392 while {1} { 1393 set paths [ShortestsPathsByBFS $residualG $s paths] 1394 1395 if { ($paths == {}) || (![dict exists $paths $t]) } break 1396 1397 set path [dict get $paths $t] 1398 #setting the path from source to sink 1399 1400 #adding sink to path 1401 lappend path $t 1402 1403 #finding the throughput of path p - the smallest value of c(f) among 1404 #edges that are contained in the path 1405 set maxThroughput Inf 1406 1407 foreach u [lrange $path 0 end-1] v [lrange $path 1 end] { 1408 set pathEdgeFlow [$residualG arc set [list $u $v] throughput] 1409 if { $maxThroughput > $pathEdgeFlow } { 1410 set maxThroughput $pathEdgeFlow 1411 } 1412 } 1413 1414 #increase of throughput using the path p, with value equal to maxThroughput 1415 foreach u [lrange $path 0 end-1] v [lrange $path 1 end] { 1416 1417 #if maximum throughput that was found for the path p (maxThroughput) is bigger than current throughput 1418 #at the edge not contained in the path p (for current pair of nodes u and v), then we add to the edge 1419 #which is contained into path p the maxThroughput value decreased by the value of throughput at 1420 #the second edge (not contained in path). That second edge's throughtput value is set to 0. 1421 1422 set f_uv [dict get $f [list $u $v]] 1423 set f_vu [dict get $f [list $v $u]] 1424 if { $maxThroughput >= $f_vu } { 1425 dict set f [list $u $v] [ expr { $f_uv + $maxThroughput - $f_vu } ] 1426 dict set f [list $v $u] 0 1427 } else { 1428 1429 #if maxThroughput is not greater than current throughput at the edge not contained in path p (here - v->u), 1430 #we add a difference between those values to edge contained in the path p (here u->v) and substract that 1431 #difference from edge not contained in the path p. 1432 set difference [ expr { $f_vu - $maxThroughput } ] 1433 dict set f [list $u $v] [ expr { $f_uv + $difference } ] 1434 dict set f [list $v $u] $maxThroughput 1435 } 1436 } 1437 1438 #when the current throughput for the graph is updated, we generate new residual graph 1439 #for new values of throughput 1440 $residualG destroy 1441 set residualG [createResidualGraph $G $f] 1442 1443 foreach e [$residualG arcs] { 1444 if { [$residualG arc set $e throughput] == 0 } { 1445 $residualG arc delete $e 1446 } 1447 } 1448 } 1449 1450 $residualG destroy 1451 1452 #removing 0-weighted edges from solution 1453 foreach e [dict keys $f] { 1454 if { [dict get $f $e] == 0 } { 1455 set f [dict remove $f $e] 1456 } 1457 } 1458 1459 return $f 1460} 1461 1462#subprocedure for FordFulkerson's algorithm, which creates 1463#for input graph G and given throughput f residual graph 1464#for further operations to find maximum flow in flow network 1465 1466proc ::struct::graph::op::createResidualGraph {G f} { 1467 1468 #initialization 1469 set residualG [struct::graph] 1470 1471 foreach v [$G nodes] { 1472 $residualG node insert $v 1473 } 1474 1475 foreach e [$G arcs] { 1476 set u [$G arc source $e] 1477 set v [$G arc target $e] 1478 dict set GF [list $u $v] [$G arc set $e throughput] 1479 } 1480 1481 foreach e [dict keys $GF] { 1482 1483 lassign $e u v 1484 1485 set c_uv [dict get $GF $e] 1486 set flow_uv [dict get $f $e] 1487 set flow_vu [dict get $f [list $v $u]] 1488 1489 if { ![$residualG arc exists $e] } { 1490 $residualG arc insert $u $v $e 1491 } 1492 1493 if { ![$residualG arc exists [list $v $u]] } { 1494 $residualG arc insert $v $u [list $v $u] 1495 } 1496 1497 #new value of c_f(u,v) for residual Graph is a max flow value for this edge 1498 #minus current flow on that edge 1499 if { ![$residualG arc keyexists $e throughput] } { 1500 if { [dict exists $GF [list $v $u]] } { 1501 $residualG arc set [list $u $v] throughput [ expr { $c_uv - $flow_uv + $flow_vu } ] 1502 } else { 1503 $residualG arc set $e throughput [ expr { $c_uv - $flow_uv } ] 1504 } 1505 } 1506 1507 if { [dict exists $GF [list $v $u]] } { 1508 #when double arcs in graph G (u->v , v->u) 1509 #so, x/y i w/z y-x+w 1510 set c_vu [dict get $GF [list $v $u]] 1511 if { ![$residualG arc keyexists [list $v $u] throughput] } { 1512 $residualG arc set [list $v $u] throughput [ expr { $c_vu - $flow_vu + $flow_uv} ] 1513 } 1514 } else { 1515 $residualG arc set [list $v $u] throughput $flow_uv 1516 } 1517 } 1518 1519 #setting all weights at edges to 1 for proper usage of shortest paths finding procedures 1520 $residualG arc setunweighted 1 1521 1522 return $residualG 1523} 1524 1525#Subprocedure for Busacker Gowen algorithm 1526# 1527#Input: 1528#graph G - flow network. Graph G has two attributes for each edge: 1529#cost and throughput. Each arc must have it's attribute value assigned. 1530#dictionary f - some flow for network G. Keys represent edges and values 1531#are flows at those edges 1532#path - set of nodes for which we transform the network 1533# 1534#Subprocedure checks 6 vital conditions and for them updates the network 1535#(let values with * be updates values for network). So, let edge (u,v) be 1536#the non-zero flow for network G, c(u,v) throughput of edge (u,v) and 1537#d(u,v) non-negative cost of edge (u,v): 1538#1. c*(v,u) = f(u,v) --- adding apparent arc 1539#2. d*(v,u) = -d(u,v) 1540#3. c*(u,v) = c(u,v) - f(u,v) --- if f(v,u) = 0 and c(u,v) > f(u,v) 1541#4. d*(u,v) = d(u,v) --- if f(v,u) = 0 and c(u,v) > f(u,v) 1542#5. c*(u,v) = 0 --- if f(v,u) = 0 and c(u,v) = f(u,v) 1543#6. d*(u,v) = Inf --- if f(v,u) = 0 and c(u,v) = f(u,v) 1544 1545proc ::struct::graph::op::createAugmentingNetwork {G f path} { 1546 1547 set Gf [struct::graph] 1548 1549 #setting the Gf graph 1550 foreach v [$G nodes] { 1551 $Gf node insert $v 1552 } 1553 1554 foreach e [$G arcs] { 1555 set u [$G arc source $e] 1556 set v [$G arc target $e] 1557 1558 $Gf arc insert $u $v [list $u $v] 1559 1560 $Gf arc set [list $u $v] throughput [$G arc set $e throughput] 1561 $Gf arc set [list $u $v] cost [$G arc set $e cost] 1562 } 1563 1564 #we set new values for each edge contained in the path from input 1565 foreach u [lrange $path 0 end-1] v [lrange $path 1 end] { 1566 1567 set f_uv [dict get $f [list $u $v]] 1568 set f_vu [dict get $f [list $v $u]] 1569 set c_uv [$G arc get [list $u $v] throughput] 1570 set d_uv [$G arc get [list $u $v] cost] 1571 1572 #adding apparent arcs 1573 if { ![$Gf arc exists [list $v $u]] } { 1574 $Gf arc insert $v $u [list $v $u] 1575 #1. 1576 $Gf arc set [list $v $u] throughput $f_uv 1577 #2. 1578 $Gf arc set [list $v $u] cost [ expr { -1 * $d_uv } ] 1579 } else { 1580 #1. 1581 $Gf arc set [list $v $u] throughput $f_uv 1582 #2. 1583 $Gf arc set [list $v $u] cost [ expr { -1 * $d_uv } ] 1584 $Gf arc set [list $u $v] cost Inf 1585 $Gf arc set [list $u $v] throughput 0 1586 } 1587 1588 if { ($f_vu == 0 ) && ( $c_uv > $f_uv ) } { 1589 #3. 1590 $Gf arc set [list $u $v] throughput [ expr { $c_uv - $f_uv } ] 1591 #4. 1592 $Gf arc set [list $u $v] cost $d_uv 1593 } 1594 1595 if { ($f_vu == 0 ) && ( $c_uv == $f_uv) } { 1596 #5. 1597 $Gf arc set [list $u $v] throughput 0 1598 #6. 1599 $Gf arc set [list $u $v] cost Inf 1600 } 1601 } 1602 1603 return $Gf 1604} 1605 1606#Busacker Gowen's algorithm - computing minimum cost maximum flow in a flow network 1607#------------------------------------------------------------------------------------- 1608# 1609#The goal is to find a flow, whose max value can be d, from source node to 1610#sink node in given flow network. That network except throughputs at edges has 1611#also defined a non-negative cost on each edge - cost of using that edge when 1612#directing flow with that edge ( it can illustrate e.g. fuel usage, time or 1613#any other measure dependent on usages ). 1614# 1615#Input: 1616#graph G - flow network, weights at edges are costs of using particular edge 1617#desiredFlow - max value of the flow for that network 1618#dictionary c - throughputs for all edges 1619#node s - the source node for graph G 1620#node t - the sink node for graph G 1621# 1622#Output: 1623#f - dictionary containing values of used throughputs for each edge ( key ) 1624#found by algorithm. 1625# 1626#Reference: http://en.wikipedia.org/wiki/Minimum_cost_flow_problem 1627# 1628 1629proc ::struct::graph::op::BusackerGowen {G desiredFlow s t} { 1630 1631 #checking if nodes s and t are in graph G 1632 if { !([$G node exists $s] && [$G node exists $t]) } { 1633 return -code error "Nodes \"$s\" and \"$t\" should be contained in graph's G set of nodes" 1634 } 1635 1636 if { $desiredFlow <= 0 } { 1637 return -code error "The \"desiredFlow\" value must be an positive integer." 1638 } 1639 1640 #checking if all attributes for input network are set well ( costs and throughputs ) 1641 foreach e [$G arcs] { 1642 if { !([$G arc keyexists $e throughput] && [$G arc keyexists $e cost]) } { 1643 return -code error "The input network doesn't have all attributes set correctly... Please, check again attributes: \"throughput\" and \"cost\" for input graph." 1644 } 1645 } 1646 1647 set Gf [struct::graph] 1648 1649 #initialization of Augmenting Network 1650 foreach v [$G nodes] { 1651 $Gf node insert $v 1652 } 1653 1654 foreach e [$G arcs] { 1655 set u [$G arc source $e] 1656 set v [$G arc target $e] 1657 $Gf arc insert $u $v [list $u $v] 1658 1659 $Gf arc set [list $u $v] throughput [$G arc set $e throughput] 1660 $Gf arc set [list $u $v] cost [$G arc set $e cost] 1661 } 1662 1663 #initialization of f 1664 foreach e [$G arcs] { 1665 set u [$G arc source $e] 1666 set v [$G arc target $e] 1667 dict set f [list $u $v] 0 1668 dict set f [list $v $u] 0 1669 } 1670 1671 set currentFlow 0 1672 1673 #main loop - it ends when we reach desired flow value or there is no path in Gf 1674 #leading from source node s to sink t 1675 1676 while { $currentFlow < $desiredFlow } { 1677 1678 #preparing correct values for pathfinding 1679 foreach edge [$Gf arcs] { 1680 $Gf arc setweight $edge [$Gf arc get $edge cost] 1681 } 1682 1683 #setting the path 'p' from 's' to 't' 1684 set paths [ShortestsPathsByBFS $Gf $s paths] 1685 1686 #if there are no more paths, the search has ended 1687 if { ($paths == {}) || (![dict exists $paths $t]) } break 1688 1689 set path [dict get $paths $t] 1690 lappend path $t 1691 1692 #counting max throughput that is availiable to send 1693 #using path 'p' 1694 set maxThroughput Inf 1695 foreach u [lrange $path 0 end-1] v [lrange $path 1 end] { 1696 set uv_throughput [$Gf arc set [list $u $v] throughput] 1697 if { $maxThroughput > $uv_throughput } { 1698 set maxThroughput $uv_throughput 1699 } 1700 } 1701 1702 #if max throughput that was found will cause exceeding the desired 1703 #flow, send as much as it's possible 1704 if { ( $currentFlow + $maxThroughput ) <= $desiredFlow } { 1705 set fAdd $maxThroughput 1706 set currentFlow [ expr { $currentFlow + $fAdd } ] 1707 } else { 1708 set fAdd [ expr { $desiredFlow - $currentFlow } ] 1709 set currentFlow $desiredFlow 1710 } 1711 1712 #update the throuputs on edges 1713 foreach v [lrange $path 0 end-1] u [lrange $path 1 end] { 1714 if { [dict get $f [list $u $v]] >= $fAdd } { 1715 dict set f [list $u $v] [ expr { [dict get $f [list $u $v]] - $fAdd } ] 1716 } 1717 1718 if { ( [dict get $f [list $u $v]] < $fAdd ) && ( [dict get $f [list $u $v]] > 0 ) } { 1719 dict set f [list $v $u] [ expr { $fAdd - [dict get $f [list $u $v]] } ] 1720 dict set f [list $u $v] 0 1721 } 1722 1723 if { [dict get $f [list $u $v]] == 0 } { 1724 dict set f [list $v $u] [ expr { [dict get $f [list $v $u]] + $fAdd } ] 1725 } 1726 } 1727 1728 #create new Augemnting Network 1729 1730 set Gfnew [createAugmentingNetwork $Gf $f $path] 1731 $Gf destroy 1732 set Gf $Gfnew 1733 } 1734 1735 set f [dict filter $f script {flow flowvalue} {expr {$flowvalue != 0}}] 1736 1737 $Gf destroy 1738 return $f 1739} 1740 1741# 1742proc ::struct::graph::op::ShortestsPathsByBFS {G s outputFormat} { 1743 1744 switch -exact -- $outputFormat { 1745 distances { 1746 set outputMode distances 1747 } 1748 paths { 1749 set outputMode paths 1750 } 1751 default { 1752 return -code error "Unknown output format \"$outputFormat\", expected distances, or paths." 1753 } 1754 } 1755 1756 set queue [list $s] 1757 set result {} 1758 1759 #initialization of marked nodes, distances and predecessors 1760 foreach v [$G nodes] { 1761 dict set marked $v 0 1762 dict set distances $v Inf 1763 dict set pred $v -1 1764 } 1765 1766 #the s node is initially marked and has 0 distance to itself 1767 dict set marked $s 1 1768 dict set distances $s 0 1769 1770 #the main loop 1771 while { [llength $queue] != 0 } { 1772 1773 #removing top element from the queue 1774 set v [lindex $queue 0] 1775 lremove queue $v 1776 1777 #for each arc that begins in v 1778 foreach arc [$G arcs -out $v] { 1779 1780 set u [$G arc target $arc] 1781 set newlabel [ expr { [dict get $distances $v] + [$G arc getweight $arc] } ] 1782 1783 if { $newlabel < [dict get $distances $u] } { 1784 1785 dict set distances $u $newlabel 1786 dict set pred $u $v 1787 1788 #case when current node wasn't placed in a queue yet - 1789 #we set u at the end of the queue 1790 if { [dict get $marked $u] == 0 } { 1791 lappend queue $u 1792 dict set marked $u 1 1793 } else { 1794 1795 #case when current node u was in queue before but it is not in it now - 1796 #we set u at the beginning of the queue 1797 if { [lsearch $queue $u] < 0 } { 1798 set queue [linsert $queue 0 $u] 1799 } 1800 } 1801 } 1802 } 1803 } 1804 1805 #if the outputformat is paths, we travel back to find shorests paths 1806 #to return sets of nodes for each node, which are their paths between 1807 #s and particular node 1808 dict set paths nopaths 1 1809 if { $outputMode eq "paths" } { 1810 foreach node [$G nodes] { 1811 1812 set path {} 1813 set lastNode $node 1814 1815 while { $lastNode != -1 } { 1816 set currentNode [dict get $pred $lastNode] 1817 if { $currentNode != -1 } { 1818 lappend path $currentNode 1819 } 1820 set lastNode $currentNode 1821 } 1822 1823 set path [lreverse $path] 1824 1825 if { [llength $path] != 0 } { 1826 dict set paths $node $path 1827 dict unset paths nopaths 1828 } 1829 } 1830 1831 if { ![dict exists $paths nopaths] } { 1832 return $paths 1833 } else { 1834 return {} 1835 } 1836 1837 #returning dictionary containing distance from start node to each other node (key) 1838 } else { 1839 return $distances 1840 } 1841 1842} 1843 1844# 1845proc ::struct::graph::op::BFS {G s outputFormat} { 1846 1847 set queue [list $s] 1848 1849 switch -exact -- $outputFormat { 1850 graph { 1851 set outputMode graph 1852 } 1853 tree { 1854 set outputMode tree 1855 } 1856 default { 1857 return -code error "Unknown output format \"$outputFormat\", expected graph, or tree." 1858 } 1859 } 1860 1861 if { $outputMode eq "graph" } { 1862 #graph initializing 1863 set BFSGraph [struct::graph] 1864 foreach v [$G nodes] { 1865 $BFSGraph node insert $v 1866 } 1867 } else { 1868 #tree initializing 1869 set BFSTree [struct::tree] 1870 $BFSTree set root name $s 1871 $BFSTree rename root $s 1872 } 1873 1874 #initilization of marked nodes 1875 foreach v [$G nodes] { 1876 dict set marked $v 0 1877 } 1878 1879 #start node is marked from the beginning 1880 dict set marked $s 1 1881 1882 #the main loop 1883 while { [llength $queue] != 0 } { 1884 #removing top element from the queue 1885 1886 set v [lindex $queue 0] 1887 lremove queue $v 1888 1889 foreach x [$G nodes -adj $v] { 1890 if { ![dict get $marked $x] } { 1891 dict set marked $x 1 1892 lappend queue $x 1893 1894 if { $outputMode eq "graph" } { 1895 $BFSGraph arc insert $v $x [list $v $x] 1896 } else { 1897 $BFSTree insert $v end $x 1898 } 1899 } 1900 } 1901 } 1902 1903 if { $outputMode eq "graph" } { 1904 return $BFSGraph 1905 } else { 1906 return $BFSTree 1907 } 1908} 1909 1910#Minimum Diameter Spanning Tree - MDST 1911#------------------------------------------------------------------------------------- 1912# 1913#The goal is to find for input graph G, the spanning tree that 1914#has the minimum diameter worth. 1915# 1916#General idea of algorithm is to run BFS over all vertices in graph 1917#G. If the diameter "d" of the tree is odd, then we are sure that tree 1918#given by BFS is minimum (considering diameter value). When, diameter "d" 1919#is even, then optimal tree can have minimum diameter equal to "d" or 1920#"d-1". 1921# 1922#In that case, what algorithm does is rebuilding the tree given by BFS, by 1923#adding a vertice between root node and root's child node (nodes), such that 1924#subtree created with child node as root node is the greatest one (has the 1925#greatests height). In the next step for such rebuilded tree, we run again BFS 1926#with new node as root node. If the height of the tree didn't changed, we have found 1927#a better solution. 1928 1929proc ::struct::graph::op::MinimumDiameterSpanningTree {G} { 1930 1931 set min_diameter Inf 1932 set best_Tree [struct::graph] 1933 1934 foreach v [$G nodes] { 1935 1936 #BFS Tree 1937 set T [BFS $G $v tree] 1938 #BFS Graph 1939 set TGraph [BFS $G $v graph] 1940 1941 #Setting all arcs to 1 for diameter procedure 1942 $TGraph arc setunweighted 1 1943 1944 #setting values for current Tree 1945 set diam [diameter $TGraph] 1946 set subtreeHeight [ expr { $diam / 2 - 1} ] 1947 1948 ############################################## 1949 #case when diameter found for tree found by BFS is even: 1950 #it's possible to decrease the diameter by one. 1951 if { ( $diam % 2 ) == 0 } { 1952 1953 #for each child u that current root node v has, we search 1954 #for the greatest subtree(subtrees) with the root in child u. 1955 # 1956 foreach u [$TGraph nodes -adj $v] { 1957 set u_depth 1 ;#[$T depth $u] 1958 set d_depth 0 1959 1960 set descendants [$T descendants $u] 1961 1962 foreach d $descendants { 1963 if { $d_depth < [$T depth $d] } { 1964 set d_depth [$T depth $d] 1965 } 1966 } 1967 1968 #depth of the current subtree 1969 set depth [ expr { $d_depth - $u_depth } ] 1970 1971 #proceed if found subtree is the greatest one 1972 if { $depth >= $subtreeHeight } { 1973 1974 #temporary Graph for holding potential better values 1975 set tempGraph [struct::graph] 1976 1977 foreach node [$TGraph nodes] { 1978 $tempGraph node insert $node 1979 } 1980 1981 #zmienic nazwy zmiennych zeby sie nie mylily 1982 foreach arc [$TGraph arcs] { 1983 set _u [$TGraph arc source $arc] 1984 set _v [$TGraph arc target $arc] 1985 $tempGraph arc insert $_u $_v [list $_u $_v] 1986 } 1987 1988 if { [$tempGraph arc exists [list $u $v]] } { 1989 $tempGraph arc delete [list $u $v] 1990 } else { 1991 $tempGraph arc delete [list $v $u] 1992 } 1993 1994 #for nodes u and v, we add a node between them 1995 #to again start BFS with root in new node to check 1996 #if it's possible to decrease the diameter in solution 1997 set node [$tempGraph node insert] 1998 $tempGraph arc insert $node $v [list $node $v] 1999 $tempGraph arc insert $node $u [list $node $u] 2000 2001 set newtempGraph [BFS $tempGraph $node graph] 2002 $tempGraph destroy 2003 set tempGraph $newtempGraph 2004 2005 $tempGraph node delete $node 2006 $tempGraph arc insert $u $v [list $u $v] 2007 $tempGraph arc setunweighted 1 2008 2009 set tempDiam [diameter $tempGraph] 2010 2011 #if better tree is found (that any that were already found) 2012 #replace it 2013 if { $min_diameter > $tempDiam } { 2014 set $min_diameter [diameter $tempGraph ] 2015 $best_Tree destroy 2016 set best_Tree $tempGraph 2017 } else { 2018 $tempGraph destroy 2019 } 2020 } 2021 2022 } 2023 } 2024 ################################################################ 2025 2026 set currentTreeDiameter $diam 2027 2028 if { $min_diameter > $currentTreeDiameter } { 2029 set min_diameter $currentTreeDiameter 2030 $best_Tree destroy 2031 set best_Tree $TGraph 2032 } else { 2033 $TGraph destroy 2034 } 2035 2036 $T destroy 2037 } 2038 2039 return $best_Tree 2040} 2041 2042#Minimum Degree Spanning Tree 2043#------------------------------------------------------------------------------------- 2044# 2045#In graph theory, minimum degree spanning tree (or degree-constrained spanning tree) 2046#is a spanning tree where the maximum vertex degree is as small as possible (or is 2047#limited to a certain constant k). The minimum degree spanning tree problem is to 2048#determine whether a particular graph has such a spanning tree for a particular k. 2049# 2050#Algorithm for input undirected graph G finds its spanning tree with the smallest 2051#possible degree. Algorithm is a 2-approximation, so it doesn't assure that optimal 2052#solution will be found. 2053# 2054#Reference: http://en.wikipedia.org/wiki/Degree-constrained_spanning_tree 2055 2056proc ::struct::graph::op::MinimumDegreeSpanningTree {G} { 2057 2058 #initialization of spanning tree for G 2059 set MST [struct::graph] 2060 2061 foreach v [$G nodes] { 2062 $MST node insert $v 2063 } 2064 2065 #forcing all arcs to be 1-weighted 2066 foreach e [$G arcs] { 2067 $G arc setweight $e 1 2068 } 2069 2070 foreach e [kruskal $G] { 2071 set u [$G arc source $e] 2072 set v [$G arc target $e] 2073 2074 $MST arc insert $u $v [list $u $v] 2075 } 2076 2077 #main loop 2078 foreach e [$G arcs] { 2079 2080 set u [$G arc source $e] 2081 set v [$G arc target $e] 2082 2083 #if nodes u and v are neighbours, proceed to next iteration 2084 if { ![$MST arc exists [list $u $v]] && ![$MST arc exists [list $v $u]] } { 2085 2086 $MST arc setunweighted 1 2087 2088 #setting the path between nodes u and v in Spanning Tree MST 2089 set path [dict get [dijkstra $MST $u] $v] 2090 lappend path $v 2091 2092 #search for the node in the path, such that its degree is greater than degree of any of nodes 2093 #u or v increased by one 2094 foreach node $path { 2095 if { [$MST node degree $node] > ([Max [$MST node degree $u] [$MST node degree $v]] + 1) } { 2096 2097 #if such node is found add the arc between nodes u and v 2098 $MST arc insert $u $v [list $u $v] 2099 2100 #then to hold MST being a spanning tree, delete any arc that is in the path 2101 #that is adjacent to found node 2102 foreach n [$MST nodes -adj $node] { 2103 if { $n in $path } { 2104 if { [$MST arc exists [list $node $n]] } { 2105 $MST arc delete [list $node $n] 2106 } else { 2107 $MST arc delete [list $n $node] 2108 } 2109 break 2110 } 2111 } 2112 2113 # Node found, stop processing the path 2114 break 2115 } 2116 } 2117 } 2118 } 2119 2120 return $MST 2121} 2122 2123#Dinic algorithm for finding maximum flow in flow network 2124#------------------------------------------------------------------------------------- 2125# 2126#Reference: http://en.wikipedia.org/wiki/Dinic's_algorithm 2127# 2128proc ::struct::graph::op::MaximumFlowByDinic {G s t blockingFlowAlg} { 2129 2130 if { !($blockingFlowAlg eq "dinic" || $blockingFlowAlg eq "mkm") } { 2131 return -code error "Uncorrect name of blocking flow algorithm. Choose \"mkm\" for Malhotra, Kumar and Maheshwari algorithm and \"dinic\" for Dinic algorithm." 2132 } 2133 2134 foreach arc [$G arcs] { 2135 set u [$G arc source $arc] 2136 set v [$G arc target $arc] 2137 2138 dict set f [list $u $v] 0 2139 dict set f [list $v $u] 0 2140 } 2141 2142 while {1} { 2143 set residualG [createResidualGraph $G $f] 2144 if { $blockingFlowAlg == "mkm" } { 2145 set blockingFlow [BlockingFlowByMKM $residualG $s $t] 2146 } else { 2147 set blockingFlow [BlockingFlowByDinic $residualG $s $t] 2148 } 2149 $residualG destroy 2150 2151 if { $blockingFlow == {} } break 2152 2153 foreach key [dict keys $blockingFlow] { 2154 dict set f $key [ expr { [dict get $f $key] + [dict get $blockingFlow $key] } ] 2155 } 2156 } 2157 2158 set f [dict filter $f script {flow flowvalue} {expr {$flowvalue != 0}}] 2159 2160 return $f 2161} 2162 2163#Dinic algorithm for finding blocking flow 2164#------------------------------------------------------------------------------------- 2165# 2166#Algorithm for given network G with source s and sink t, finds a blocking 2167#flow, which can be used to obtain a maximum flow for that network G. 2168# 2169#Some steps that algorithm takes: 2170#1. constructing the level graph from network G 2171#2. until there are edges in level graph: 2172# 3. find the path between s and t nodes in level graph 2173# 4. for each edge in path update current throughputs at those edges and... 2174# 5. ...deleting nodes from which there are no residual edges 2175#6. return the dictionary containing the blocking flow 2176 2177proc ::struct::graph::op::BlockingFlowByDinic {G s t} { 2178 2179 #initializing blocking flow dictionary 2180 foreach edge [$G arcs] { 2181 set u [$G arc source $edge] 2182 set v [$G arc target $edge] 2183 2184 dict set b [list $u $v] 0 2185 } 2186 2187 #1. 2188 set LevelGraph [createLevelGraph $G $s] 2189 2190 #2. the main loop 2191 while { [llength [$LevelGraph arcs]] > 0 } { 2192 2193 if { ![$LevelGraph node exists $s] || ![$LevelGraph node exists $t] } break 2194 2195 #3. 2196 set paths [ShortestsPathsByBFS $LevelGraph $s paths] 2197 2198 if { $paths == {} } break 2199 if { ![dict exists $paths $t] } break 2200 2201 set path [dict get $paths $t] 2202 lappend path $t 2203 2204 #setting the max throughput to go with the path found one step before 2205 set maxThroughput Inf 2206 foreach u [lrange $path 0 end-1] v [lrange $path 1 end] { 2207 2208 set uv_throughput [$LevelGraph arc get [list $u $v] throughput] 2209 2210 if { $maxThroughput > $uv_throughput } { 2211 set maxThroughput $uv_throughput 2212 } 2213 } 2214 2215 #4. updating throughputs and blocking flow 2216 foreach u [lrange $path 0 end-1] v [lrange $path 1 end] { 2217 2218 set uv_throughput [$LevelGraph arc get [list $u $v] throughput] 2219 #decreasing the throughputs contained in the path by max flow value 2220 $LevelGraph arc set [list $u $v] throughput [ expr { $uv_throughput - $maxThroughput } ] 2221 2222 #updating blocking flows 2223 dict set b [list $u $v] [ expr { [dict get $b [list $u $v]] + $maxThroughput } ] 2224 #dict set b [list $v $u] [ expr { -1 * [dict get $b [list $u $v]] } ] 2225 2226 #5. deleting the arcs, whose throughput is completely used 2227 if { [$LevelGraph arc get [list $u $v] throughput] == 0 } { 2228 $LevelGraph arc delete [list $u $v] 2229 } 2230 2231 #deleting the node, if it hasn't any outgoing arcs 2232 if { ($u != $s) && ( ![llength [$LevelGraph nodes -out $u]] || ![llength [$LevelGraph nodes -in $u]] ) } { 2233 $LevelGraph node delete $u 2234 } 2235 } 2236 2237 } 2238 2239 set b [dict filter $b script {flow flowvalue} {expr {$flowvalue != 0}}] 2240 2241 $LevelGraph destroy 2242 2243 #6. 2244 return $b 2245} 2246 2247#Malhotra, Kumar and Maheshwari Algorithm for finding blocking flow 2248#------------------------------------------------------------------------------------- 2249# 2250#Algorithm for given network G with source s and sink t, finds a blocking 2251#flow, which can be used to obtain a maximum flow for that network G. 2252# 2253#For given node v, Let c(v) be the min{ a, b }, where a is the sum of all incoming 2254#throughputs and b is the sum of all outcoming throughputs from the node v. 2255# 2256#Some steps that algorithm takes: 2257#1. constructing the level graph from network G 2258#2. until there are edges in level graph: 2259# 3. finding the node with the minimum c(v) 2260# 4. sending c(v) units of throughput by incoming arcs of v 2261# 5. sending c(v) units of throughput by outcoming arcs of v 2262# 6. 4 and 5 steps can cause excess or deficiency of throughputs at nodes, so we 2263# send exceeds forward choosing arcs greedily and... 2264# 7. ...the same with deficiencies but we send those backward. 2265# 8. delete the v node from level graph 2266# 9. upgrade the c values for all nodes 2267# 2268#10. if no other edges left in level graph, return b - found blocking flow 2269# 2270 2271proc ::struct::graph::op::BlockingFlowByMKM {G s t} { 2272 2273 #initializing blocking flow dictionary 2274 foreach edge [$G arcs] { 2275 set u [$G arc source $edge] 2276 set v [$G arc target $edge] 2277 2278 dict set b [list $u $v] 0 2279 } 2280 2281 #1. setting the level graph 2282 set LevelGraph [createLevelGraph $G $s] 2283 2284 #setting the in/out throughputs for each node 2285 set c [countThroughputsAtNodes $LevelGraph $s $t] 2286 2287 #2. the main loop 2288 while { [llength [$LevelGraph nodes]] > 2 } { 2289 2290 #if there is no path between s and t nodes, end the procedure and 2291 #return current blocking flow 2292 set distances [ShortestsPathsByBFS $LevelGraph $s distances] 2293 if { [dict get $distances $t] == "Inf" } { 2294 $LevelGraph destroy 2295 set b [dict filter $b script {flow flowvalue} {expr {$flowvalue != 0}}] 2296 return $b 2297 } 2298 2299 #3. finding the node with minimum value of c(v) 2300 set min_cv Inf 2301 2302 dict for {node cv} $c { 2303 if { $min_cv > $cv } { 2304 set min_cv $cv 2305 set minCv_node $node 2306 } 2307 } 2308 2309 #4. sending c(v) by all incoming arcs of node with minimum c(v) 2310 set _min_cv $min_cv 2311 foreach arc [$LevelGraph arcs -in $minCv_node] { 2312 2313 set t_arc [$LevelGraph arc get $arc throughput] 2314 set u [$LevelGraph arc source $arc] 2315 set v [$LevelGraph arc target $arc] 2316 set b_uv [dict get $b [list $u $v]] 2317 2318 if { $t_arc >= $min_cv } { 2319 $LevelGraph arc set $arc throughput [ expr { $t_arc - $min_cv } ] 2320 dict set b [list $u $v] [ expr { $b_uv + $min_cv } ] 2321 break 2322 } else { 2323 set difference [ expr { $min_cv - $t_arc } ] 2324 set min_cv $difference 2325 dict set b [list $u $v] [ expr { $b_uv + $difference } ] 2326 $LevelGraph arc set $arc throughput 0 2327 } 2328 } 2329 2330 #5. sending c(v) by all outcoming arcs of node with minimum c(v) 2331 foreach arc [$LevelGraph arcs -out $minCv_node] { 2332 2333 set t_arc [$LevelGraph arc get $arc throughput] 2334 set u [$LevelGraph arc source $arc] 2335 set v [$LevelGraph arc target $arc] 2336 set b_uv [dict get $b [list $u $v]] 2337 2338 if { $t_arc >= $min_cv } { 2339 $LevelGraph arc set $arc throughput [ expr { $t_arc - $_min_cv } ] 2340 dict set b [list $u $v] [ expr { $b_uv + $_min_cv } ] 2341 break 2342 } else { 2343 set difference [ expr { $_min_cv - $t_arc } ] 2344 set _min_cv $difference 2345 dict set b [list $u $v] [ expr { $b_uv + $difference } ] 2346 $LevelGraph arc set $arc throughput 0 2347 } 2348 } 2349 2350 #find exceeds and if any, send them forward or backwards 2351 set distances [ShortestsPathsByBFS $LevelGraph $s distances] 2352 2353 #6. 2354 for {set i [ expr {[dict get $distances $minCv_node] + 1}] } { $i < [llength [$G nodes]] } { incr i } { 2355 foreach w [$LevelGraph nodes] { 2356 if { [dict get $distances $w] == $i } { 2357 set excess [findExcess $LevelGraph $w $b] 2358 if { $excess > 0 } { 2359 set b [sendForward $LevelGraph $w $b $excess] 2360 } 2361 } 2362 } 2363 } 2364 2365 #7. 2366 for { set i [ expr { [dict get $distances $minCv_node] - 1} ] } { $i > 0 } { incr i -1 } { 2367 foreach w [$LevelGraph nodes] { 2368 if { [dict get $distances $w] == $i } { 2369 set excess [findExcess $LevelGraph $w $b] 2370 if { $excess < 0 } { 2371 set b [sendBack $LevelGraph $w $b [ expr { (-1) * $excess } ]] 2372 } 2373 } 2374 } 2375 } 2376 2377 #8. delete current node from the network 2378 $LevelGraph node delete $minCv_node 2379 2380 #9. correctingg the in/out throughputs for each node after 2381 #deleting one of the nodes in network 2382 set c [countThroughputsAtNodes $LevelGraph $s $t] 2383 2384 #if node has no availiable outcoming or incoming throughput 2385 #delete that node from the graph 2386 dict for {key val} $c { 2387 if { $val == 0 } { 2388 $LevelGraph node delete $key 2389 dict unset c $key 2390 } 2391 } 2392 } 2393 2394 set b [dict filter $b script {flow flowvalue} {expr {$flowvalue != 0}}] 2395 2396 $LevelGraph destroy 2397 #10. 2398 return $b 2399} 2400 2401#Subprocedure for algorithms that find blocking-flows. 2402#It's creating a level graph from the residual network. 2403proc ::struct::graph::op::createLevelGraph {Gf s} { 2404 2405 set LevelGraph [struct::graph] 2406 2407 $Gf arc setunweighted 1 2408 2409 #deleting arcs with 0 throughputs for proper pathfinding 2410 foreach arc [$Gf arcs] { 2411 if { [$Gf arc get $arc throughput] == 0 } { 2412 $Gf arc delete $arc 2413 } 2414 } 2415 2416 set distances [ShortestsPathsByBFS $Gf $s distances] 2417 2418 foreach v [$Gf nodes] { 2419 $LevelGraph node insert $v 2420 $LevelGraph node set $v distance [dict get $distances $v] 2421 } 2422 2423 foreach e [$Gf arcs] { 2424 set u [$Gf arc source $e] 2425 set v [$Gf arc target $e] 2426 2427 if { ([$LevelGraph node get $u distance] + 1) == [$LevelGraph node get $v distance]} { 2428 $LevelGraph arc insert $u $v [list $u $v] 2429 $LevelGraph arc set [list $u $v] throughput [$Gf arc get $e throughput] 2430 } 2431 } 2432 2433 $LevelGraph arc setunweighted 1 2434 return $LevelGraph 2435} 2436 2437#Subprocedure for blocking flow finding by MKM algorithm 2438# 2439#It computes for graph G and each of his nodes the throughput value - 2440#for node v: from the sum of availiable throughputs from incoming arcs and 2441#the sum of availiable throughputs from outcoming arcs chooses lesser and sets 2442#as the throughput of the node. 2443# 2444#Throughputs of nodes are returned in the dictionary. 2445# 2446proc ::struct::graph::op::countThroughputsAtNodes {G s t} { 2447 2448 set c {} 2449 foreach v [$G nodes] { 2450 2451 if { ($v eq $t) || ($v eq $s) } continue 2452 2453 set outcoming [$G arcs -out $v] 2454 set incoming [$G arcs -in $v] 2455 2456 set outsum 0 2457 set insum 0 2458 2459 foreach o $outcoming i $incoming { 2460 2461 if { [llength $o] > 0 } { 2462 set outsum [ expr { $outsum + [$G arc get $o throughput] } ] 2463 } 2464 2465 if { [llength $i] > 0 } { 2466 set insum [ expr { $insum + [$G arc get $i throughput] } ] 2467 } 2468 2469 set value [Min $outsum $insum] 2470 } 2471 2472 dict set c $v $value 2473 } 2474 2475 return $c 2476} 2477 2478#Subprocedure for blocking-flow finding algorithm by MKM 2479# 2480#If for a given input node, outcoming flow is bigger than incoming, then that deficiency 2481#has to be send back by that subprocedure. 2482proc ::struct::graph::op::sendBack {G node b value} { 2483 2484 foreach arc [$G arcs -in $node] { 2485 set u [$G arc source $arc] 2486 set v [$G arc target $arc] 2487 2488 if { $value > [$G arc get $arc throughput] } { 2489 set value [ expr { $value - [$G arc get $arc throughput] } ] 2490 dict set b [list $u $v] [ expr { [dict get $b [list $u $v]] + [$G arc get $arc throughput] } ] 2491 $G arc set $arc throughput 0 2492 } else { 2493 $G arc set $arc throughput [ expr { [$G arc get $arc throughput] - $value } ] 2494 dict set b [list $u $v] [ expr { [dict get $b [list $u $v]] + $value } ] 2495 set value 0 2496 break 2497 } 2498 } 2499 2500 return $b 2501} 2502 2503#Subprocedure for blocking-flow finding algorithm by MKM 2504# 2505#If for a given input node, incoming flow is bigger than outcoming, then that exceed 2506#has to be send forward by that sub procedure. 2507proc ::struct::graph::op::sendForward {G node b value} { 2508 2509 foreach arc [$G arcs -out $node] { 2510 2511 set u [$G arc source $arc] 2512 set v [$G arc target $arc] 2513 2514 if { $value > [$G arc get $arc throughput] } { 2515 set value [ expr { $value - [$G arc get $arc throughput] } ] 2516 dict set b [list $u $v] [ expr { [dict get $b [list $u $v]] + [$G arc get $arc throughput] } ] 2517 $G arc set $arc throughput 0 2518 } else { 2519 $G arc set $arc throughput [ expr { [$G arc get $arc throughput] - $value } ] 2520 dict set b [list $u $v] [ expr { [dict get $b [list $u $v]] + $value } ] 2521 2522 set value 0 2523 break 2524 } 2525 } 2526 2527 return $b 2528} 2529 2530#Subprocedure for blocking-flow finding algorithm by MKM 2531# 2532#It checks for graph G if node given at input has a exceed 2533#or deficiency of throughput. 2534# 2535#For exceed the positive value of exceed is returned, for deficiency 2536#procedure returns negative value. If the incoming throughput 2537#is the same as outcoming, procedure returns 0. 2538# 2539proc ::struct::graph::op::findExcess {G node b} { 2540 2541 set incoming 0 2542 set outcoming 0 2543 2544 foreach key [dict keys $b] { 2545 2546 lassign $key u v 2547 if { $u eq $node } { 2548 set outcoming [ expr { $outcoming + [dict get $b $key] } ] 2549 } 2550 if { $v eq $node } { 2551 set incoming [ expr { $incoming + [dict get $b $key] } ] 2552 } 2553 } 2554 2555 return [ expr { $incoming - $outcoming } ] 2556} 2557 2558#Travelling Salesman Problem - Heuristic of local searching 2559#2 - approximation Algorithm 2560#------------------------------------------------------------------------------------- 2561# 2562 2563proc ::struct::graph::op::TSPLocalSearching {G C} { 2564 2565 foreach arc $C { 2566 if { ![$G arc exists $arc] } { 2567 return -code error "Given cycle has arcs not included in graph G." 2568 } 2569 } 2570 2571 #initialization 2572 set CGraph [struct::graph] 2573 set GCopy [struct::graph] 2574 set w 0 2575 2576 foreach node [$G nodes] { 2577 $CGraph node insert $node 2578 $GCopy node insert $node 2579 } 2580 2581 foreach arc [$G arcs] { 2582 set u [$G arc source $arc] 2583 set v [$G arc target $arc] 2584 $GCopy arc insert $u $v [list $u $v] 2585 $GCopy arc set [list $u $v] weight [$G arc get $arc weight] 2586 } 2587 2588 foreach arc $C { 2589 2590 set u [$G arc source $arc] 2591 set v [$G arc target $arc] 2592 set arcWeight [$G arc get $arc weight] 2593 2594 $CGraph arc insert $u $v [list $u $v] 2595 $CGraph arc set [list $u $v] weight $arcWeight 2596 2597 set w [ expr { $w + $arcWeight } ] 2598 } 2599 2600 set reductionDone 1 2601 2602 while { $reductionDone } { 2603 2604 set queue {} 2605 set reductionDone 0 2606 2607 #double foreach loop goes through all pairs of arcs 2608 foreach i [$CGraph arcs] { 2609 2610 #source and target nodes of first arc 2611 set iu [$CGraph arc source $i] 2612 set iv [$CGraph arc target $i] 2613 2614 #second arc 2615 foreach j [$CGraph arcs] { 2616 2617 #if pair of arcs already was considered, continue with next pair of arcs 2618 if { [list $j $i] ni $queue } { 2619 2620 #add current arc to queue to mark that it was used 2621 lappend queue [list $i $j] 2622 2623 set ju [$CGraph arc source $j] 2624 set jv [$CGraph arc target $j] 2625 2626 #we consider only arcs that are not adjacent 2627 if { !($iu eq $ju) && !($iu eq $jv) && !($iv eq $ju) && !($iv eq $jv) } { 2628 2629 #set the current cycle 2630 set CPrim [copyGraph $CGraph] 2631 2632 #transform the current cycle: 2633 #1. 2634 $CPrim arc delete $i 2635 $CPrim arc delete $j 2636 2637 2638 set param 0 2639 2640 #adding new edges instead of erased ones 2641 if { !([$CPrim arc exists [list $iu $ju]] || [$CPrim arc exists [list $iv $jv]] || [$CPrim arc exists [list $ju $iu]] || [$CPrim arc exists [list $jv $iv]] ) } { 2642 2643 $CPrim arc insert $iu $ju [list $iu $ju] 2644 $CPrim arc insert $iv $jv [list $iv $jv] 2645 2646 if { [$GCopy arc exists [list $iu $ju]] } { 2647 $CPrim arc set [list $iu $ju] weight [$GCopy arc get [list $iu $ju] weight] 2648 } else { 2649 $CPrim arc set [list $iu $ju] weight [$GCopy arc get [list $ju $iu] weight] 2650 } 2651 2652 if { [$GCopy arc exists [list $iv $jv]] } { 2653 $CPrim arc set [list $iv $jv] weight [$GCopy arc get [list $iv $jv] weight] 2654 } else { 2655 $CPrim arc set [list $iv $jv] weight [$GCopy arc get [list $jv $iv] weight] 2656 } 2657 } else { 2658 set param 1 2659 } 2660 2661 $CPrim arc setunweighted 1 2662 2663 #check if it's still a cycle or if any arcs were added instead those erased 2664 if { !([struct::graph::op::distance $CPrim $iu $ju] > 0 ) || $param } { 2665 2666 #deleting new edges if they were added before in current iteration 2667 if { !$param } { 2668 $CPrim arc delete [list $iu $ju] 2669 } 2670 2671 if { !$param } { 2672 $CPrim arc delete [list $iv $jv] 2673 } 2674 2675 #adding new ones that will assure the graph is still a cycle 2676 $CPrim arc insert $iu $jv [list $iu $jv] 2677 $CPrim arc insert $iv $ju [list $iv $ju] 2678 2679 if { [$GCopy arc exists [list $iu $jv]] } { 2680 $CPrim arc set [list $iu $jv] weight [$GCopy arc get [list $iu $jv] weight] 2681 } else { 2682 $CPrim arc set [list $iu $jv] weight [$GCopy arc get [list $jv $iu] weight] 2683 } 2684 2685 if { [$GCopy arc exists [list $iv $ju]] } { 2686 $CPrim arc set [list $iv $ju] weight [$GCopy arc get [list $iv $ju] weight] 2687 } else { 2688 $CPrim arc set [list $iv $ju] weight [$GCopy arc get [list $ju $iv] weight] 2689 } 2690 } 2691 2692 #count current value of cycle 2693 set cycleWeight [countCycleWeight $CPrim] 2694 2695 #if we found cycle with lesser sum of weights, we set is as a result and 2696 #marked that reduction was successful 2697 if { $w > $cycleWeight } { 2698 set w $cycleWeight 2699 set reductionDone 1 2700 set C [$CPrim arcs] 2701 } 2702 2703 $CPrim destroy 2704 } 2705 } 2706 } 2707 } 2708 2709 #setting the new current cycle if the reduction was successful 2710 if { $reductionDone } { 2711 foreach arc [$CGraph arcs] { 2712 $CGraph arc delete $arc 2713 } 2714 for {set i 0} { $i < [llength $C] } { incr i } { 2715 lset C $i [lsort [lindex $C $i]] 2716 } 2717 2718 foreach arc [$GCopy arcs] { 2719 if { [lsort $arc] in $C } { 2720 set u [$GCopy arc source $arc] 2721 set v [$GCopy arc target $arc] 2722 $CGraph arc insert $u $v [list $u $v] 2723 $CGraph arc set $arc weight [$GCopy arc get $arc weight] 2724 } 2725 } 2726 } 2727 } 2728 2729 $GCopy destroy 2730 $CGraph destroy 2731 2732 return $C 2733} 2734 2735proc ::struct::graph::op::copyGraph {G} { 2736 2737 set newGraph [struct::graph] 2738 2739 foreach node [$G nodes] { 2740 $newGraph node insert $node 2741 } 2742 foreach arc [$G arcs] { 2743 set u [$G arc source $arc] 2744 set v [$G arc target $arc] 2745 $newGraph arc insert $u $v $arc 2746 $newGraph arc set $arc weight [$G arc get $arc weight] 2747 } 2748 2749 return $newGraph 2750} 2751 2752proc ::struct::graph::op::countCycleWeight {G} { 2753 2754 set result 0 2755 2756 foreach arc [$G arcs] { 2757 set result [ expr { $result + [$G arc get $arc weight] } ] 2758 } 2759 2760 return $result 2761} 2762 2763# ### ### ### ######### ######### ######### 2764## 2765 2766# This command finds a minimum spanning tree/forest (MST) of the graph 2767# argument, using the algorithm developed by Joseph Kruskal. The 2768# result is a set (as list) containing the names of the arcs in the 2769# MST. The set of nodes of the MST is implied by set of arcs, and thus 2770# not given explicitly. The algorithm does not consider arc 2771# directions. Note that unconnected nodes are left out of the result. 2772 2773# Reference: http://en.wikipedia.org/wiki/Kruskal%27s_algorithm 2774 2775proc ::struct::graph::op::kruskal {g} { 2776 # Check graph argument for proper configuration. 2777 2778 VerifyWeightsAreOk $g 2779 2780 # Transient helper data structures. A priority queue for the arcs 2781 # under consideration, using their weights as priority, and a 2782 # disjoint-set to keep track of the forest of partial minimum 2783 # spanning trees we are working with. 2784 2785 set consider [::struct::prioqueue -dictionary consider] 2786 set forest [::struct::disjointset forest] 2787 2788 # Start with all nodes in the graph each in their partition. 2789 2790 foreach n [$g nodes] { 2791 $forest add-partition $n 2792 } 2793 2794 # Then fill the queue with all arcs, using their weight to 2795 # prioritize. The weight is the cost of the arc. The lesser the 2796 # better. 2797 2798 foreach {arc weight} [$g arc weights] { 2799 $consider put $arc $weight 2800 } 2801 2802 # And now we can construct the tree. This is done greedily. In 2803 # each round we add the arc with the smallest weight to the 2804 # minimum spanning tree, except if doing so would violate the tree 2805 # condition. 2806 2807 set result {} 2808 2809 while {[$consider size]} { 2810 set minarc [$consider get] 2811 set origin [$g arc source $minarc] 2812 set destin [$g arc target $minarc] 2813 2814 # Ignore the arc if both ends are in the same partition. Using 2815 # it would add a cycle to the result, i.e. it would not be a 2816 # tree anymore. 2817 2818 if {[$forest equal $origin $destin]} continue 2819 2820 # Take the arc for the result, and merge the trees both ends 2821 # are in into a single tree. 2822 2823 lappend result $minarc 2824 $forest merge $origin $destin 2825 } 2826 2827 # We are done. Get rid of the transient helper structures and 2828 # return our result. 2829 2830 $forest destroy 2831 $consider destroy 2832 2833 return $result 2834} 2835 2836# ### ### ### ######### ######### ######### 2837## 2838 2839# This command finds a minimum spanning tree/forest (MST) of the graph 2840# argument, using the algorithm developed by Prim. The result is a 2841# set (as list) containing the names of the arcs in the MST. The set 2842# of nodes of the MST is implied by set of arcs, and thus not given 2843# explicitly. The algorithm does not consider arc directions. 2844 2845# Reference: http://en.wikipedia.org/wiki/Prim%27s_algorithm 2846 2847proc ::struct::graph::op::prim {g} { 2848 VerifyWeightsAreOk $g 2849 2850 # Fill an array with all nodes, to track which nodes have been 2851 # visited at least once. When the inner loop runs out of nodes and 2852 # we still have some left over we restart using one of the 2853 # leftover as new starting point. In this manner we get the MST of 2854 # the whole graph minus unconnected nodes, instead of only the MST 2855 # for the component the initial starting node is in. 2856 2857 array set unvisited {} 2858 foreach n [$g nodes] { set unvisited($n) . } 2859 2860 # Transient helper data structure. A priority queue for the nodes 2861 # and arcs under consideration for inclusion into the MST. Each 2862 # element of the queue is a list containing node name, a flag bit, 2863 # and arc name, in this order. The associated priority is the 2864 # weight of the arc. The flag bit is set for the initial queue 2865 # entry only, containing a fake (empty) arc, to trigger special 2866 # handling. 2867 2868 set consider [::struct::prioqueue -dictionary consider] 2869 2870 # More data structures, the result arrays. 2871 array set weightmap {} ; # maps nodes to min arc weight seen so 2872 # far. This is the threshold other arcs 2873 # on this node will have to beat to be 2874 # added to the MST. 2875 array set arcmap {} ; # maps arcs to nothing, these are the 2876 # arcs in the MST. 2877 2878 while {[array size unvisited]} { 2879 # Choose a 'random' node as the starting point for the inner 2880 # loop, prim's algorithm, and put it on the queue for 2881 # consideration. Then we iterate until we have considered all 2882 # nodes in the its component. 2883 2884 set startnode [lindex [array names unvisited] 0] 2885 $consider put [list $startnode 1 {}] 0 2886 2887 while {[$consider size] > 0} { 2888 # Pull the next minimum weight to look for. This is the 2889 # priority of the next item we can get from the queue. And the 2890 # associated node/decision/arc data. 2891 2892 set arcweight [$consider peekpriority 1] 2893 2894 foreach {v arcundefined arc} [$consider get] break 2895 #8.5: lassign [$consider get] v arcundefined arc 2896 2897 # Two cases to consider: The node v is already part of the 2898 # MST, or not. If yes we check if the new arcweight is better 2899 # than what we have stored already, and update accordingly. 2900 2901 if {[info exists weightmap($v)]} { 2902 set currentweight $weightmap($v) 2903 if {$arcweight < $currentweight} { 2904 # The new weight is better, update to use it as 2905 # the new threshold. Note that this fill not touch 2906 # any other arcs found for this node, as these are 2907 # still minimal. 2908 2909 set weightmap($v) $arcweight 2910 set arcmap($arc) . 2911 } 2912 } else { 2913 # Node not yet present. Save weight and arc. The 2914 # latter if and only the arc is actually defined. For 2915 # the first, initial queue entry, it is not. Then we 2916 # add all the arcs adjacent to the current node to the 2917 # queue to consider them in the next rounds. 2918 2919 set weightmap($v) $arcweight 2920 if {!$arcundefined} { 2921 set arcmap($arc) . 2922 } 2923 foreach adjacentarc [$g arcs -adj $v] { 2924 set weight [$g arc getweight $adjacentarc] 2925 set neighbour [$g node opposite $v $adjacentarc] 2926 $consider put [list $neighbour 0 $adjacentarc] $weight 2927 } 2928 } 2929 2930 # Mark the node as visited, belonging to the current 2931 # component. Future iterations will ignore it. 2932 unset -nocomplain unvisited($v) 2933 } 2934 } 2935 2936 # We are done. Get rid of the transient helper structure and 2937 # return our result. 2938 2939 $consider destroy 2940 2941 return [array names arcmap] 2942} 2943 2944# ### ### ### ######### ######### ######### 2945## 2946 2947# This command checks whether the graph argument is bi-partite or not, 2948# and returns the result as a boolean value, true for a bi-partite 2949# graph, and false otherwise. A variable can be provided to store the 2950# bi-partition into. 2951# 2952# Reference: http://en.wikipedia.org/wiki/Bipartite_graph 2953 2954proc ::struct::graph::op::isBipartite? {g {bipartitionvar {}}} { 2955 2956 # Handle the special cases of empty graphs, or one without arcs 2957 # quickly. Both are bi-partite. 2958 2959 if {$bipartitionvar ne ""} { 2960 upvar 1 $bipartitionvar bipartitions 2961 } 2962 if {![llength [$g nodes]]} { 2963 set bipartitions {{} {}} 2964 return 1 2965 } elseif {![llength [$g arcs]]} { 2966 if {$bipartitionvar ne ""} { 2967 set bipartitions [list [$g nodes] {}] 2968 } 2969 return 1 2970 } 2971 2972 # Transient helper data structure, a queue of the nodes waiting 2973 # for processing. 2974 2975 set pending [struct::queue pending] 2976 set nodes [$g nodes] 2977 2978 # Another structure, a map from node names to their 'color', 2979 # indicating which of the two partitions a node belngs to. All 2980 # nodes start out as undefined (0). Traversing the arcs we 2981 # set and flip them as needed (1,2). 2982 2983 array set color {} 2984 foreach node $nodes { 2985 set color($node) 0 2986 } 2987 2988 # Iterating over all nodes we use their connections to traverse 2989 # the components and assign colors. We abort when encountering 2990 # paradox, as that means that the graph is not bi-partite. 2991 2992 foreach node $nodes { 2993 # Ignore nodes already in the second partition. 2994 if {$color($node)} continue 2995 2996 # Flip the color, then travel the component and check for 2997 # conflicts with the neighbours. 2998 2999 set color($node) 1 3000 3001 $pending put $node 3002 while {[$pending size]} { 3003 set current [$pending get] 3004 foreach neighbour [$g nodes -adj $current] { 3005 if {!$color($neighbour)} { 3006 # Exchange the color between current and previous 3007 # nodes, and remember the neighbour for further 3008 # processing. 3009 set color($neighbour) [expr {3 - $color($current)}] 3010 $pending put $neighbour 3011 } elseif {$color($neighbour) == $color($current)} { 3012 # Color conflict between adjacent nodes, should be 3013 # different. This graph is not bi-partite. Kill 3014 # the data structure and abort. 3015 3016 $pending destroy 3017 return 0 3018 } 3019 } 3020 } 3021 } 3022 3023 # The graph is bi-partite. Kill the transient data structure, and 3024 # move the partitions into the provided variable, if there is any. 3025 3026 $pending destroy 3027 3028 if {$bipartitionvar ne ""} { 3029 # Build bipartition, then set the data into the variable 3030 # passed as argument to this command. 3031 3032 set X {} 3033 set Y {} 3034 3035 foreach {node partition} [array get color] { 3036 if {$partition == 1} { 3037 lappend X $node 3038 } else { 3039 lappend Y $node 3040 } 3041 } 3042 set bipartitions [list $X $Y] 3043 } 3044 3045 return 1 3046} 3047 3048# ### ### ### ######### ######### ######### 3049## 3050 3051# This command computes a maximal matching, if it exists, for the 3052# graph argument G and its bi-partition as specified through the node 3053# sets X and Y. As is implied, this method requires that the graph is 3054# bi-partite. Use the command 'isBipartite?' to check for this 3055# property, and to obtain the bi-partition. 3056if 0 { 3057 proc ::struct::graph::op::maxMatching {g X Y} { 3058 return -code error "not implemented yet" 3059 }} 3060 3061# ### ### ### ######### ######### ######### 3062## 3063 3064# This command computes the strongly connected components (SCCs) of 3065# the graph argument G. The result is a list of node-sets, each set 3066# containing the nodes of one SCC of G. In any SCC there is a directed 3067# path between any two nodes U, V from U to V. If all SCCs contain 3068# only a single node the graph is acyclic. 3069 3070proc ::struct::graph::op::tarjan {g} { 3071 set all [$g nodes] 3072 3073 # Quick bailout for simple special cases, i.e. graphs without 3074 # nodes or arcs. 3075 if {![llength $all]} { 3076 # No nodes => no SCCs 3077 return {} 3078 } elseif {![llength [$g arcs]]} { 3079 # Have nodes, but no arcs => each node is its own SCC. 3080 set r {} ; foreach a $all { lappend r [list $a] } 3081 return $r 3082 } 3083 3084 # Transient data structures. Stack of nodes to consider, the 3085 # result, and various state arrays. TarjanSub upvar's all them 3086 # into its scope. 3087 3088 set pending [::struct::stack pending] 3089 set result {} 3090 3091 array set index {} 3092 array set lowlink {} 3093 array set instack {} 3094 3095 # Invoke the main search system while we have unvisited 3096 # nodes. TarjanSub will remove all visited nodes from 'all', 3097 # ensuring termination. 3098 3099 while {[llength $all]} { 3100 TarjanSub [lindex $all 0] 0 3101 } 3102 3103 # Release the transient structures and return result. 3104 $pending destroy 3105 return $result 3106} 3107 3108proc ::struct::graph::op::TarjanSub {start counter} { 3109 # Import the tracer state from our caller. 3110 upvar 1 g g index index lowlink lowlink instack instack result result pending pending all all 3111 3112 struct::set subtract all $start 3113 3114 set component {} 3115 set index($start) $counter 3116 set lowlink($start) $counter 3117 incr counter 3118 3119 $pending push $start 3120 set instack($start) 1 3121 3122 foreach outarc [$g arcs -out $start] { 3123 set neighbour [$g arc target $outarc] 3124 3125 if {![info exists index($neighbour)]} { 3126 # depth-first-search of reachable nodes from the neighbour 3127 # node. Original from the chosen startnode. 3128 TarjanSub $neighbour $counter 3129 set lowlink($start) [Min $lowlink($start) $lowlink($neighbour)] 3130 3131 } elseif {[info exists instack($neighbour)]} { 3132 set lowlink($start) [Min $lowlink($start) $lowlink($neighbour)] 3133 } 3134 } 3135 3136 # Check if the 'start' node on this recursion level is the root 3137 # node of a SCC, and collect the component if yes. 3138 3139 if {$lowlink($start) == $index($start)} { 3140 while {1} { 3141 set v [$pending pop] 3142 unset instack($v) 3143 lappend component $v 3144 if {$v eq $start} break 3145 } 3146 lappend result $component 3147 } 3148 3149 return 3150} 3151 3152# ### ### ### ######### ######### ######### 3153## 3154 3155# This command computes the connected components (CCs) of the graph 3156# argument G. The result is a list of node-sets, each set containing 3157# the nodes of one CC of G. In any CC there is UN-directed path 3158# between any two nodes U, V. 3159 3160proc ::struct::graph::op::connectedComponents {g} { 3161 set all [$g nodes] 3162 3163 # Quick bailout for simple special cases, i.e. graphs without 3164 # nodes or arcs. 3165 if {![llength $all]} { 3166 # No nodes => no CCs 3167 return {} 3168 } elseif {![llength [$g arcs]]} { 3169 # Have nodes, but no arcs => each node is its own CC. 3170 set r {} ; foreach a $all { lappend r [list $a] } 3171 return $r 3172 } 3173 3174 # Invoke the main search system while we have unvisited 3175 # nodes. 3176 3177 set result {} 3178 while {[llength $all]} { 3179 set component [ComponentOf $g [lindex $all 0]] 3180 lappend result $component 3181 # all = all - component 3182 struct::set subtract all $component 3183 } 3184 return $result 3185} 3186 3187# A derivative command which computes the connected component (CC) of 3188# the graph argument G containing the node N. The result is a node-set 3189# containing the nodes of the CC of N in G. 3190 3191proc ::struct::graph::op::connectedComponentOf {g n} { 3192 # Quick bailout for simple special cases 3193 if {![$g node exists $n]} { 3194 return -code error "node \"$n\" does not exist in graph \"$g\"" 3195 } elseif {![llength [$g arcs -adj $n]]} { 3196 # The chosen node has no neighbours, so is its own CC. 3197 return [list $n] 3198 } 3199 3200 # Invoke the main search system for the chosen node. 3201 3202 return [ComponentOf $g $n] 3203} 3204 3205# Internal helper for finding connected components. 3206 3207proc ::struct::graph::op::ComponentOf {g start} { 3208 set pending [::struct::queue pending] 3209 $pending put $start 3210 3211 array set visited {} 3212 set visited($start) . 3213 3214 while {[$pending size]} { 3215 set current [$pending get 1] 3216 foreach neighbour [$g nodes -adj $current] { 3217 if {[info exists visited($neighbour)]} continue 3218 $pending put $neighbour 3219 set visited($neighbour) 1 3220 } 3221 } 3222 $pending destroy 3223 return [array names visited] 3224} 3225 3226# ### ### ### ######### ######### ######### 3227## 3228 3229# This command determines if the specified arc A in the graph G is a 3230# bridge, i.e. if its removal will split the connected component its 3231# end nodes belong to, into two. The result is a boolean value. Uses 3232# the 'ComponentOf' helper command. 3233 3234proc ::struct::graph::op::isBridge? {g arc} { 3235 if {![$g arc exists $arc]} { 3236 return -code error "arc \"$arc\" does not exist in graph \"$g\"" 3237 } 3238 3239 # Note: We could avoid the need for a copy of the graph if we were 3240 # willing to modify G (*). As we are not willing using a copy is 3241 # the easiest way to allow us a trivial modification. For the 3242 # future consider the creation of a graph class which represents 3243 # virtual graphs over a source, generated by deleting nodes and/or 3244 # arcs. without actually modifying the source. 3245 # 3246 # (Ad *): Create a new unnamed helper node X. Move the arc 3247 # destination to X. Recompute the component and ignore 3248 # X. Then move the arc target back to its original node 3249 # and remove X again. 3250 3251 set src [$g arc source $arc] 3252 set compBefore [ComponentOf $g $src] 3253 if {[llength $compBefore] == 1} { 3254 # Special case, the arc is a loop on an otherwise unconnected 3255 # node. The component will not split, this is not a bridge. 3256 return 0 3257 } 3258 3259 set copy [struct::graph BridgeCopy = $g] 3260 $copy arc delete $arc 3261 set compAfter [ComponentOf $copy $src] 3262 $copy destroy 3263 3264 return [expr {[llength $compBefore] != [llength $compAfter]}] 3265} 3266 3267# This command determines if the specified node N in the graph G is a 3268# cut vertex, i.e. if its removal will split the connected component 3269# it belongs to into two. The result is a boolean value. Uses the 3270# 'ComponentOf' helper command. 3271 3272proc ::struct::graph::op::isCutVertex? {g n} { 3273 if {![$g node exists $n]} { 3274 return -code error "node \"$n\" does not exist in graph \"$g\"" 3275 } 3276 3277 # Note: We could avoid the need for a copy of the graph if we were 3278 # willing to modify G (*). As we are not willing using a copy is 3279 # the easiest way to allow us a trivial modification. For the 3280 # future consider the creation of a graph class which represents 3281 # virtual graphs over a source, generated by deleting nodes and/or 3282 # arcs. without actually modifying the source. 3283 # 3284 # (Ad *): Create two new unnamed helper nodes X and Y. Move the 3285 # icoming and outgoing arcs to these helpers. Recompute 3286 # the component and ignore the helpers. Then move the arcs 3287 # back to their original nodes and remove the helpers 3288 # again. 3289 3290 set compBefore [ComponentOf $g $n] 3291 3292 if {[llength $compBefore] == 1} { 3293 # Special case. The node is unconnected. Its removal will 3294 # cause no changes. Therefore not a cutvertex. 3295 return 0 3296 } 3297 3298 # We remove the node from the original component, so that we can 3299 # select a new start node without fear of hitting on the 3300 # cut-vertex candidate. Also makes the comparison later easier 3301 # (straight ==). 3302 struct::set subtract compBefore $n 3303 3304 set copy [struct::graph CutVertexCopy = $g] 3305 $copy node delete $n 3306 set compAfter [ComponentOf $copy [lindex $compBefore 0]] 3307 $copy destroy 3308 3309 return [expr {[llength $compBefore] != [llength $compAfter]}] 3310} 3311 3312# This command determines if the graph G is connected. 3313 3314proc ::struct::graph::op::isConnected? {g} { 3315 return [expr { [llength [connectedComponents $g]] == 1 }] 3316} 3317 3318# ### ### ### ######### ######### ######### 3319## 3320 3321# This command determines if the specified graph G has an eulerian 3322# cycle (aka euler tour, <=> g is eulerian) or not. If yes, it can 3323# return the cycle through the named variable, as a list of arcs 3324# traversed. 3325# 3326# Note that for a graph to be eulerian all nodes have to have an even 3327# degree, and the graph has to be connected. And if more than two 3328# nodes have an odd degree the graph is not even semi-eulerian (cannot 3329# even have an euler path). 3330 3331proc ::struct::graph::op::isEulerian? {g {eulervar {}} {tourstart {}}} { 3332 set nodes [$g nodes] 3333 if {![llength $nodes] || ![llength [$g arcs]]} { 3334 # Quick bailout for special cases. No nodes, or no arcs imply 3335 # that no euler cycle is present. 3336 return 0 3337 } 3338 3339 # Check the condition regarding even degree nodes, then 3340 # connected-ness. 3341 3342 foreach n $nodes { 3343 if {([$g node degree $n] % 2) == 0} continue 3344 # Odd degree node found, not eulerian. 3345 return 0 3346 } 3347 3348 if {![isConnected? $g]} { 3349 return 0 3350 } 3351 3352 # At this point the graph is connected, with all nodes of even 3353 # degree. As per Carl Hierholzer the graph has to have an euler 3354 # tour. If the user doesn't request it we do not waste the time to 3355 # actually compute one. 3356 3357 if {$tourstart ne ""} { 3358 upvar 1 $tourstart start 3359 } 3360 3361 # We start the tour at an arbitrary node. 3362 set start [lindex $nodes 0] 3363 3364 if {$eulervar eq ""} { 3365 return 1 3366 } 3367 3368 upvar 1 $eulervar tour 3369 Fleury $g $start tour 3370 return 1 3371} 3372 3373# This command determines if the specified graph G has an eulerian 3374# path (<=> g is semi-eulerian) or not. If yes, it can return the 3375# path through the named variable, as a list of arcs traversed. 3376# 3377# (*) Aka euler tour. 3378# 3379# Note that for a graph to be semi-eulerian at most two nodes are 3380# allowed to have an odd degree, all others have to be of even degree, 3381# and the graph has to be connected. 3382 3383proc ::struct::graph::op::isSemiEulerian? {g {eulervar {}}} { 3384 set nodes [$g nodes] 3385 if {![llength $nodes] || ![llength [$g arcs]]} { 3386 # Quick bailout for special cases. No nodes, or no arcs imply 3387 # that no euler path is present. 3388 return 0 3389 } 3390 3391 # Check the condition regarding oddd/even degree nodes, then 3392 # connected-ness. 3393 3394 set odd 0 3395 foreach n $nodes { 3396 if {([$g node degree $n] % 2) == 0} continue 3397 incr odd 3398 set lastodd $n 3399 } 3400 if {($odd > 2) || ![isConnected? $g]} { 3401 return 0 3402 } 3403 3404 # At this point the graph is connected, with the node degrees 3405 # supporting existence of an euler path. If the user doesn't 3406 # request it we do not waste the time to actually compute one. 3407 3408 if {$eulervar eq ""} { 3409 return 1 3410 } 3411 3412 upvar 1 $eulervar path 3413 3414 # We start at either an odd-degree node, or any node, if there are 3415 # no odd-degree ones. In the last case we are actually 3416 # constructing an euler tour, i.e. a closed path. 3417 3418 if {$odd} { 3419 set start $lastodd 3420 } else { 3421 set start [lindex $nodes 0] 3422 } 3423 3424 Fleury $g $start path 3425 return 1 3426} 3427 3428proc ::struct::graph::op::Fleury {g start eulervar} { 3429 upvar 1 $eulervar path 3430 3431 # We start at the chosen node. 3432 3433 set copy [struct::graph FleuryCopy = $g] 3434 set path {} 3435 3436 # Edges are chosen per Fleury's algorithm. That is easy, 3437 # especially as we already have a command to determine whether an 3438 # arc is a bridge or not. 3439 3440 set arcs [$copy arcs] 3441 while {![struct::set empty $arcs]} { 3442 set adjacent [$copy arcs -adj $start] 3443 3444 if {[llength $adjacent] == 1} { 3445 # No choice in what arc to traverse. 3446 set arc [lindex $adjacent 0] 3447 } else { 3448 # Choose first non-bridge arcs. The euler conditions force 3449 # that at least two such are present. 3450 3451 set has 0 3452 foreach arc $adjacent { 3453 if {[isBridge? $copy $arc]} { 3454 continue 3455 } 3456 set has 1 3457 break 3458 } 3459 if {!$has} { 3460 $copy destroy 3461 return -code error {Internal error} 3462 } 3463 } 3464 3465 set start [$copy node opposite $start $arc] 3466 $copy arc delete $arc 3467 struct::set exclude arcs $arc 3468 lappend path $arc 3469 } 3470 3471 $copy destroy 3472 return 3473} 3474 3475# ### ### ### ######### ######### ######### 3476## 3477 3478# This command uses dijkstra's algorithm to find all shortest paths in 3479# the graph G starting at node N. The operation can be configured to 3480# traverse arcs directed and undirected, and the format of the result. 3481 3482proc ::struct::graph::op::dijkstra {g node args} { 3483 # Default traversal is undirected. 3484 # Default output format is tree. 3485 3486 set arcTraversal undirected 3487 set resultFormat tree 3488 3489 # Process options to override the defaults, if any. 3490 foreach {option param} $args { 3491 switch -exact -- $option { 3492 -arcmode { 3493 switch -exact -- $param { 3494 directed - 3495 undirected { 3496 set arcTraversal $param 3497 } 3498 default { 3499 return -code error "Bad value for -arcmode, expected one of \"directed\" or \"undirected\"" 3500 } 3501 } 3502 } 3503 -outputformat { 3504 switch -exact -- $param { 3505 tree - 3506 distances { 3507 set resultFormat $param 3508 } 3509 default { 3510 return -code error "Bad value for -outputformat, expected one of \"distances\" or \"tree\"" 3511 } 3512 } 3513 } 3514 default { 3515 return -code error "Bad option \"$option\", expected one of \"-arcmode\" or \"-outputformat\"" 3516 } 3517 } 3518 } 3519 3520 # We expect that all arcs of g are given a weight. 3521 VerifyWeightsAreOk $g 3522 3523 # And the start node has to belong to the graph too, of course. 3524 if {![$g node exists $node]} { 3525 return -code error "node \"$node\" does not exist in graph \"$g\"" 3526 } 3527 3528 # TODO: Quick bailout for special cases (no arcs). 3529 3530 # Transient and other data structures for the core algorithm. 3531 set pending [::struct::prioqueue -dictionary DijkstraQueue] 3532 array set distance {} ; # array: node -> distance to 'n' 3533 array set previous {} ; # array: node -> parent in shortest path to 'n'. 3534 array set visited {} ; # array: node -> bool, true when node processed 3535 3536 # Initialize the data structures. 3537 foreach n [$g nodes] { 3538 set distance($n) Inf 3539 set previous($n) undefined 3540 set visited($n) 0 3541 } 3542 3543 # Compute the distances ... 3544 $pending put $node 0 3545 set distance($node) 0 3546 set previous($node) none 3547 3548 while {[$pending size]} { 3549 set current [$pending get] 3550 set visited($current) 1 3551 3552 # Traversal to neighbours according to the chosen mode. 3553 if {$arcTraversal eq "undirected"} { 3554 set arcNeighbours [$g arcs -adj $current] 3555 } else { 3556 set arcNeighbours [$g arcs -out $current] 3557 } 3558 3559 # Compute distances, record newly discovered nodes, minimize 3560 # distances for nodes reachable through multiple paths. 3561 foreach arcNeighbour $arcNeighbours { 3562 set cost [$g arc getweight $arcNeighbour] 3563 set neighbour [$g node opposite $current $arcNeighbour] 3564 set delta [expr {$distance($current) + $cost}] 3565 3566 if { 3567 ($distance($neighbour) eq "Inf") || 3568 ($delta < $distance($neighbour)) 3569 } { 3570 # First path, or better path to the node folund, 3571 # update our records. 3572 3573 set distance($neighbour) $delta 3574 set previous($neighbour) $current 3575 if {!$visited($neighbour)} { 3576 $pending put $neighbour $delta 3577 } 3578 } 3579 } 3580 } 3581 3582 $pending destroy 3583 3584 # Now generate the result based on the chosen format. 3585 if {$resultFormat eq "distances"} { 3586 return [array get distance] 3587 } else { 3588 array set listofprevious {} 3589 foreach n [$g nodes] { 3590 set current $n 3591 while {1} { 3592 if {$current eq "undefined"} break 3593 if {$current eq $node} { 3594 lappend listofprevious($n) $current 3595 break 3596 } 3597 if {$current ne $n} { 3598 lappend listofprevious($n) $current 3599 } 3600 set current $previous($current) 3601 } 3602 } 3603 return [array get listofprevious] 3604 } 3605} 3606 3607# This convenience command is a wrapper around dijkstra's algorithm to 3608# find the (un)directed distance between two nodes in the graph G. 3609 3610proc ::struct::graph::op::distance {g origin destination args} { 3611 if {![$g node exists $origin]} { 3612 return -code error "node \"$origin\" does not exist in graph \"$g\"" 3613 } 3614 if {![$g node exists $destination]} { 3615 return -code error "node \"$destination\" does not exist in graph \"$g\"" 3616 } 3617 3618 set arcTraversal undirected 3619 3620 # Process options to override the defaults, if any. 3621 foreach {option param} $args { 3622 switch -exact -- $option { 3623 -arcmode { 3624 switch -exact -- $param { 3625 directed - 3626 undirected { 3627 set arcTraversal $param 3628 } 3629 default { 3630 return -code error "Bad value for -arcmode, expected one of \"directed\" or \"undirected\"" 3631 } 3632 } 3633 } 3634 default { 3635 return -code error "Bad option \"$option\", expected \"-arcmode\"" 3636 } 3637 } 3638 } 3639 3640 # Quick bailout for special case: the distance from a node to 3641 # itself is zero 3642 3643 if {$origin eq $destination} { 3644 return 0 3645 } 3646 3647 # Compute all distances, then pick and return the one we are 3648 # interested in. 3649 array set distance [dijkstra $g $origin -outputformat distances -arcmode $arcTraversal] 3650 return $distance($destination) 3651} 3652 3653# This convenience command is a wrapper around dijkstra's algorithm to 3654# find the (un)directed eccentricity of the node N in the graph G. The 3655# eccentricity is the maximal distance to any other node in the graph. 3656 3657proc ::struct::graph::op::eccentricity {g node args} { 3658 if {![$g node exists $node]} { 3659 return -code error "node \"$node\" does not exist in graph \"$g\"" 3660 } 3661 3662 set arcTraversal undirected 3663 3664 # Process options to override the defaults, if any. 3665 foreach {option param} $args { 3666 switch -exact -- $option { 3667 -arcmode { 3668 switch -exact -- $param { 3669 directed - 3670 undirected { 3671 set arcTraversal $param 3672 } 3673 default { 3674 return -code error "Bad value for -arcmode, expected one of \"directed\" or \"undirected\"" 3675 } 3676 } 3677 } 3678 default { 3679 return -code error "Bad option \"$option\", expected \"-arcmode\"" 3680 } 3681 } 3682 } 3683 3684 # Compute all distances, then pick out the max 3685 3686 set ecc 0 3687 foreach {n distance} [dijkstra $g $node -outputformat distances -arcmode $arcTraversal] { 3688 if {$distance eq "Inf"} { return Inf } 3689 if {$distance > $ecc} { set ecc $distance } 3690 } 3691 3692 return $ecc 3693} 3694 3695# This convenience command is a wrapper around eccentricity to find 3696# the (un)directed radius of the graph G. The radius is the minimal 3697# eccentricity over all nodes in the graph. 3698 3699proc ::struct::graph::op::radius {g args} { 3700 return [lindex [RD $g $args] 0] 3701} 3702 3703# This convenience command is a wrapper around eccentricity to find 3704# the (un)directed diameter of the graph G. The diameter is the 3705# maximal eccentricity over all nodes in the graph. 3706 3707proc ::struct::graph::op::diameter {g args} { 3708 return [lindex [RD $g $args] 1] 3709} 3710 3711proc ::struct::graph::op::RD {g options} { 3712 set arcTraversal undirected 3713 3714 # Process options to override the defaults, if any. 3715 foreach {option param} $options { 3716 switch -exact -- $option { 3717 -arcmode { 3718 switch -exact -- $param { 3719 directed - 3720 undirected { 3721 set arcTraversal $param 3722 } 3723 default { 3724 return -code error "Bad value for -arcmode, expected one of \"directed\" or \"undirected\"" 3725 } 3726 } 3727 } 3728 default { 3729 return -code error "Bad option \"$option\", expected \"-arcmode\"" 3730 } 3731 } 3732 } 3733 3734 set radius Inf 3735 set diameter 0 3736 foreach n [$g nodes] { 3737 set e [eccentricity $g $n -arcmode $arcTraversal] 3738 #puts "$n ==> ($e)" 3739 if {($e eq "Inf") || ($e > $diameter)} { 3740 set diameter $e 3741 } 3742 if {($radius eq "Inf") || ($e < $radius)} { 3743 set radius $e 3744 } 3745 } 3746 3747 return [list $radius $diameter] 3748} 3749 3750# 3751## place holder for operations to come 3752# 3753 3754# ### ### ### ######### ######### ######### 3755## Internal helpers 3756 3757proc ::struct::graph::op::Min {first second} { 3758 if {$first > $second} { 3759 return $second 3760 } else { 3761 return $first 3762 } 3763} 3764 3765proc ::struct::graph::op::Max {first second} { 3766 if {$first < $second} { 3767 return $second 3768 } else { 3769 return $first 3770 } 3771} 3772 3773# This method verifies that every arc on the graph has a weight 3774# assigned to it. This is required for some algorithms. 3775proc ::struct::graph::op::VerifyWeightsAreOk {g} { 3776 if {![llength [$g arc getunweighted]]} return 3777 return -code error "Operation invalid for graph with unweighted arcs." 3778} 3779 3780# ### ### ### ######### ######### ######### 3781## Ready 3782 3783namespace eval ::struct::graph::op { 3784 #namespace export ... 3785} 3786 3787package provide struct::graph::op 0.11.3 3788